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 = null;
        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};
        double sum=0;
        ARMA armaMM = null;
        ARMA armaLS = null;
        ARMAMaxLikelihood maxArma = null;
        double[][] printOutput    = null;
        double[][] printOutput2   = null;
        String title="ARMA(p=2, q=1)";
        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);
       
    }
    
}