Example 2 (Logging and Forecasting): ARAutoUnivariate

Using the Canadian Lynx data included in TIMSAC-78, ARAutoUnivariate is used to find the minimum AIC autoregressive model using a maximum number of lags of maxlag =20.

This example compares the three different methods for estimating the autoregressive coefficients, and it illustrates the relationship between these estimates and those calculated within the TIMSAC UNIMAR code. As illustrated, the UNIMAR code estimates the coefficients and innovation variance using only the last N-maxlag values in the time series. The other estimation methods use all N-k values, where k is the number of lags with minimum AIC selected by ARAutoUnivariate .

This example also illustrates how to generate forecasts for the observed series values and beyond by setting the backward orgin for the forecasts.

import java.text.*;
import com.imsl.stat.*;
import com.imsl.math.*;
import com.imsl.math.PrintMatrix;
import java.util.logging.*;

public class ARAutoUnivariateEx2 {
    public static void main(String args[]) throws Exception {
        /* THE CANDIAN LYNX DATA AS USED IN TIMSAC 1821-1934 */
        double[] y =    
        {0.24300e01, 0.25060e01, 0.27670e01, 0.29400e01,
         0.31690e01, 0.34500e01, 0.35940e01, 0.37740e01,
         0.36950e01, 0.34110e01, 0.27180e01, 0.19910e01,
         0.22650e01, 0.24460e01, 0.26120e01, 0.33590e01,
         0.34290e01, 0.35330e01, 0.32610e01, 0.26120e01,
         0.21790e01, 0.16530e01, 0.18320e01, 0.23280e01,
         0.27370e01, 0.30140e01, 0.33280e01, 0.34040e01,
         0.29810e01, 0.25570e01, 0.25760e01, 0.23520e01,
         0.25560e01, 0.28640e01, 0.32140e01, 0.34350e01,
         0.34580e01, 0.33260e01, 0.28350e01, 0.24760e01,
         0.23730e01, 0.23890e01, 0.27420e01, 0.32100e01,
         0.35200e01, 0.38280e01, 0.36280e01, 0.28370e01,
         0.24060e01, 0.26750e01, 0.25540e01, 0.28940e01,
         0.32020e01, 0.32240e01, 0.33520e01, 0.31540e01,
         0.28780e01, 0.24760e01, 0.23030e01, 0.23600e01,
         0.26710e01, 0.28670e01, 0.33100e01, 0.34490e01,
         0.36460e01, 0.34000e01, 0.25900e01, 0.18630e01,
         0.15810e01, 0.16900e01, 0.17710e01, 0.22740e01,
         0.25760e01, 0.31110e01, 0.36050e01, 0.35430e01,
         0.27690e01, 0.20210e01, 0.21850e01, 0.25880e01,
         0.28800e01, 0.31150e01, 0.35400e01, 0.38450e01,
         0.38000e01, 0.35790e01, 0.32640e01, 0.25380e01,
         0.25820e01, 0.29070e01, 0.31420e01, 0.34330e01,
         0.35800e01, 0.34900e01, 0.34750e01, 0.35790e01,
         0.28290e01, 0.19090e01, 0.19030e01, 0.20330e01,
         0.23600e01, 0.26010e01, 0.30540e01, 0.33860e01,
         0.35530e01, 0.34680e01, 0.31870e01, 0.27230e01,
         0.26860e01, 0.28210e01, 0.30000e01, 0.32010e01,
         0.34240e01, 0.35310e01};
        double[][] printOutput=null;
        double timsacAR[], mmAR[], mleAR[], lsAR[];
        double forecasts[], residuals[];
        double timsacConstant, mmConstant, mleConstant, lsConstant;
        double timsacVar, timsacEquivalentVar, mmVar, mleVar, lsVar;
        int    maxlag = 20;
        String[] colLabels = {"TIMSAC", "Method of Moments", "Least Squares", 
        "Maximum Likelihood"};
        String[] colLabels2= {"Observed", "Forecast", "Residual"};
        PrintMatrixFormat pmf = new PrintMatrixFormat();
        PrintMatrix pm        = new PrintMatrix();
        NumberFormat nf = NumberFormat.getNumberInstance();
        pm.setColumnSpacing(4);
        nf.setMinimumFractionDigits(4);
        pmf.setNumberFormat(nf);
        pmf.setColumnLabels(colLabels);
        System.out.println("Automatic Selection of Minimum AIC AR Model");
        System.out.println("");
        
        ARAutoUnivariate autoAR = new ARAutoUnivariate(maxlag, y);

        // logging to console
        Logger logger = autoAR.getLogger();
        ConsoleHandler ch = new ConsoleHandler();
        ch.setLevel(Level.ALL);         // default ConsoleHandler Level is INFO
        logger.setLevel(Level.FINE);
        logger.addHandler(ch);
        ch.setFormatter(new com.imsl.IMSLFormatter());

        autoAR.compute();
        int orderSelected = autoAR.getOrder();
        System.out.println("Minimum AIC Selected="+autoAR.getAIC()+
                           " with an optimum lag of k= "+ autoAR.getOrder());
        System.out.println("");
        
        timsacAR       = autoAR.getTimsacAR();
        timsacConstant = autoAR.getTimsacConstant();
        timsacVar      = autoAR.getTimsacVariance();
        lsAR           = autoAR.getAR();
        lsConstant     = autoAR.getConstant();
        lsVar          = autoAR.getInnovationVariance();
        
        autoAR.setEstimationMethod(ARAutoUnivariate.METHOD_OF_MOMENTS);
        autoAR.compute();
        mmAR       = autoAR.getAR();
        mmConstant = autoAR.getConstant();
        mmVar      = autoAR.getInnovationVariance();
        
        autoAR.setEstimationMethod(ARAutoUnivariate.MAX_LIKELIHOOD);
        autoAR.compute();
        mleAR       = autoAR.getAR();
        mleConstant = autoAR.getConstant();
        mleVar      = autoAR.getInnovationVariance();
        
        printOutput = new double[orderSelected+1][4];
        printOutput[0][0] = timsacConstant;
        for(int i=0; i<orderSelected; i++) printOutput[i+1][0] = timsacAR[i];
        printOutput[0][1] = mmConstant;
        for(int i=0; i<orderSelected; i++) printOutput[i+1][1] = mmAR[i];
        printOutput[0][2] = lsConstant;
        for(int i=0; i<orderSelected; i++) printOutput[i+1][2] = lsAR[i];
        printOutput[0][3] = mleConstant;
        for(int i=0; i<orderSelected; i++) printOutput[i+1][3] = mleAR[i];
        pm.setTitle("Comparison of AR Estimates");
        pm.print(pmf, printOutput);
       
        /* calculation of equivalent innovation variance using TIMSAC 
          coefficients.  The Timsac innovation variance is calculated using
          only N-maxlag observations in the series.  The following code 
          calculates the innovation variance using N-k observations in the 
          series with the Timsac coefficient.  This illustrates that the
          least squares Timsac coefficients will not have the least value for
          the sum of squared residuals, which is calculated using all N-k
          observations. */
        ARMA armaLS = new ARMA(orderSelected,0,y);
        armaLS.setArmaInfo(autoAR.getTimsacConstant(), autoAR.getTimsacAR(), 
                new double[0], autoAR.getTimsacVariance());
        armaLS.setBackwardOrigin(y.length-orderSelected);
        forecasts = armaLS.getForecast(1);
        double sumResiduals = 0.0;
        for(int i=0; i<y.length-orderSelected; i++){
            sumResiduals += (y[i+orderSelected]-forecasts[i])*
                    (y[i+orderSelected]-forecasts[i]);
        }
        timsacEquivalentVar = sumResiduals / (y.length-orderSelected-1);
        printOutput = new double[1][4];
        printOutput[0][0] = timsacEquivalentVar; 
               /* the method of moments variance */
        printOutput[0][1] = mmVar;  
               /* the least squares variance */
        printOutput[0][2] = lsVar;  
               /* the maximum likelihood estimate of the variance */
        printOutput[0][3] = mleVar; 
        nf.setMinimumFractionDigits(5);
        pmf.setNumberFormat(nf);
        pm.setTitle("Comparison of Equivalent Innovation Variances");
        pm.print(pmf, printOutput);
        
        /* FORECASTING - An example of forecasting using the maximum
         * likelihood estimates for the minimum AIC AR model.  In this example,
         * forecasts are returned for the last 10 values in the series followed
         * by the forecasts for the next 5 values.
         */
        autoAR.setBackwardOrigin(10);
        forecasts = autoAR.getForecast(15);
        residuals = autoAR.getResiduals();
        printOutput = new double[15][3];
        for(int i=0; i<10; i++){
            printOutput[i][0] = y[y.length-10+i];
            printOutput[i][1] = forecasts[i];
            printOutput[i][2] = residuals[i];
        }
        for(int i=10; i<15; i++){
            printOutput[i][0] = Double.NaN;
            printOutput[i][1] = forecasts[i];
            printOutput[i][2] = Double.NaN;
        }
        nf.setMaximumFractionDigits(3);
        pmf.setFirstRowNumber(105);
        pmf.setNumberFormat(nf);
        pmf.setColumnLabels(colLabels2);
        pm.setTitle("Maximum Likelihood Forecasts of Last 10 Observations & " +
                "the Next 5");
        pm.print(pmf, printOutput);
    }
}

