- Main Page
- Related Pages
- Classes
- Files

# work/git/mixpp/bdm/stat/emix.h

Go to the documentation of this file.

    00001
    00013 #ifndef MX_H
    00014 #define MX_H
    00015
    00016 #include "libBM.h"
    00017 #include "libEF.h"
    00018 //#include <std>
    00019
    00020 using namespace itpp;
    00021
    00032 class emix : public epdf {
    00033 protected:
    00035         vec w;
    00037         Array<epdf*> Coms;
    00038 public:
    00040         emix ( RV &rv ) : epdf ( rv ) {};
    00042         void set_parameters ( const vec &w, const Array<epdf*> &Coms );
    00043
    00044         vec sample() const;
    00045         vec mean() const {
    00046                 int i; vec mu = zeros ( rv.count() );
    00047                 for ( i = 0;i < w.length();i++ ) {mu += w ( i ) * Coms ( i )->mean(); }
    00048                 return mu;
    00049         }
    00050         double evalpdflog ( const vec &val ) const {
    00051                 int i;
    00052                 double sum = 0.0;
    00053                 for ( i = 0;i < w.length();i++ ) {sum += w ( i ) * Coms ( i )->evalpdflog ( val );}
    00054                 return log ( sum );
    00055         };
    00056
    00057 //Access methods
    00059         vec& _w() {return w;}
    00060 };
    00061
    00070 class mprod: public mpdf {
    00071 protected:
    00072         //
    00073         int n;
    00074         // pointers to epdfs
    00075         Array<epdf*> epdfs;
    00076         Array<mpdf*> mpdfs;
    00077         //
    00079         Array<ivec> rvinds;
    00081         Array<ivec> rvcinrv;
    00083         Array<ivec> rvcinds;
    00084 //      //! Indicate independence of its factors
    00085 //      bool independent;
    00086 //      //! Indicate internal creation of mpdfs which must be destroyed
    00087 //      bool intermpdfs;
    00088 public:
    00092         mprod ( Array<mpdf*> mFacs, bool overlap=false );
    00093
    00094         double evalpdflog ( const vec &val ) const {
    00095                 int i;
    00096                 double res = 0.0;
    00097                 for ( i = n - 1;i > 0;i++ ) {
    00098                         if ( rvcinds ( i ).length() > 0 )
    00099                                 {mpdfs ( i )->condition ( val ( rvcinds ( i ) ) );}
    00100                         // add logarithms
    00101                         res += epdfs ( i )->evalpdflog ( val ( rvinds ( i ) ) );
    00102                 }
    00103                 return res;
    00104         }
    00105         vec samplecond ( const vec &cond, double &ll ) {
    00106                 vec smp=zeros ( rv.count() );
    00107                 vec condi;
    00108                 vec smpi;
    00109                 ll = 0;
    00110                 for ( int i = ( n - 1 );i >= 0;i-- ) {
    00111                         if ( rvcinds ( i ).length() > 0 ) {
    00112                                 condi = zeros ( rvcinrv.length() + rvcinds.length() );
    00113                                 // copy data in condition
    00114                                 set_subvector ( condi,rvcinds ( i ), cond );
    00115                                 // copy data in already generated sample
    00116                                 set_subvector ( condi,rvcinrv ( i ), smp );
    00117
    00118                                 mpdfs ( i )->condition ( condi );
    00119                         }
    00120                         smpi = epdfs ( i )->sample();
    00121                         // copy contribution of this pdf into smp
    00122                         set_subvector ( smp,rvinds ( i ), smpi );
    00123                         // add ith likelihood contribution
    00124                         ll+=epdfs ( i )->evalpdflog ( smpi );
    00125                 }
    00126                 return smp;
    00127         }
    00128         mat samplecond ( const vec &cond, vec &ll, int N ) {
    00129                 mat Smp ( rv.count(),N );
    00130                 for ( int i=0;i<N;i++ ) {Smp.set_col ( i,samplecond ( cond,ll ( i ) ) );}
    00131                 return Smp;
    00132         }
    00133
    00134         ~mprod() {};
    00135 };
    00136
    00138 class eprod: public epdf {
    00139 protected:
    00141         Array<const epdf*> epdfs;
    00143         Array<ivec> rvinds;
    00144 public:
    00145         eprod ( const Array<const epdf*> epdfs0 ) : epdf ( RV() ),epdfs ( epdfs0 ),rvinds ( epdfs.length() ) {
    00146                 bool independent=true;
    00147                 for ( int i=0;i<epdfs.length();i++ ) {
    00148                         independent=rv.add ( epdfs ( i )->_rv() );
    00149                         it_assert_debug ( independent==true, "eprod:: given components are not independent ." );
    00150                         rvinds ( i ) = ( epdfs ( i )->_rv() ).dataind ( rv );
    00151                 }
    00152         }
    00153
    00154         vec mean() const {
    00155                 vec tmp ( rv.count() );
    00156                 for ( int i=0;i<epdfs.length();i++ ) {
    00157                         vec pom = epdfs ( i )->mean();
    00158                         set_subvector ( tmp,rvinds ( i ), pom );
    00159                 }
    00160                 return tmp;
    00161         }
    00162         vec sample() const {
    00163                 vec tmp ( rv.count() );
    00164                 for ( int i=0;i<epdfs.length();i++ ) {
    00165                         vec pom = epdfs ( i )->sample();
    00166                         set_subvector ( tmp,rvinds ( i ), pom );
    00167                 }
    00168                 return tmp;
    00169         }
    00170         double evalpdflog ( const vec &val ) const {
    00171                 double tmp=0;
    00172                 for ( int i=0;i<epdfs.length();i++ ) {
    00173                         tmp+=epdfs(i)->evalpdflog(val(rvinds(i)));
    00174                 }
    00175                 return tmp;
    00176         }
    00178         const epdf* operator () (int i) const {it_assert_debug(i<epdfs.length(),"wrong index");return epdfs(i);}
    00179 };
    00180
    00181
    00185 class mmix : public mpdf {
    00186 protected:
    00188         Array<mpdf*> Coms;
    00190         emix Epdf;
    00191 public:
    00193         mmix ( RV &rv, RV &rvc ) : mpdf ( rv, rvc ), Epdf ( rv ) {ep = &Epdf;};
    00195         void set_parameters ( const vec &w, const Array<mpdf*> &Coms ) {
    00196                 Array<epdf*> Eps ( Coms.length() );
    00197
    00198                 for ( int i = 0;i < Coms.length();i++ ) {
    00199                         Eps ( i ) = & ( Coms ( i )->_epdf() );
    00200                 }
    00201                 Epdf.set_parameters ( w, Eps );
    00202         };
    00203
    00204         void condition ( const vec &cond ) {
    00205                 for ( int i = 0;i < Coms.length();i++ ) {Coms ( i )->condition ( cond );}
    00206         };
    00207 };
    00208 #endif //MX_H