/*!
  \file
  \brief Bayesian Filtering for generalized autoregressive (ARX) model
  \author Vaclav Smidl.

  -----------------------------------
  BDM++ - C++ library for Bayesian Decision Making under Uncertainty

  Using IT++ for numerical operations
  -----------------------------------
*/

#ifndef AR_H
#define AR_H

#include "../math/functions.h"
#include "../stat/exp_family.h"
#include "../base/user_info.h"

namespace bdm {

/*!
* \brief Linear Autoregressive model with Gaussian noise

Regression of the following kind:
\f[
y_t = \theta_1 \psi_1 + \theta_2 + \psi_2 +\ldots + \theta_n \psi_n + r e_t
\f]
where unknown parameters \c rv are \f$[\theta r]\f$, regression vector \f$\psi=\psi(y_{1:t},u_{1:t})\f$ is a known function of past outputs and exogeneous variables \f$u_t\f$. Distrubances \f$e_t\f$ are supposed to be normally distributed:
\f[
e_t \sim \mathcal{N}(0,1).
\f]

See \ref tut_arx for mathematical treatment.

The easiest way how to use the class is:
\include arx_simple.cpp

	\todo sort out constant terms - bayes should accept vec without additional 1s
*/
class ARX: public BMEF {
protected:
	//!size of output variable (needed in regressors)
	int dimx;
	//!description of modelled data \f$ y_t \f$ in the likelihood function
	//! Do NOT access directly, only via \c get_yrv().
	RV _yrv;
	//! Posterior estimate of \f$\theta,r\f$ in the form of Normal-inverse Wishart density
	egiw est;
	//! cached value of est.V
	ldmat &V;
	//! cached value of est.nu
	double &nu;
public:
	//! \name Constructors
	//!@{
	ARX ( const double frg0 = 1.0 ) : BMEF ( frg0 ), est (), V ( est._V() ), nu ( est._nu() ) {};
	ARX ( const ARX &A0 ) : BMEF (), est (), V ( est._V() ), nu ( est._nu() ) {
		set_statistics ( A0.dimx, A0.V, A0.nu );
		set_parameters ( A0.frg );
	};
	ARX* _copy_() const;
	void set_parameters ( double frg0 ) {
		frg = frg0;
	}
	void set_statistics ( int dimx0, const ldmat V0, double nu0 = -1.0 ) {
		est.set_parameters ( dimx0, V0, nu0 );
		last_lognc = est.lognc();
		dimx = dimx0;
	}
	//!@}

//	//! Set parameters given by moments, \c mu (mean of theta), \c R (mean of R) and \c C (variance of theta)
//	void set_parameters ( const vec &mu, const mat &R, const mat &C, double dfm){};
	//! Set sufficient statistics
	void set_statistics ( const BMEF* BM0 );
//	//! Returns sufficient statistics
//	void get_parameters ( mat &V0, double &nu0 ) {V0=est._V().to_mat(); nu0=est._nu();}
	//!\name Mathematical operations
	//!@{

	//! Weighted Bayes \f$ dt = [y_t psi_t] \f$.
	void bayes ( const vec &dt, const double w );
	void bayes ( const vec &dt ) {
		bayes ( dt, 1.0 );
	};
	double logpred ( const vec &dt ) const;
	void flatten ( const BMEF* B ) {
		const ARX* A = dynamic_cast<const ARX*> ( B );
		// nu should be equal to B.nu
		est.pow ( A->nu / nu );
		if ( evalll ) {
			last_lognc = est.lognc();
		}
	}
	//! Conditioned version of the predictor
	enorm<ldmat>* epredictor ( const vec &rgr ) const;
	//! Predictor for empty regressor
	enorm<ldmat>* epredictor() const {
		bdm_assert_debug ( dimx == V.rows() - 1, "Regressor is not only 1" );
		return epredictor ( vec_1 ( 1.0 ) );
	}
	//! conditional version of the predictor
	mlnorm<ldmat>* predictor() const;
	mlstudent* predictor_student() const;
	//! Brute force structure estimation.\return indeces of accepted regressors.
	ivec structure_est ( egiw Eg0 );
	//! Smarter structure estimation by Ludvik Tesar.\return indeces of accepted regressors.
	ivec structure_est_LT ( egiw Eg0 );
	//!@}

	//!\name Access attributes
	//!@{
	const egiw& posterior() const {
		return est;
	}
	//!@}

	//!\name Connection
	//!@{
	void set_drv ( const RV &drv0 ) {
		drv = drv0;
	}

	RV& get_yrv() {
		//if yrv is not ready create it
		if ( _yrv._dsize() != dimx ) {
			int i = 0;
			while ( _yrv._dsize() < dimx ) {
				_yrv.add ( drv ( vec_1 ( i ) ) );
				i++;
			}
		}
		//yrv should be ready by now
		bdm_assert_debug ( _yrv._dsize() == dimx, "incompatible drv" );
		return _yrv;
	}
	//!@}

	// TODO dokumentace - aktualizovat
	/*! UI for ARX estimator

	The ARX is constructed from a structure with fields:
	\code
	estimator = {
		type = "ARX";
		y = {type="rv", ...}   // description of output variables
		rgr = {type="rv", ...} // description of regressor variables
		constant = true;       // boolean switch if the constant term is modelled or not

		//optional fields
		dV0 = [1e-3, 1e-5, 1e-5, 1e-5];
							   // default: 1e-3 for y, 1e-5 for rgr
		nu0 = 6;               // default: rgrlen + 2
		frg = 1.0;             // forgetting, default frg=1.0
	};
	\endcode

	The estimator will assign names of the posterior in the form ["theta_i" and "r_i"]
	*/
	void from_setting ( const Setting &set );

};

UIREGISTER ( ARX );
SHAREDPTR ( ARX );

}

#endif // AR_H


