/*!
  \file
  \brief Bayesian Filtering for mixtures of exponential family (EF) members
  \author Vaclav Smidl.

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

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

#ifndef MIXTURES_H
#define MIXTURES_H


#include "../math/functions.h"
#include "../stat/exp_family.h"
#include "../stat/emix.h"
#include "arx.h"

namespace bdm {

//! enum switch for internal approximation used in MixEF
enum MixEF_METHOD { EM = 0, QB = 1};

/*!
* \brief Mixture of Exponential Family Densities

An approximate estimation method for models with latent discrete variable, such as
mixture models of the following kind:
\f[
f(y_t|\psi_t, \Theta) = \sum_{i=1}^{n} w_i f(y_t|\psi_t, \theta_i)
\f]
where  \f$\psi\f$ is a known function of past outputs, \f$w=[w_1,\ldots,w_n]\f$ are component weights, and component parameters \f$\theta_i\f$ are assumed to be mutually independent. \f$\Theta\f$ is an aggregation af all component parameters and weights, i.e. \f$\Theta = [\theta_1,\ldots,\theta_n,w]\f$.

The characteristic feature of this model is that if the exact values of the latent variable were known, estimation of the parameters can be handled by a single model. For example, for the case of mixture models, posterior density for each component parameters would be a BayesianModel from Exponential Family.

This class uses EM-style type algorithms for estimation of its parameters. Under this simplification, the posterior density is a product of exponential family members, hence under EM-style approximate estimation this class itself belongs to the exponential family.

TODO: Extend BM to use rvc.
*/
class MixEF: public BMEF {
protected:
	//!Number of components
	int n;
	//! Models for Components of \f$\theta_i\f$
	Array<BMEF*> Coms;
	//! Statistics for weights
	multiBM weights;
	//aux
	friend class eprod_mix;
	//!Posterior on component parameters
	class eprod_mix: public eprod_base {
		protected:
		const MixEF &mix; // pointer to parent.n
		public:
		eprod_mix(const MixEF &m):mix(m){}
		const epdf* factor(int i) const {return (i==(mix.n-1)) ? &mix.weights.posterior() : &mix.Coms(i)->posterior();}
		const int no_factors()const {return mix.n+1;}
	} est;
	////!indices of component rvc in common rvc

	//! Flag for a method that is used in the inference
	MixEF_METHOD method;

public:
	//! Full constructor
	MixEF ( const Array<BMEF*> &Coms0, const vec &alpha0 ) :
			BMEF ( ), n ( Coms0.length() ), Coms ( n ),
			weights (), est(*this), method ( QB ) {
		for ( int i = 0; i < n; i++ ) {
			Coms ( i ) = ( BMEF* ) Coms0 ( i )->_copy();
		}
		weights.set_parameters(alpha0);
		weights.validate();
	}

	//! Constructor of empty mixture
	MixEF () :
			BMEF ( ), n ( 0 ), Coms ( 0 ),
			weights (), est(*this), method ( QB ) {
	}
	//! Copy constructor
	MixEF ( const MixEF &M2 ) : BMEF ( ), n ( M2.n ), Coms ( n ),
			weights ( M2.weights ), est(*this), method ( M2.method ) {
		for ( int i = 0; i < n; i++ ) {
			Coms ( i ) = (BMEF*) M2.Coms ( i )->_copy();
		}
	}

	//! Initializing the mixture by a random pick of centroids from data
	//! \param Com0 Initial component - necessary to determine its type.
	//! \param Data Data on which the initialization will be done
	//! \param c Initial number of components, default=5
	void init ( BMEF* Com0, const mat &Data, const int c = 5 );
	//Destructor
	//! Recursive EM-like algorithm (QB-variant), see Karny et. al, 2006
	void bayes ( const vec &yt, const vec &cond );
	//! EM algorithm
	void bayes ( const mat &yt, const vec &cond );
	//! batch weighted Bayes rule
	void bayes_batch ( const mat &yt, const mat &cond, const vec &wData );
	double logpred ( const vec &yt ) const;
	//! return correctly typed posterior (covariant return)
	const eprod_mix& posterior() const {
		return est;
	}

	emix* epredictor(const vec &cond=vec()) const;
	//! Flatten the density as if it was not estimated from the data
	void flatten ( const BMEF* M2 );
	//! Access function
	BMEF* _Coms ( int i ) {
		return Coms ( i );
	}

	//!Set which method is to be used
	void set_method ( MixEF_METHOD M ) {
		method = M;
	}

	void to_setting ( Setting &set ) const {
		BMEF::to_setting( set );
		UI::save ( Coms, set, "Coms" );
		UI::save ( &weights, set, "weights" );
	}
};
UIREGISTER ( MixEF );

class ARXprod: public ProdBMBase {
	Array<shared_ptr<ARX> > arxs;
	public:
		ARX* bm(int i){return arxs(i).get();}
		int no_bms() {return arxs.length();}
};
UIREGISTER(ARXprod);

}
#endif // MIXTURES_H
