ARMA Model

Example: ARMA Model in TSMT

Autoregressive models take the general form:

 

 

 

 

The time series module provides a number of routines for performing pre-estimation data analysis, model parameter estimation, and post-estimation diagnosis of autoregressive time series.

Preliminary Autocorrelations and Partial Autocorrelations

Suppose we have stationary univariate time series data but are uncertain about the order of the regressive relationship. You may wish to first use the sample autocorrelation function can help identify the general structure of the model. This can be achieved in Gauss using the autocormt and autocovmt functions.

The autocor function computes the autocorrelations for each column of a matrix. It assumes that all data has been demeaned and produces a output matrix containing the autocorrelations for each column of an input matrix. The autocor function requires three inputs: (1) an N x K datamatrix, x, (2) a scalar, f, denoting the first autocorrelation to be computed, and (3) a scalar, l, denoting the last autocorrelation to be computed.

Consider an autoregressive time series randomly generated using the Gauss ARMA data simulation function. For example, consider a time series that follows the purely autoregressive data generation process given by

 

 

 

 

Realizations of any ARMA data generation process can be randomly created using the function simarmamt ,

//SPECIFY ARMA PARAMETERS
b_sim =  {0.5,-0.8};

//SPECIFY MA ORDER
q_sim = 0;

//SPECIFY AR ORDER
p_sim = 2;

//SPECIFY DETERMINISTIC FEATURES (CONSTANT AND TREND)
const_sim=1.5;

//NUMBER OF OBSERVATIONS
n_sim=200;

//NUMBER OF SERIES
k_sim=4;

//STANDARD DEVIATION OF ERROR TERMS
std_sim=1;

seed_sim = 5012;

y_sim = simarmamt(b_sim,p_sim,q_sim,const_sim,n_sim,k_sim,std_sim,seed_sim);

//VISUALIZE DATA
x_plot=seqa(1,1,rows(y_sim));

//PLOT Y OBS
plotXY(x_plot,y_sim);

The code above will produce a 200X4 matrix of data that follows the AR(2) data generating process specified. Each column of the matrix represents a different random realization of the process. The code graphs the first realization:

Though the data generating process of the data above is known, it is generally more realistic the one will wish to estimate the parameters of an unknown data generating process behind a realization of data. To begin consider finding the ACF and the PACF the autocor and pacf  Gauss functions ,

//FIND AUTOCORRELATIONS FOR 0 to 10 lags

//DEMEAN DATA

plotSetTitle(&myPlot, "Demeaned Data" );

y = y_sim - meanc(y_sim)';

plotXY(myPlot,x_plot,y_sim);

//FIND AUTOCORRELATIONS
a = autocor(y_sim, 0, 10);

a_lab = seqa(0,1,rows(a));

format /rd 12,2;

//DISPLAY AUTOCORRELATIONS

let a_col_lab = { "Lags" "ACF" };

print /rd $a_col_lab;;

print a_lab~a;

This code calculates and lists the sample autocorrelations for lags 0 through 10 for each sample of data. All columns of the list have a precision of two decimal points. Column labels “Lags” and “ACF” are included at the top of the list. This list will change depending on the data that has been randomly generated but should look similar to:

Lags          ACF

0.00         1.00         1.00         1.00         1.00

1.00         0.28         0.30         0.28         0.25

2.00        -0.71        -0.64        -0.63        -0.68

3.00        -0.59        -0.61        -0.54        -0.47

4.00         0.31         0.15         0.22         0.31

5.00         0.64         0.58         0.54         0.48

6.00         0.05         0.24         0.11        -0.05

7.00        -0.50        -0.31        -0.37        -0.38

8.00        -0.25        -0.39        -0.30        -0.05

9.00         0.33        -0.01         0.08         0.27

10.00        0.38         0.32         0.28         0.11

As an alternative, a bar graph visually presents the autocorrelation information. For example,

//DISPLAY CORRELOGRAM

//Initialize plotControl STRUCTURE
struct plotControl myPlot;

myPlot = plotGetDefaults("bar");

//ADD A TITLE
plotSetTitle(&myPlot, "Correlogram" );

//ADD AXIS LABLES
plotSetYLabel(&myPlot, "ACF");
plotSetXLabel(&myPlot, "Lag");

//PLOT ACF
x_lag=seqa(0, 1, 11);
plotBar(myPlot, x_lag, a[1:rows(a),1]);

This produces the following output:

The pacf function supplements the diagnostic information in the ACF and finds the partial autocorrelation function. The previous example is easily extended to find the PACF for the same randomly generated data. The pacf function requires three inputs: (1) an N x 1 data vector, x, (2) a scalar, k, denoting the maximum number of partial autocorrelation to compute, and (3) a scalar, d, denoting the order of differencing.  Note that in GAUSS the pacf function always begins with the minimum of one lag. Hence, if k=10, pacf will find the partial autocorrelations for lags 1 through 10. For example,

//FIND PACF
b = pacf(y[1:rows(y),1], 10, 0);

format /rd 12,2;

b_lab = seqa(0, 1, rows(a));

//DISPLAY AUTOCORRELATIONS

let b_col_lab = { "Lags" "ACF" };