Output

Automatic Selection of Minimum AIC AR Model

ORDER SELECTED = 11
AIC = -122.36968839314528


Minimum AIC Selected=-296.13013263562374 with an optimum lag of k= 11

ORDER SELECTED = 11
AIC = -122.36968839314528


ORDER SELECTED = 11
AIC = -122.36968839314528


                         Comparison of AR Estimates
       TIMSAC     Method of Moments    Least Squares    Maximum Likelihood   
 0     1.0427          1.1679             1.1144             1.1187          
 1     1.1813          1.1381             1.1481             1.1664          
 2    -0.5516         -0.5061            -0.5331            -0.5419          
 3     0.2314          0.2098             0.2757             0.2624          
 4    -0.1780         -0.2672            -0.3263            -0.3052          
 5     0.0199          0.1112             0.1685             0.1519          
 6    -0.0626         -0.1246            -0.1643            -0.1460          
 7     0.0286          0.0693             0.0728             0.0581          
 8    -0.0507         -0.0419            -0.0305            -0.0310          
 9     0.1999          0.1366             0.1509             0.1379          
10     0.1618          0.1828             0.1935             0.1995          
11    -0.3391         -0.3101            -0.3414            -0.3375          

               Comparison of Equivalent Innovation Variances
      TIMSAC     Method of Moments    Least Squares    Maximum Likelihood   
0    0.03769         0.04274            0.03687            0.03618          

Maximum Likelihood Forecasts of Last 10 Observations & the Next 5
        Observed    Forecast    Residual   
105     3.553       3.439        0.114     
106     3.468       3.480       -0.012     
107     3.187       2.924        0.263     
108     2.723       2.703        0.020     
109     2.686       2.556        0.130     
110     2.821       2.785        0.036     
111     3.000       2.949        0.051     
112     3.201       3.186        0.015     
113     3.424       3.386        0.038     
114     3.531       3.527        0.004     
115     ?           3.446        ?         
116     ?           3.195        ?         
117     ?           2.829        ?         
118     ?           2.492        ?         
119     ?           2.414        ?         

Link to Java source.