SARIMAX model fitting too slow in statsmodels

Question

I am trying a grid search to perform model selection by fitting SARIMAX(p, d, q)x(P, D, Q, s) models using SARIMAX() method in statsmodels. I do set d and D to 1 and s to 7 and iterate over values of p in {0, 1}, q in {0, 1, 2}, P in {0, 1}, Q in {0, 1}, and trend in {None, 'c'}, which makes for a total of 48 iterations. During the model fitting phase, if any combination of the parameters leads to a non-stationary or non-invertible model, I move to the next combination of parameters.

I have a set of time-series, each one representing the performance of an agent over time and consisting of 83 (daily) measurements with a weekly seasonality. I keep 90% of the data for model fitting, and the last 10% for forecasting/testing purposes.

What I find is that model fitting during the grid search takes a very long time, about 11 minutes, for a couple of agents, whereas the same 48 iterations take much less time, less than 10 seconds, for others.

However, if, before performing my grid search, I log-transform the data corresponding to the agents whose analyses take a very long time, the same 48 iterations take about 15 seconds! However, as much as I love the speed-up factor, the final forecast turns out to be poorer compared to the case where the original (that is, not log-transformed) data was used. So, I'd rather keep the data in its original format.

My questions are the following: What causes such slow down for certain time-serires? And is there a way to speed-up the model fitting by giving SARIMAX() or SARIMAX.fit() certain arguments? I have tried simple_differencing = True which, constructing a smaller model in the state-space, reduced the time from 11 minutes to 6 minutes, but that's still too long.

I'd appreciate any help.

Seventeen

If the slowdown is only for certain time series, then they are most likely difficult to optimize. You could try to find better starting values for the optimization, i.e. start_params. For example, when trying different models it's often possible to warm start, by adding zeros for extra parameters or dropping parameters. – Tello

@user333700 Thanks a lot! I have tried start_params, by first performing a burn-in .fit() and then fitting the final model using .fit(params = <burn_in_model>.params). The problem is that I get a speed-up IF I set enforce_stationarity and enforce_invertibility to False for the burn-in phase (the final model fit is performed with these two parameters set to True). I don't know if (A) my final estimates are biased by not enforcing stationarity and invertiblity in the burn-in phase and if (B) the burn-in parameters would lead to the best fit possible (given the data and the model). – Seventeen

@user333700 But I don't quite understand the very last part of your comment, "by adding zeros for extra parameters or dropping parameters." Could you elaborate a bit please? To what exactly should I be adding zeros or which parameters should I be dropping? Maybe a one-liner example could help me. Much appreciated :) – Seventeen

suppose you have ARMA(1,1) and you want to estimate ARMA(2,1), then it has one additional parameter. So we could take the ARMA(1,1) parameter vector and insert a zero at the right place for the second AR coefficient. It requires to check how params combines the different coefficients. (I don't remember for SARIMAX) – Tello

@user333700 Ah, I see! Thank you. I went back to using the burn-in phase without enforcing stationarity and invertibility while using simple_differencing = True. It cut down the run time from ~ 11 minutes down to 14 seconds! The best-fit models are largely the same and have similar AICc scores. I'll stick with this for now, but will try your suggestion when I get a bit more free time to implement it. Cheers! – Seventeen

UPDATE: Not enforcing stationarity and invertibility within the a burn-in step sometimes leads to the optimizer not being able to find/converge to a value in the final .fit(). I'd get Warning: Desired error not necessarily achieved due to precision loss. when method = 'bfgs' is used, and the nm algorithm exceeds the maximum number of iterations. – Seventeen

I realise it's been more than a year, but I got the same errors as you when I tried these same two methods. I was able to get nm to converge by increasing the max iter to a ridiculously high value. – Wallsend

this library is just useless. I'll just use R for time series. – Cerebration

With using simple_differencing = True, Forecasts and predictions will be about the differenced data. Any work around for this @Seventeen ? – Target