print /rd $b_col_lab;;
print b_lab~a;

//PLOT PACF

//INITIALIZE plotControl STRUCTURE
struct plotControl myPlot;

myPlot = plotGetDefaults("bar");

//ADD A TITLE
plotSetTitle(&myPlot, "Partial Autocorrelation Function");

//ADD AXIS LABLES
plotSetYLabel(&myPlot, "PACF");
plotSetXLabel(&myPlot, "Lag");

//PLOT ACF
plotBar(myPlot, b_lab,b);

This code calculates, lists, and graphs the partial autocorrelations for lags 0 through 10.  This list and graph will change depending on the data that has been randomly generated but should look similar to:

Lags         PACF

1.00         0.28

2.00        -0.85

3.00         0.05

4.00         0.00

5.00        -0.06

6.00         0.00

7.00         0.04

8.00         0.09

9.00         0.09

10.00        0.07

Estimating the Autoregressive Model

The arimamt function is a convenient tool for estimating the parameters of any ARIMA model, including ARMA models, purely AR models, and purely MA models.  The arimamt function takes two arguments, (1) an arimamt control structure, (2) N x 1 data vector, y. Note that this function can only be used to estimate models with a specified number of sequential ordered lags. For example, it will estimate the parameters for the first through fourth lag but will not estimate a model that includes ONLY the first and fourth lag. In addition, note that this arimamt will produce an error if the data is not stationary and invertible.

The arimamt Control Structure

This structure handles the matrices and strings that control the estimation process such as setting specifying whether to include a constant, output appearance, and estimation parameters.  Continuing with the example above,

//CREATE AND SPECIFY ALL STRUCTURES
struct arimamtControl amc;
struct arimamtOut amo;

amc = arimamtControlCreate;

//SPECIFIES TO INCLUDE CONSTANT (Default)

amc.ar = 2;
amc.diff = 0;
amc.ma = 0;
amc.const = 1;

This appropriately passes to arimamt that the model is an AR(2) model and includes a constants.  The remaining controls elements, set at the default values, specify preferences for output display and estimation.

The Command File

Finally we put it all together in the command file:

//Aptech Systems, Inc
//ARMA MODEL ANALYSIS 
//Objective: Determining model parameters

library tsmt;
#include tsmt.sdf

//STEP ONE: RANDOMLY GENERATE DATA

//SPECIFY ARMA PARAMETERS
b_sim =  { 0.5, -0.8 };
q_sim = 0;
p_sim = 2;
const_sim = 1.5;

n_sim=200;
k_sim=4;
std_sim=1;
seed_sim = 5012;

y_sim = simarmamt(b_sim,p_sim,q_sim,const_sim,n_sim,k_sim,std_sim,seed_sim);

//VISUALIZE DATA
x_plot=seqa(1, 1, rows(y_sim));
//PLOT Y OBS
plotXY(x_plot, y_sim);

//CREATE AND SPECIFY ALL STRUCTURES
struct arimamtControl amc;
struct arimamtOut amo;

amc = arimamtControlCreate;

amc.ar = 2;
amc.diff = 0;
amc.ma = 0;
amc.const = 1;

//DEFINE OUTPUT FILE
output file = arimamt1.out reset;

//CALL ARIMA ROUTINE
amo = arimamt(amc, y_sim[1:rows(y_sim),1]);

//FORECAST USING RESULTS
f = tsforecastmt(amc, amo.b, z_sim,amo.e,10);

{alpha, tstat, zcrit} = vmadfmt(z_sim, -1, 3);

output off;

This program produces the following output:

Model:  ARIMA(2,0,0)

Final Results:

Iterations Until Convergence:  18

Log Likelihood:   -286.868570         Number of Residuals: 200

AIC           :    579.737139         Error Variance     : 1.032660354

SBC           :    589.632091         Standard Error     : 1.016198974

DF: 197     Adj. SSE: 206.257625544        SSE: 203.434089775

Coefficients     Std. Errors      T-Ratio    Approx. Prob.

AR1           0.67379713      1.27040000      0.53038          0.59644

AR2          -0.85412220      0.24522685     -3.48299          0.00061

Mean          1.49350894      0.07172501     20.82271          0.00000

Constant:   1.76282603

Total Computation Time: 0.03 (seconds)

AR Roots and Moduli:

Real :    0.39444   0.39444

Imag.:    1.00758  -1.00758

Mod. :    1.08203   1.08203

Forecasts for ARIMA(2,0,0) Model.   95% Confidence Interval Computed.

Period        LCL        Forecasts        UCL     Forecast Std. Err.

201     -0.710854      1.280860      3.272573      1.016199

202      0.684241      3.085889      5.487537      1.225353

203      0.217663      2.748079      5.278494      1.291052

204     -2.060362      0.978747      4.017856      1.550594

205     -2.997646      0.075108      3.147862      1.567761

206     -2.297213      0.977464      4.252141      1.670784

207     -1.113494      2.357287      5.828068      1.770839

208     -0.959802      2.516286      5.992373      1.773547

209     -2.204658      1.444881      5.094420      1.862044

210     -3.108973      0.587167      4.283308      1.885821