- Main Page
- Related Pages
- Modules
- Classes
- Files

- File List
- File Members

# merger.h

Go to the documentation of this file.

    00001
    00013 #ifndef MER_H
    00014 #define MER_H
    00015
    00016
    00017 #include "mixef.h"
    00018
    00019 namespace bdm
    00020 {
    00021         using std::string;
    00022
    00029         class merger : public compositepdf, public epdf
    00030         {
    00031                 protected:
    00033                         MixEF Mix;
    00035                         Array<datalink_m2e*> dls;
    00037                         Array<RV> rvzs;
    00039                         Array<datalink_m2e*> zdls;
    00040
    00042                         int Ns;
    00044                         int Nc;
    00046                         double beta;
    00048                         eEmp eSmp;
    00050                         double effss_coef;
    00051
    00053                         bool DBG;
    00055                         it_file* dbg;
    00057                         bool fix_smp;
    00058                 public:
    00060                         merger ( const Array<mpdf*> &S ) :
    00061                                         compositepdf ( S ), epdf ( ),
    00062                                         Mix ( Array<BMEF*> ( 0 ),vec ( 0 ) ), dls ( n ), rvzs ( n ), zdls ( n ), eSmp()
    00063                         {
    00064                                 RV ztmp;
    00065                                 rv = getrv ( false );
    00066                                 RV rvc; setrvc ( rv,rvc ); // Extend rv by rvc!
    00067                                 rv.add ( rvc );
    00068                                 // get dimension
    00069                                 dim = rv._dsize();
    00070
    00071                                 for ( int i=0;i<n;i++ )
    00072                                 {
    00073                                         //Establich connection between mpdfs and merger
    00074                                         dls ( i ) = new datalink_m2e;dls ( i )->set_connection ( mpdfs ( i )->_rv(), mpdfs ( i )->_rvc(), rv );
    00075                                         // find out what is missing in each mpdf
    00076                                         ztmp= mpdfs ( i )->_rv();
    00077                                         ztmp.add ( mpdfs ( i )->_rvc() );
    00078                                         rvzs ( i ) =rv.subt ( ztmp );
    00079                                         zdls ( i ) = new datalink_m2e; zdls ( i )->set_connection ( rvzs ( i ), ztmp, rv ) ;
    00080                                 };
    00081                                 //Set Default values of parameters
    00082                                 beta=2.0;
    00083                                 Ns=100;
    00084                                 Nc=10;
    00085                                 Mix.set_method ( EM );
    00086                                 DBG = false;
    00087                                 fix_smp = false;
    00088                         }
    00090                         void debug_file ( const string fname ) { if ( DBG ) delete dbg; dbg = new it_file ( fname ); if ( dbg ) DBG=true;}
    00092                         void set_parameters ( double beta0, int Ns0, int Nc0, double effss_coef0=0.5 ) {beta=beta0;
    00093                                 Ns=Ns0;
    00094                                 Nc=Nc0;
    00095                                 effss_coef=effss_coef0;
    00096                                 eSmp.set_parameters ( Ns0,false );
    00097                         }
    00098                         void set_grid ( Array<vec> &XYZ )
    00099                         {
    00100                                 int dim=XYZ.length(); ivec szs ( dim );
    00101                                 for(int i=0; i<dim;i++){szs=XYZ(i).length();}
    00102                                 Ns=prod(szs);
    00103                                 eSmp.set_parameters(Ns,false);
    00104                                 Array<vec> &samples=eSmp._samples();
    00105                                 eSmp._w()=ones(Ns)/Ns;
    00106
    00107                                 //set samples
    00108                                 ivec is=zeros_i(dim);//indeces of dimensions in for cycle;
    00109                                 vec smpi(dim);
    00110                                 for(int i=0; i<Ns; i++){
    00111                                         for(int j=0; j<dim; j++){smpi(j)=XYZ(j)(is(j)); /* jty vector*/ }
    00112                                         samples(i)=smpi;
    00113                                         // shift indeces
    00114                                         for (int j=0;j<dim;j++){
    00115                                                 if (is(j)==szs(j)-1) { //j-th index is full
    00116                                                         is(j)=0; //shift back
    00117                                                         is(j+1)++; //increase th next dimension;
    00118                                                         if (is(j+1)<szs(j+1)-1) break;
    00119                                                 } else {
    00120                                                         is(j)++; break;
    00121                                                 }
    00122                                         }
    00123                                 }
    00124
    00125                                 fix_smp = true;
    00126                         }
    00128                         void init()
    00129                         {
    00130                                 Array<vec> Smps ( n );
    00131                                 //Gibbs sampling
    00132                                 for ( int i=0;i<n;i++ ) {Smps ( i ) =zeros ( 0 );}
    00133                         }
    00135                         void merge ( const epdf* g0 );
    00137                         void merge () {merge ( & ( Mix.posterior() ) );};
    00138
    00140                         vec lognorm_merge ( mat &lW );
    00143                         vec sample ( ) const { return Mix.posterior().sample();}
    00144                         double evallog ( const vec &dt ) const
    00145                         {
    00146                                 vec dtf=ones ( dt.length() +1 );
    00147                                 dtf.set_subvector ( 0,dt );
    00148                                 return Mix.logpred ( dtf );
    00149                         }
    00150                         vec mean() const
    00151                         {
    00152                                 const Vec<double> &w = eSmp._w();
    00153                                 const Array<vec> &S = eSmp._samples();
    00154                                 vec tmp=zeros ( dim );
    00155                                 for ( int i=0; i<Ns; i++ )
    00156                                 {
    00157                                         tmp+=w ( i ) *S ( i );
    00158                                 }
    00159                                 return tmp;
    00160                         }
    00161                         mat covariance() const
    00162                         {
    00163                                 const vec &w = eSmp._w();
    00164                                 const Array<vec> &S = eSmp._samples();
    00165
    00166                                 vec mea = mean();
    00167
    00168                                 cout << sum ( w ) << "," << w*w <<endl;
    00169
    00170                                 mat Tmp=zeros ( dim, dim );
    00171                                 for ( int i=0; i<Ns; i++ )
    00172                                 {
    00173                                         Tmp+=w ( i ) *outer_product ( S ( i ), S ( i ) );
    00174                                 }
    00175                                 return Tmp-outer_product ( mea,mea );
    00176                         }
    00177                         vec variance() const
    00178                         {
    00179                                 const vec &w = eSmp._w();
    00180                                 const Array<vec> &S = eSmp._samples();
    00181
    00182                                 vec tmp=zeros ( dim );
    00183                                 for ( int i=0; i<Ns; i++ )
    00184                                 {
    00185                                         tmp+=w ( i ) *pow ( S ( i ),2 );
    00186                                 }
    00187                                 return tmp-pow ( mean(),2 );
    00188                         }
    00190                         virtual ~merger()
    00191                         {
    00192                                 for ( int i=0; i<n; i++ )
    00193                                 {
    00194                                         delete dls ( i );
    00195                                         delete zdls ( i );
    00196                                 }
    00197                                 if ( DBG ) delete dbg;
    00198                         };
    00199
    00201                         MixEF& _Mix() {return Mix;}
    00203                         emix* proposal() {emix* tmp=Mix.epredictor(); tmp->set_rv(rv); return tmp;}
    00205                         eEmp& _Smp() {return eSmp;}
    00206         };
    00207
    00208 }
    00209
    00210 #endif // MER_H