package com.imsl.test.example.stat;

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

/**
 * <p>
 * Computes the filtered estimates and the one-step-ahead estimates using the
 * Kalman filter.
 * </p>
 * <p>
 * In this example, {@link KalmanFilter} is used to compute the filtered 
 * estimates and one-step-ahead
 * estimates for a scalar problem discussed by Harvey (1981, pages 116-117). The
 * observation equation and state equation are given by</p>
 * <p>
 * $$y_k = b_k + e_k$$</p>
 * <p>
 * $$b_{k+1} = b_k + w_{k+1}$$</p>
 * <p>
 * for \(k = 1,2,3,4\). The \(e_k\)s are identically and independently distributed normal with
 * mean 0 and variance \(\sigma ^2\), the \(w_k\)s are identically and
 * independently distributed normal with mean 0 and variance \(4 \sigma ^2\),
 * and \(b_1\) is distributed normal with mean 4 and variance \(16 \sigma ^2\).
 * Two <code>KalmanFilter</code> objects are needed for each time point in order to compute
 * the filtered estimate and the one-step-ahead estimate. The first object does
 * not use the methods <code>SetTransitionMatrix</code> and <code>setQ</code> so
 * that the prediction equations are skipped in the computations. The update
 * equations are skipped in the computations in the second object.</p>
 * <p>
 * This example also computes the one-step-ahead prediction errors. Note that
 * Harvey (1981, page 117) contains a misprint for the value \(v_4\) that he
 * gives as \(1.197\). The correct value of \(v_4 = 1.003\) is computed by
 * <code>KalmanFilter</code>.</p>
 *
 * @see <a href="KalmanFilterEx1.java">Code</a>
 * @see <a href="doc-files/KalmanFilterEx1.out">Output</a>
 */
public class KalmanFilterEx1 {

    static private final MessageFormat mf
            = new MessageFormat("{0}/{1}\t{2}\t{3}\t{4}\t{5}\t{6}\t{7}\t{8}");

    public static void main(String args[]) {
        int nobs = 4;
        int rank = 0;
        double logDeterminant = 0.0;
        double ss = 0.0;
        double[] b = {4};
        double[][] covb = {{16}};
        double[][] q = {{4}};
        double[][] r = {{1}};
        double[][] t = {{1}};
        double[][] z = {{1}};
        double[] ydata = {4.4, 4.0, 3.5, 4.6};

        Object argFormat[] = {
            "k", "j", "b", "cov(b)", "rank", "ss", "ln(det)", "v", "cov(v)"
        };
        System.out.println(mf.format(argFormat));

        KalmanFilter kalman
                = new KalmanFilter(b, covb, rank, ss, logDeterminant);

        for (int i = 0; i < nobs; i++) {
            double y[] = {ydata[i]};
            kalman.update(y, z, r);
            kalman.filter();
            double v[] = kalman.getPredictionError();
            double covv[][] = kalman.getCovV();
            argFormat[0] = new Integer(i);
            argFormat[1] = new Integer(i);
            argFormat[2] = new Double(kalman.getStateVector()[0]);
            argFormat[3] = new Double(kalman.getCovB()[0][0]);
            argFormat[4] = new Integer(kalman.getRank());
            argFormat[5] = new Double(kalman.getSumOfSquares());
            argFormat[6] = new Double(kalman.getLogDeterminant());
            argFormat[7] = new Double(v[0]);
            argFormat[8] = new Double(covv[0][0]);
            System.out.println(mf.format(argFormat));
            kalman.resetUpdate();

            kalman.setTransitionMatrix(t);
            kalman.setQ(q);
            kalman.filter();
            argFormat[0] = new Integer(i + 1);
            argFormat[1] = new Integer(i);
            argFormat[2] = new Double(kalman.getStateVector()[0]);
            argFormat[3] = new Double(kalman.getCovB()[0][0]);
            argFormat[4] = new Integer(kalman.getRank());
            argFormat[5] = new Double(kalman.getSumOfSquares());
            argFormat[6] = new Double(kalman.getLogDeterminant());
            argFormat[7] = new Double(v[0]);
            argFormat[8] = new Double(covv[0][0]);
            System.out.println(mf.format(argFormat));
            kalman.resetTransitionMatrix();
            kalman.resetQ();
        }
    }
}