package com.imsl.test.example.stat; import java.text.*; import com.imsl.stat.*; import com.imsl.math.*; /** *
* Fits an \(\text{ARMA}(2,1)\) to the Wolfer * sunspot data using the method of maximum likelihood.
** 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 \(ARMA(2, 1)\) model can be expressed as
** $$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 *
* $${\rm{\hat \theta }}_{\rm{0}} {\rm{,\hat \phi }}_{\rm{1}} {\rm{,\hat \phi * }}_{\rm{2}} {\rm{, and \, \hat \theta }}_{\rm{1}}$$
* fromARMAMaxLikelihood
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.
*
* @see Code
* @see Output
*/
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);
}
}