root/bdm/estim/mixef.h @ 254

/*!
  \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 MEF_H
#define MEF_H

#include <itpp/itbase.h>
#include "../stat/libFN.h"
#include "../stat/libEF.h"
#include "../stat/emix.h"

namespace bdm{

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;
        //!Posterior on component parameters
        eprod* est;
        ////!Indeces of component rvc in common rvc
        
        //! Flag for a method that is used in the inference
        MixEF_METHOD method;
        
        //! Auxiliary function for use in constructors
        void build_est() {
                Array<const epdf*> epdfs ( n+1 );
                for ( int i=0;i<Coms.length();i++ ) {
//                      it_assert_debug(!x,"MixEF::MixEF : Incompatible components");
                        epdfs ( i ) =& ( Coms ( i )->_epdf() );
                }
                // last in the product is the weight
                epdfs ( n ) =& ( weights._epdf() );
                est = new eprod ( epdfs );
        }

public:
        //! Full constructor
        MixEF ( const Array<BMEF*> &Coms0, const vec &alpha0 ) :
                        BMEF ( RV() ), n ( Coms0.length() ), Coms ( n ),
                        weights ( RV ( "{w }", vec_1 ( n ) ),alpha0 ), method(QB) {
        //      it_assert_debug ( n>0,"MixEF::MixEF : Empty Component list" );

                for ( int i=0;i<n;i++ ) {Coms ( i ) = ( BMEF* ) Coms0 ( i )->_copy_();}
                build_est();
        };
        //! Constructor of empty mixture
        MixEF () :
                        BMEF ( RV() ), n ( 0 ), Coms ( 0 ),
                        weights ( RV ( "{w }", vec_1 ( 0 ) ),vec ( 0 ) ),method(QB) {build_est();}
        //! Copy constructor
        MixEF(const MixEF &M2): BMEF ( RV() ), n ( M2.n ), Coms ( n ),
                   weights ( M2.weights ), method(M2.method) {
        //      it_assert_debug ( n>0,"MixEF::MixEF : Empty Component list" );

                           for ( int i=0;i<n;i++ ) {Coms ( i ) = M2.Coms ( i )->_copy_();}
                           build_est();
                   }
        //! 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, int c=5 );
        //Destructor
        ~MixEF() {
                delete est;
                for ( int i=0;i<n;i++ ) {delete Coms ( i );}
        }
        //! Recursive EM-like algorithm (QB-variant), see Karny et. al, 2006
        void bayes ( const vec &dt );
        //! EM algorithm
        void bayes ( const mat &dt );
        void bayesB ( const mat &dt, const vec &wData );
        double logpred ( const vec &dt ) const;
        const epdf& _epdf() const {return *est;}
        const eprod* _e() const {return est;}
        emix* predictor(const RV &rv) 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;}
};

}
#endif // MEF_H