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
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. However, the least squares method does not converge using the default relative error convergence criteria. As a result, a warning is produced and the estimates are very different from the method of moments and maximum likelihood estimates.
using System; using Imsl.Stat; using Imsl.Math; public class ARMAMaxLikelihoodEx1 { public static void Main(String[] args) { 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 = null; double[] arMM, maMM; double constantMM; double[] arLS, maLS; double constantLS; double[] ar, ma; double constant; double[,] forecastMM, forecastLS, forecast; double[] deviations; double likelihood, var, varMM, varLS; double[] avgDev = new double[]{0.0, 0.0, 0.0}; ARMA armaMM = null; ARMA armaLS = null; ARMAMaxLikelihood maxArma = null; double[,] printOutput = null; double[,] printOutput2 = null; 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(); pm.SetColumnSpacing(3); pmf.SetNoRowLabels(); 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.Method = Imsl.Stat.ARMA.ParamEstimation.MethodOfMoments; armaMM.Compute(); armaMM.BackwardOrigin = backwardOrigin; arMM = armaMM.GetAR(); maMM = armaMM.GetMA(); constantMM = armaMM.Constant; forecastMM = armaMM.Forecast(n_forecast); varMM = armaMM.InnovationVariance; /* Least Squares ARMA(2,1) Estimation */ armaLS = new ARMA(2, 1, z); armaLS.Method = Imsl.Stat.ARMA.ParamEstimation.LeastSquares; armaLS.Compute(); armaLS.BackwardOrigin = backwardOrigin; arLS = armaLS.GetAR(); maLS = armaLS.GetMA(); constantLS = armaLS.Constant; varLS = armaLS.InnovationVariance; forecastLS = armaLS.Forecast(n_forecast); /* Maximum Likelihood ARMA(2,1) Estimation */ maxArma = new ARMAMaxLikelihood(2, 1, z); maxArma.Compute(); maxArma.BackwardOrigin = backwardOrigin; ar = maxArma.GetAR(); ma = maxArma.GetMA(); constant = maxArma.Constant; likelihood = maxArma.Likelihood; var = maxArma.InnovationVariance; maxArma.Confidence = 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); Console.Out.WriteLine("INNOVATION VARIANCE:"); Console.Out.WriteLine("Method of Moments " + varMM); Console.Out.WriteLine("Least Squares " + varLS); Console.Out.WriteLine("Maximum Likelihood " + var); Console.Out.WriteLine(""); 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"); pmf.SetColumnLabels(colLabels2); pmf.FirstRowNumber = 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]; } pmf.SetColumnLabels(colLabels3); pmf.FirstRowNumber = 1860; pm.SetTitle("SUNSPOT MAX. LIKELIHOOD 90% CONFIDENCE INTERVALS"); pm.Print(pmf, printOutput2); } }
ARMA(2,1) - Constant Term Method of Moments Least Squares Maximum Likelihood 15.5439819435637 17.9316722242554 15.7580039183081 ARMA(2,1) - AR(1) Coefficient Method of Moments Least Squares Maximum Likelihood 1.24425777984372 1.53128221858841 1.22506595122835 ARMA(2,1) - AR(2) Coefficient Method of Moments Least Squares Maximum Likelihood -0.575149766040151 -0.894433318954122 -0.56051392292258 ARMA(2,1) - MA(1) Coefficient Method of Moments Least Squares Maximum Likelihood -0.124089747872598 -0.132018760615022 -0.382891621968495 INNOVATION VARIANCE: Method of Moments 287.242403738122 Least Squares 239.68797223859 Maximum Likelihood 214.508788425629 SUNSPOT FORECASTS FOR 1860-1869 Observed Sunspots Method of Moments Least Squares Maximum Likelihood 1860 95.7 87.6861734610048 98.0265069437306 87.9339016395945 1861 77.2 82.1446177467774 101.939296986813 82.0608538716595 1862 59.1 67.3203794962267 86.351331124405 66.9997857593059 1863 44 52.0624301952585 58.982026390731 51.8409090696475 1864 47 41.6037652345956 31.0142927589859 41.7122237493372 1865 30.5 37.3660959612162 12.6678176248844 37.8006776731389 1866 16.3 38.1086417046268 9.58945929412723 38.6860449014629 1867 7.3 41.469854513895 21.2853225650988 41.9631541830784 1868 37.3 45.2249946909406 41.9483762816881 45.4815685240649 1869 73.9 47.9641563097715 63.1281732161496 47.9549927562368 SUNSPOT MAX. LIKELIHOOD 90% CONFIDENCE INTERVALS Lower Confidence Limit Forecast Upper Confidence Limit 1860 63.8431804305709 87.9339016395945 112.024622848618 1861 36.4438944595113 82.0608538716595 127.677813283808 1862 10.134623371631 66.9997857593059 123.864948146981 1863 0 51.8409090696475 112.081240519376 1864 0 41.7122237493372 102.18743849516 1865 0 37.8006776731389 98.4521257842839 1866 0 38.6860449014629 99.9504700070846 1867 0 41.9631541830784 103.74778523107 1868 0 45.4815685240649 107.46464011018 1869 0 47.9549927562368 109.958371498445 Imsl.Stat.ARMA: Relative function convergence - Both the scaled actual and predicted reductions in the function are less than or equal to the relative function convergence tolerance "ConvergenceTolerance" = 5.00000041370186E-10. Imsl.Stat.ARMA: Least squares estimation of the parameters has failed to converge. Increase "length" and/or "tolerance" and/or "convergence_tolerance". The estimates of the parameters at the last iteration may be used as new starting values.Link to C# source.