How to Tune ARIMA Parameters in Python
There are many parameters to consider when configuring an ARIMA model with Statsmodels in Python.
In this tutorial, we take a look at a few key parameters (other than the order parameter) that you may be curious about.
Specifically, after completing this tutorial, you will know:
- How to suppress noisy output from the underlying mathematical libraries when fitting an ARIMA model.
- The effect of enabling or disabling a trend term in your ARIMA model.
- The influence of using different mathematical solvers to fit coefficients to your training data.
Note, if you are interested in tuning the order parameter, see the post:
Let’s get started.
Shampoo Sales Dataset
This dataset describes the monthly number of sales of shampoo over a 3 year period.
The units are a sales count and there are 36 observations. The original dataset is credited to Makridakis, Wheelwright, and Hyndman (1998).
The example below loads and creates a plot of the loaded dataset.
| 12345678910111213 | # load and plot datasetfrom pandas import read_csvfrom pandas import datetimefrom matplotlib import pyplot# load datasetdef parser(x):returndatetime.strptime('190'+x,'%Y-%m')series=read_csv('shampoo-sales.csv',header=0,parse_dates=[0],index_col=0,squeeze=True,date_parser=parser)# summarize first few rowsprint(series.head())# line plotseries.plot()pyplot.show() |
Running the example loads the dataset as a Pandas Series and prints the first 5 rows.
| 1234567 | Month1901-01-01 266.01901-02-01 145.91901-03-01 183.11901-04-01 119.31901-05-01 180.3Name: Sales, dtype: float64 |
A line plot of the series is then created showing a clear increasing trend.
Line Plot of Monthly Shampoo Sales Dataset
Experimental Test-Setup
It is important to evaluate time series forecasting models consistently.
In this section, we will define how we will evaluate the three forecast models in this tutorial.
First, we will hold the last one year of data back and evaluate forecasts on this data. Given the data is monthly, this means that the last 12 observations will be used as test data.
We will use a walk-forward validation method to evaluate model performance. This means that each time step in the test dataset will be enumerated, a model constructed on history data, and the forecast compared to the expected value. The observation will then be added to the training dataset and the process repeated.
Walk-forward validation is a realistic way to evaluate time series forecast models as one would expect models to be updated as new observations are made available.
Finally, forecasts will be evaluated using root mean squared error, or RMSE. The benefit of RMSE is that it penalizes large errors and the scores are in the same units as the forecast values (car sales per month).
An ARIMA(4,1,0) forecast model will be used as the baseline to explore the additional parameters of the model. This may not be the optimal model for the problem, but is generally skillful against some other hand tested configurations.
In summary, the test harness involves:
- The last 2 years of data used a test set.
- Walk-forward validation for model evaluation.
- Root mean squared error used to report model skill.
- An ARIMA(4,1,0) model will be used as a baseline.
The complete example is listed below.
| 1234567891011121314151617181920212223242526272829303132 | from pandas import read_csvfrom pandas import datetimefrom matplotlib import pyplotfrom statsmodels.tsa.arima_model import ARIMAfrom sklearn.metrics import mean_squared_errorfrom math import sqrt# load datasetdef parser(x):returndatetime.strptime('190'+x,'%Y-%m')series=read_csv('shampoo-sales.csv',header=0,parse_dates=[0],index_col=0,squeeze=True,date_parser=parser)# split into train and test setsX=series.valuestrain,test=X[0:-12],X[-12:]history=[xforxintrain]predictions=list()# walk-forward validationfortinrange(len(test)):# fit modelmodel=ARIMA(history,order=(4,1,0))model_fit=model.fit()# one step forecastyhat=model_fit.forecast()[0]# store forecast and obpredictions.append(yhat)history.append(test[t])# evaluate forecastsrmse=sqrt(mean_squared_error(test,predictions))print('Test RMSE: %.3f'%rmse)# plot forecasts against actual outcomespyplot.plot(test)pyplot.plot(predictions,color='red')pyplot.show() |
Running the example spews a lot of convergence information and finishes with an RMSE score of 84.832 monthly shampoo sales.
| 123456789101112131415161718192021222324 | ...Tit = total number of iterationsTnf = total number of function evaluationsTnint = total number of segments explored during Cauchy searchesSkip = number of BFGS updates skippedNact = number of active bounds at final generalized Cauchy pointProjg = norm of the final projected gradientF = final function value * * * N Tit Tnf Tnint Skip Nact Projg F 5 15 20 1 0 0 8.882D-08 5.597D+00 F = 5.5972342395324288CONVERGENCE: REL_REDUCTION_OF_F_<=_FACTR*EPSMCH Cauchy time 0.000E+00 seconds. Subspace minimization time 0.000E+00 seconds. Line search time 0.000E+00 seconds. Total User time 0.000E+00 seconds.Test RMSE: 84.832 |
A plot of the forecast vs the actual observations in the test harness is created to give some context for the model we are working with.
ARIMA Forecast for Monthly Shampoo Sales Dataset
Now let’s dive into some of the other ARIMA parameters.
The “disp” Parameter
The first parameter we will look at is the disp parameter.
This is described as follows:
If True, convergence information is printed. For the default l_bfgs_b solver, disp controls the frequency of the output during the iterations. disp < 0 means no output in this case.
By default, this parameter is set to 1, which shows output.
We are dealing with this first because it is critical in removing all of the convergence output when evaluating the ARIMA model using walk-forward validation.
Setting it to False turns off all of this noise.
The complete example is listed below.
| 1234567891011121314151617181920212223242526272829 | from pandas import read_csvfrom pandas import datetimefrom matplotlib import pyplotfrom statsmodels.tsa.arima_model import ARIMAfrom sklearn.metrics import mean_squared_errorfrom math import sqrt# load datasetdef parser(x):returndatetime.strptime('190'+x,'%Y-%m')series=read_csv('shampoo-sales.csv',header=0,parse_dates=[0],index_col=0,squeeze=True,date_parser=parser)# split into train and test setsX=series.valuessize=int(len(X)*0.66)train,test=X[0:size],X[size:len(X)]history=[xforxintrain]predictions=list()# walk-forward validationfortinrange(len(test)):# fit modelmodel=ARIMA(history,order=(4,1,0))model_fit=model.fit(disp=False)# one step forecastyhat=model_fit.forecast()[0]# store forecast and obpredictions.append(yhat)history.append(test[t])# evaluate forecastsrmse=sqrt(mean_squared_error(test,predictions))print('Test RMSE: %.3f'%rmse) |
Running this example not only produces a cleaner output, but also is much faster to execute.
| 1 | Test RMSE: 81.545 |
We will leave disp=False on all following examples.
The “transparams” Parameter
This parameter controls whether or not to perform a transform on AR parameters.
Specifically, it is described as:
Whether or not to transform the parameters to ensure stationarity. Uses the transformation suggested in Jones (1980). If False, no checking for stationarity or invertibility is done.
By default, transparams is set to True, meaning this transform is performed.
This parameter is also used on the R version of the ARIMA implementation (see docs) and I expect this is why it is here in statsmodels.
The statsmodels doco is weak on this, but you can learn more about the transform in the paper:
The example below demonstrates turning this parameter off.
| 1234567891011121314151617181920212223242526272829 | from pandas import read_csvfrom pandas import datetimefrom matplotlib import pyplotfrom statsmodels.tsa.arima_model import ARIMAfrom sklearn.metrics import mean_squared_errorfrom math import sqrt# load datasetdef parser(x):returndatetime.strptime('190'+x,'%Y-%m')series=read_csv('shampoo-sales.csv',header=0,parse_dates=[0],index_col=0,squeeze=True,date_parser=parser)# split into train and test setsX=series.valuessize=int(len(X)*0.66)train,test=X[0:size],X[size:len(X)]history=[xforxintrain]predictions=list()# walk-forward validationfortinrange(len(test)):# fit modelmodel=ARIMA(history,order=(4,1,0))model_fit=model.fit(disp=False,transparams=False |