root/library/tests/arx_elem_test.cpp @ 679

#include "estim/arx.h"
using namespace bdm;

int main() {
        // Setup model : ARX for 1D Gaussian
        //Test constructor
        mat V0 = 0.00001 * eye ( 2 );
        V0 ( 0, 0 ) = 0.1; //
        ARX Ar;
        Ar.set_statistics ( 1, V0, -1.0 );
        Ar.set_constant(true);
        Ar.validate();

        mat mu ( 1, 1 );
        mat R ( 1, 1 );
        Ar.posterior().mean_mat ( mu, R );
        cout << "Prior moments: mu=" << mu << ", R=" << R << endl;

        int ndat = 200;
        vec smp = randn ( ndat );
        //
        mat Smp = ones ( 2, ndat );
        Smp.set_row ( 0, smp );
        //
        Ar.bayesB ( Smp );
        // Ar is now filled with estimates of N(0,1);
        cout << "Empirical moments: mu=" << sum ( smp ) / ndat << ", R=" << sum_sqr ( smp ) / ndat - pow ( sum ( smp ) / ndat, 2 ) << endl;
        Ar.posterior().mean_mat ( mu, R );
        cout << "Posterior moments: mu=" << mu << ", R=" << R << endl;

        //////// TEST prediction
        vec x = linspace ( -3.0, 3.0, 100 );
        double xstep = 6.0 / 100.0;
        mat X ( 1, 100 );
        mat X2 ( 2, 100 );
        X.set_row ( 0, x );
        X2.set_row ( 0, x );

        mlstudent* Ap = Ar.predictor_student();
        vec Ap_x = Ap->evallogcond_m ( X, empty_vec );
        vec ll_x = Ar.logpred_m ( X2 );

        cout << "normalize : " << xstep*sum ( exp ( Ap_x ) ) << endl;
        cout << "normalize : " << xstep*sum ( exp ( ll_x ) ) << endl;

        it_file it ( "arx_elem_test.it" );
        it << Name ( "Ap_x" ) << Ap_x;
        it << Name ( "ll_x" ) << ll_x;
}