Example: ARMAMaxLikelihood

Consider the Wolfer Sunspot Data (Anderson 1971, p. 660) consisting of the number of sunspots observed each year from 1749 through 1924. The data set for this example consists of the number of sunspots observed from 1770 through 1869. The maximum likelihood estimates

{\rm{\hat \theta }}_{\rm{0}} {\rm{,\hat \phi }}_{\rm{1}} {\rm{,\hat \phi }}_{\rm{2}} {\rm{, and \, \hat \theta }}_{\rm{1}}

for the ARMA(2, 1) model

z_t = \theta _0 + \phi _1 z_{t - 1} + \phi _2 z_{t - 2} - \theta _1 A_{t - 1} + A_t

where the errors A_t are independently normally distributed with mean zero and variance

\sigma _A^2

The maximum likelihood estimates from ARMAMaxLikelihood are compared to the same estimates using the method of moments and least squares from the ARMA class. For each method, the coefficients and forecasts for the last ten years, 1860-1869, are compared. The method of moments and maximum likelihood estimates produced similar results, but the least squares estimates were very different from the other two.


import java.text.*;
import com.imsl.stat.*;
import com.imsl.math.*;

public class ARMAMaxLikelihoodEx1 {

    public static void main(String args[]) throws Exception {
        int backwardOrigin = 0, n_forecast = 10, n_series = 100;
        double[] sunspots = {
            100.8, 81.6, 66.5, 34.8, 30.6, 7, 19.8, 92.5,
            154.4, 125.9, 84.8, 68.1, 38.5, 22.8, 10.2, 24.1, 82.9,
            132, 130.9, 118.1, 89.9, 66.6, 60, 46.9, 41, 21.3, 16,
            6.4, 4.1, 6.8, 14.5, 34, 45, 43.1, 47.5, 42.2, 28.1, 10.1,
            8.1, 2.5, 0, 1.4, 5, 12.2, 13.9, 35.4, 45.8, 41.1, 30.4,
            23.9, 15.7, 6.6, 4, 1.8, 8.5, 16.6, 36.3, 49.7, 62.5,
            67, 71, 47.8, 27.5, 8.5, 13.2, 56.9, 121.5, 138.3, 103.2,
            85.8, 63.2, 36.8, 24.2, 10.7, 15, 40.1, 61.5, 98.5,
            124.3, 95.9, 66.5, 64.5, 54.2, 39, 20.6, 6.7, 4.3, 22.8,
            54.8, 93.8, 95.7, 77.2, 59.1, 44, 47, 30.5, 16.3, 7.3,
            37.3, 73.9
        };
        double[] z;
        double arMM[], maMM[], constantMM;
        double arLS[], maLS[], constantLS;
        double ar[], ma[], constant;
        double forecastMM[][], forecastLS[][], forecast[][];
        double deviations[];
        double likelihood, var, varMM, varLS;
        double[] avgDev = {0.0, 0.0, 0.0};
        ARMA armaMM, armaLS;
        ARMAMaxLikelihood maxArma;
        double[][] printOutput, printOutput2;
        String[] colLabels = {
            "Method of Moments", "Least Squares", "Maximum Likelihood"
        };
        String[] colLabels1 = {"Least Squares", "Maximum Likelihood"};
        String[] colLabels2 = {
            "Observed Sunspots", "Method of Moments",
            "Least Squares", "Maximum Likelihood"
        };
        String[] colLabels3 = {
            "Lower Confidence Limit", "Forecast", "Upper Confidence Limit"
        };
        PrintMatrixFormat pmf = new PrintMatrixFormat();
        PrintMatrix pm = new PrintMatrix();
        NumberFormat nf = NumberFormat.getNumberInstance();
        pm.setColumnSpacing(3);
        nf.setMinimumFractionDigits(3);
        pmf.setNoRowLabels();
        pmf.setNumberFormat(nf);
        pmf.setColumnLabels(colLabels);
        printOutput = new double[1][3];
        printOutput2 = new double[n_forecast][4];
        z = new double[n_series];
        for (int i = 0; i < n_series; i++) {
            z[i] = sunspots[i];
        }
        /* Method of Moments ARMA(2,1) Estimation */
        armaMM = new ARMA(2, 1, z);
        armaMM.setMethod(ARMA.METHOD_OF_MOMENTS);
        armaMM.compute();
        armaMM.setBackwardOrigin(backwardOrigin);
        arMM = armaMM.getAR();
        maMM = armaMM.getMA();
        constantMM = armaMM.getConstant();
        forecastMM = armaMM.forecast(n_forecast);
        varMM = armaMM.getInnovationVariance();

        /* Least Squares ARMA(2,1) Estimation */
        armaLS = new ARMA(2, 1, z);
        armaLS.setMethod(ARMA.LEAST_SQUARES);
        armaLS.compute();
        armaLS.setBackwardOrigin(backwardOrigin);
        arLS = armaLS.getAR();
        maLS = armaLS.getMA();
        constantLS = armaLS.getConstant();
        varLS = armaLS.getInnovationVariance();
        forecastLS = armaLS.forecast(n_forecast);

        /* Maximum Likelihood ARMA(2,1) Estimation */
        maxArma = new ARMAMaxLikelihood(2, 1, z);
        maxArma.compute();
        maxArma.setBackwardOrigin(backwardOrigin);
        ar = maxArma.getAR();
        ma = maxArma.getMA();
        constant = maxArma.getConstant();
        likelihood = maxArma.getLikelihood();
        var = maxArma.getInnovationVariance();
        maxArma.setConfidence(0.9);
        forecast = maxArma.forecast(n_forecast);
        deviations = maxArma.getDeviations();

        printOutput[0][0] = constantMM;
        printOutput[0][1] = constantLS;
        printOutput[0][2] = constant;
        pm.setTitle("ARMA(2,1) - Constant Term");
        pm.print(pmf, printOutput);
        printOutput[0][0] = arMM[0];
        printOutput[0][1] = arLS[0];
        printOutput[0][2] = ar[0];
        pm.setTitle("ARMA(2,1) - AR(1) Coefficient");
        pm.print(pmf, printOutput);
        printOutput[0][0] = arMM[1];
        printOutput[0][1] = arLS[1];
        printOutput[0][2] = ar[1];
        pm.setTitle("ARMA(2,1) - AR(2) Coefficient");
        pm.print(pmf, printOutput);
        printOutput[0][0] = maMM[0];
        printOutput[0][1] = maLS[0];
        printOutput[0][2] = ma[0];
        pm.setTitle("ARMA(2,1) - MA(1) Coefficient");
        pm.print(pmf, printOutput);
        System.out.println("INNOVATION VARIANCE:");
        System.out.println("Method of Moments  " + varMM);
        System.out.println("Least Squares      " + varLS);
        System.out.println("Maximum Likelihood " + var);
        System.out.println("");

        for (int i = 0; i < n_forecast; i++) {
            printOutput2[i][0] = sunspots[90 + i];
            printOutput2[i][1] = forecastMM[i][backwardOrigin];
            printOutput2[i][2] = forecastLS[i][backwardOrigin];
            printOutput2[i][3] = forecast[i][backwardOrigin];
        }
        pm.setTitle("SUNSPOT FORECASTS FOR 1860-1869");
        nf.setMaximumFractionDigits(0);
        pmf.setNumberFormat(nf);
        pmf.setColumnLabels(colLabels2);
        pmf.setFirstRowNumber(1860);
        pm.print(pmf, printOutput2);
        /* Get Confidence Interval Deviations */
        printOutput2 = new double[n_forecast][3];
        for (int i = 0; i < n_forecast; i++) {
            printOutput2[i][0] = Math.max(0, forecast[i][backwardOrigin]
                    - deviations[i + backwardOrigin]);
            printOutput2[i][1] = forecast[i][backwardOrigin];
            printOutput2[i][2] = forecast[i][backwardOrigin] + deviations[i];
        }
        nf.setMaximumFractionDigits(0);
        pmf.setNumberFormat(nf);
        pmf.setColumnLabels(colLabels3);
        pmf.setFirstRowNumber(1860);
        pm.setTitle("SUNSPOT MAX. LIKELIHOOD 90% CONFIDENCE INTERVALS");
        pm.print(pmf, printOutput2);
    }
}

