merger.h

00001
00013 #ifndef MERGER_H
00014 #define MERGER_H
00015 
00016
00017 #include "../estim/mixtures.h"
00018 #include "discrete.h"
00019
00020 namespace bdm {
00021 using std::string;
00022
00024 enum MERGER_METHOD {ARITHMETIC = 1, GEOMETRIC = 2, LOGNORMAL = 3};
00025
00040 class MergerBase: public root{
00041         public:
00043                 vec w;
00045                 Array<shared_ptr<epdf> > sources;
00047                 virtual void merge() NOT_IMPLEMENTED_VOID;
00050                 virtual void set_sources(const Array<shared_ptr<pdf> > &A) {
00051                         int n=A.length();
00052
00053                         bdm_assert(n>0,"merger has no sources to merge");
00054                         // check if weights are ok
00055                         if (w.length()!=n) {
00056                                 w = ones(n)/n;
00057                         }
00058
00059                         // check compatibility of sources -- no types are needed
00060                         int dim0 = A(0)->dimension();
00061                         for (int i=0; i<n; i++){
00062                                 const epdf *si=dynamic_cast<const epdf*>(A(i).get());
00063                                 if (si){
00064                                         sources(i) = shared_ptr<epdf>(const_cast<epdf*>(si));
00065                                 }
00066                                 bdm_assert(sources(i)->dimension()==dim0, "Merger: Incompatible dimensions of sources");
00067                                 if (sources(i)->isnamed()){
00068                                         //check if rv match
00069                                 }
00070                         }
00071                 }
00073                 virtual const epdf& merger()=0;
00074 };
00075
00077 template<class sq_T>
00078 class ENormMerger: public MergerBase{
00079         protected:
00081                 enorm<sq_T> iepdf;
00082         public:
00083                 ENormMerger():method(GEOMETRIC){};
00084                 MERGER_METHOD method;
00085                 void merge(){
00086                         int n = sources.length();
00087                         int dim = sources(0)->dimension();
00088                         sq_T Var(zeros(dim));
00089                         vec mea=zeros(dim);
00090                         // go through all sources
00091                         epdf * si;
00092                         for (int i=0; i<n; i++){
00093                                 si = dynamic_cast<epdf*>(sources(i).get());  // transform pdf into epdf
00094
00095                                 if (!si) {bdm_error("Can't merger this type of pdf: " + sources(i)->to_string());}
00096
00097                                 sq_T Ci = si->covariance();
00098                                 sq_T wCi = Ci; wCi*=w(i);
00099                                 vec mi = si->mean();
00100                                 switch (method) {
00101                                         case ARITHMETIC:{
00102                                                 Var += wCi;
00103                                                 Var.opupdt(mi,w(i)); // + mean * mean'
00104                                                 mea += w(i)*mi;
00105                                                 break;
00106                                         }
00107                                         case GEOMETRIC:{
00108                                                 Var += wCi;
00109                                                 sq_T iCi;
00110                                                 Ci.inv(iCi);
00111                                                 mea += iCi.to_mat()*w(i)*mi;
00112                                                 break;
00113                                         }
00114                                         default:
00115                                                 bdm_error("Method not implemneted");
00116                                 }
00117                         }
00118                         // post pocessing
00119                         switch (method) {
00120                                 case ARITHMETIC:
00121                                         iepdf._R()=Var;
00122                                         iepdf._R().opupdt(mea,-1.0); // Var - mean * mean'
00123                                         iepdf._mu()=mea;
00124                                         break;
00125                                 case GEOMETRIC:
00126                                         iepdf._R()=Var;
00127                                         iepdf._mu() = Var.to_mat()*mea; // mean=sqrt(Var)*mea
00128                                         break;
00129                                 default:
00130                                         bdm_error("Method not implemneted");
00131                         }
00132                 }
00134                 const epdf& merger(){return iepdf;}
00135 };
00136
00157 class MergerDiscrete : public MergerBase {
00158 protected:
00160         Array<shared_ptr< pdf > > part_sources;
00161
00163     Array<datalink_m2e*> dls;
00164
00166     Array<RV> rvzs;
00167
00169     Array<datalink_m2e*> zdls;
00170
00172     int Npoints;
00173
00175     int Nsources;
00176
00178     MERGER_METHOD METHOD;
00180     static const MERGER_METHOD DFLT_METHOD;
00181
00183     double beta;
00185     static const double DFLT_beta;
00186
00188     eEmp eSmp;
00189
00191     bool DBG;
00192
00194     it_file* dbg_file;
00195 public:
00198
00200     MergerDiscrete () : Npoints ( 0 ), Nsources ( 0 ), DBG ( false ), dbg_file ( 0 ) {
00201     }
00202
00204     MergerDiscrete ( const Array<shared_ptr<pdf> > &S );
00205
00207     void set_sources ( const Array<shared_ptr<pdf> > &PartSources );
00208
00210     void set_support ( rectangular_support &Sup );
00211
00213     void set_support ( discrete_support &Sup ) {
00214         Npoints = Sup.points();
00215         eSmp.set_parameters ( Sup._Spoints() );
00216         eSmp.validate();
00217     }
00218
00219
00221     void set_debug_file ( const string fname ) {
00222         if ( DBG ) delete dbg_file;
00223         dbg_file = new it_file ( fname );
00224         DBG = ( dbg_file != 0 );
00225     }
00226
00228     void set_method ( MERGER_METHOD MTH = GEOMETRIC, double beta0 = 1.2 ) {
00229         METHOD = MTH;
00230         beta = beta0;
00231     }
00233     void set_support ( const epdf &overall, int N ) {
00234         eSmp.set_statistics ( overall, N );
00235         Npoints = N;
00236     }
00237
00239     virtual ~MergerDiscrete() {
00240         for ( int i = 0; i < Nsources; i++ ) {
00241             delete dls ( i );
00242             delete zdls ( i );
00243         }
00244         if ( DBG ) delete dbg_file;
00245     };
00247
00250
00252     void merge ();
00253
00255     vec merge_points ( mat &lW );
00256
00258
00261
00263     eEmp& _Smp() {
00264         return eSmp;
00265     }
00266
00267         eEmp & merger() {return eSmp;}
00268
00282     void from_setting ( const Setting& set );
00283
00284     void to_setting  (Setting  &set) const ;
00285
00286     void validate() ;
00288 };
00289 UIREGISTER ( MergerDiscrete );
00290 SHAREDPTR ( MergerDiscrete );
00291
00293 class merger_mix : public MergerDiscrete {
00294 protected:
00296     MixEF Mix;
00298     int Ncoms;
00300     double effss_coef;
00302     int stop_niter;
00303
00305     static const int DFLT_Ncoms;
00307     static const double DFLT_effss_coef;
00308
00309 public:
00312     merger_mix () : Ncoms ( 0 ), effss_coef ( 0 ), stop_niter ( 0 ) { }
00313
00314     merger_mix ( const Array<shared_ptr<pdf> > &S ) :
00315         Ncoms ( 0 ), effss_coef ( 0 ), stop_niter ( 0 ) {
00316         set_sources ( S );
00317     }
00318
00320     void set_sources ( const Array<shared_ptr<pdf> > &S ) {
00321         MergerDiscrete::set_sources ( S );
00322         //Nsources = S.length();
00323     }
00324
00326     void set_parameters ( int Ncoms0 = DFLT_Ncoms, double effss_coef0 = DFLT_effss_coef ) {
00327         Ncoms = Ncoms0;
00328         effss_coef = effss_coef0;
00329     }
00331
00334
00336     void merge ();
00337
00339     vec sample () const {
00340         return Mix.posterior().sample();
00341     }
00343     double evallog ( const vec &yt ) const {
00344         vec dtf = ones ( yt.length() + 1 );
00345         dtf.set_subvector ( 0, yt );
00346         return Mix.logpred ( dtf, vec(0) );
00347     }
00349
00353     MixEF& _Mix() {
00354         return Mix;
00355     }
00357     emix* proposal() {
00358         emix* tmp = Mix.epredictor();
00359         tmp->set_rv ( merger()._rv() );
00360         return tmp;
00361     }
00381     void from_setting ( const Setting& set );
00382
00383     void to_setting  (Setting  &set) const;
00384     void validate();
00385
00386 //      const emix &merger(){return *Mpred;}
00388 };
00389 UIREGISTER ( merger_mix );
00390 SHAREDPTR ( merger_mix );
00391
00392 }
00393
00394 #endif // MER_H