Output

                 ARMA(2,1) - Constant Term
   Method of Moments   Least Squares   Maximum Likelihood   
     15.544            17.932             15.758         

               ARMA(2,1) - AR(1) Coefficient
   Method of Moments   Least Squares   Maximum Likelihood   
      1.244             1.531             1.225          

               ARMA(2,1) - AR(2) Coefficient
   Method of Moments   Least Squares   Maximum Likelihood   
     -0.575            -0.894             -0.561         

               ARMA(2,1) - MA(1) Coefficient
   Method of Moments   Least Squares   Maximum Likelihood   
     -0.124            -0.132             -0.383         

INNOVATION VARIANCE:
Method of Moments  287.2424037381216
Least Squares      239.68797223858988
Maximum Likelihood 214.5087884256287

                          SUNSPOT FORECASTS FOR 1860-1869
       Observed Sunspots   Method of Moments   Least Squares   Maximum Likelihood   
1860          96                  88                 98                88           
1861          77                  82                102                82           
1862          59                  67                 86                67           
1863          44                  52                 59                52           
1864          47                  42                 31                42           
1865          30                  37                 13                38           
1866          16                  38                 10                39           
1867           7                  41                 21                42           
1868          37                  45                 42                45           
1869          74                  48                 63                48           

          SUNSPOT MAX. LIKELIHOOD 90% CONFIDENCE INTERVALS
       Lower Confidence Limit   Forecast   Upper Confidence Limit   
1860             64                88               112             
1861             36                82               128             
1862             10                67               124             
1863              0                52               112             
1864              0                42               102             
1865              0                38                98             
1866              0                39               100             
1867              0                42               104             
1868              0                45               107             
1869              0                48               110             

Link to Java source.