- 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
    00024 enum MERGER_METHOD {ARITHMETIC = 0, GEOMETRIC = 1, LOGNORMAL = 2};
    00025
    00045 class merger_base : public compositepdf, public epdf
    00046 {
    00047         protected:
    00049                 Array<datalink_m2e*> dls;
    00051                 Array<RV> rvzs;
    00053                 Array<datalink_m2e*> zdls;
    00055                 int Ns;
    00056
    00058                 MERGER_METHOD METHOD;
    00060                 double beta;
    00061
    00063                 eEmp eSmp;
    00064
    00066                 bool DBG;
    00067
    00069                 it_file* dbg;
    00070         public:
    00073
    00075                 merger_base () : compositepdf() {};
    00077                 merger_base (const Array<mpdf*> &S) {set_sources (S);};
    00079                 void set_sources (const Array<mpdf*> &Sources) {
    00080                         compositepdf::set_elements (Sources);
    00081                         //set sizes
    00082                         dls.set_size (Sources.length());
    00083                         rvzs.set_size (Sources.length());
    00084                         zdls.set_size (Sources.length());
    00085
    00086                         rv = getrv (/* checkoverlap = */ false);
    00087                         RV rvc; setrvc (rv, rvc);  // Extend rv by rvc!
    00088                         // join rv and rvc - see descriprion
    00089                         rv.add (rvc);
    00090                         // get dimension
    00091                         dim = rv._dsize();
    00092
    00093                         // create links between sources and common rv
    00094                         RV xytmp;
    00095                         for (int i = 0;i < mpdfs.length();i++) {
    00096                                 //Establich connection between mpdfs and merger
    00097                                 dls (i) = new datalink_m2e;
    00098                                 dls (i)->set_connection (mpdfs (i)->_rv(), mpdfs (i)->_rvc(), rv);
    00099
    00100                                 // find out what is missing in each mpdf
    00101                                 xytmp = mpdfs (i)->_rv();
    00102                                 xytmp.add (mpdfs (i)->_rvc());
    00103                                 // z_i = common_rv-xy
    00104                                 rvzs (i) = rv.subt (xytmp);
    00105                                 //establish connection between extension (z_i|x,y)s and common rv
    00106                                 zdls (i) = new datalink_m2e; zdls (i)->set_connection (rvzs (i), xytmp, rv) ;
    00107                         };
    00108                         //Set Default values of parameters
    00109                         DBG = false;
    00110                         METHOD = LOGNORMAL;
    00111                         beta = 2.0;
    00112                 }
    00114                 void set_debug_file (const string fname) { if (DBG) delete dbg; dbg = new it_file (fname); if (dbg) DBG = true;}
    00116                 void set_parameters (double beta0, MERGER_METHOD MTH) {beta = beta0;METHOD = MTH;}
    00118                 void set_support(epdf &overall, int N){eSmp.set_statistics(&overall,N);Ns=N;}
    00119
    00121                 virtual ~merger_base() {
    00122                         for (int i = 0; i < mpdfs.length(); i++) {
    00123                                 delete dls (i);
    00124                                 delete zdls (i);
    00125                         }
    00126                         if (DBG) delete dbg;
    00127                 };
    00129
    00132
    00134                 void merge () {
    00135                         if (eSmp._w().length() ==0) {it_error("Empty support points use set_support" );}
    00136                         //check if sources overlap:
    00137                         bool OK = true;
    00138                         for (int i = 0;i < mpdfs.length(); i++) {
    00139                                 OK &= (rvzs (i)._dsize() == 0); // z_i is empty
    00140                                 OK &= (mpdfs (i)->_rvc()._dsize() == 0); // y_i is empty
    00141                         }
    00142
    00143                         if (OK) {
    00144                                 mat lW = zeros (mpdfs.length(), eSmp._w().length());
    00145
    00146                                 vec emptyvec (0);
    00147                                 for (int i = 0; i < mpdfs.length(); i++) {
    00148                                         for (int j = 0; j < eSmp._w().length(); j++) {
    00149                                                 lW (i, j) = mpdfs (i)->evallogcond (eSmp._samples() (j), emptyvec);
    00150                                         }
    00151                                 }
    00152
    00153                                 vec wtmp = exp (merge_points (lW));
    00154                                 //renormalize
    00155                                 eSmp._w() = wtmp / sum (wtmp);
    00156                         } else {
    00157                                 it_error("Sources are not compatible - use merger_mix");
    00158                         }
    00159                         ;
    00160                 };
    00161
    00162
    00164                 vec merge_points (mat &lW);
    00165
    00166
    00169                 vec mean() const {
    00170                         const Vec<double> &w = eSmp._w();
    00171                         const Array<vec> &S = eSmp._samples();
    00172                         vec tmp = zeros (dim);
    00173                         for (int i = 0; i < Ns; i++) {
    00174                                 tmp += w (i) * S (i);
    00175                         }
    00176                         return tmp;
    00177                 }
    00178                 mat covariance() const {
    00179                         const vec &w = eSmp._w();
    00180                         const Array<vec> &S = eSmp._samples();
    00181
    00182                         vec mea = mean();
    00183
    00184                         cout << sum (w) << "," << w*w << endl;
    00185
    00186                         mat Tmp = zeros (dim, dim);
    00187                         for (int i = 0; i < Ns; i++) {
    00188                                 Tmp += w (i) * outer_product (S (i), S (i));
    00189                         }
    00190                         return Tmp -outer_product (mea, mea);
    00191                 }
    00192                 vec variance() const {
    00193                         const vec &w = eSmp._w();
    00194                         const Array<vec> &S = eSmp._samples();
    00195
    00196                         vec tmp = zeros (dim);
    00197                         for (int i = 0; i < Ns; i++) {
    00198                                 tmp += w (i) * pow (S (i), 2);
    00199                         }
    00200                         return tmp -pow (mean(), 2);
    00201                 }
    00203
    00206
    00208                 eEmp& _Smp() {return eSmp;}
    00209
    00211 };
    00212
    00213 class merger_mix : public merger_base
    00214 {
    00215         protected:
    00217                 MixEF Mix;
    00219                 int Nc;
    00221                 double effss_coef;
    00222
    00223         public:
    00225                 merger_mix () {};
    00226                 merger_mix (const Array<mpdf*> &S) {set_sources(S);};
    00228                 void set_sources (const Array<mpdf*> &S) {
    00229                         merger_base::set_sources(S);
    00230                         //Set Default values of parameters
    00231                         Ns = 100;
    00232                         Nc = 10;
    00233                         Mix.set_method (EM);
    00234                 }
    00236                 void set_parameters (double beta0, int Ns0, int Nc0, double effss_coef0 = 0.5) {
    00237                         beta = beta0;
    00238                         Ns = Ns0;
    00239                         Nc = Nc0;
    00240                         effss_coef = effss_coef0;
    00241                         eSmp.set_parameters (Ns0, false);
    00242                 }
    00244                 void init() {
    00245                         Array<vec> Smps (Ns);
    00246                         //Gibbs sampling
    00247                         for (int i = 0;i < mpdfs.length();i++) {Smps (i) = zeros (0);}
    00248                 }
    00250                 void merge ();
    00251
    00254                 vec sample () const { return Mix.posterior().sample();}
    00255                 double evallog (const vec &dt) const {
    00256                         vec dtf = ones (dt.length() + 1);
    00257                         dtf.set_subvector (0, dt);
    00258                         return Mix.logpred (dtf);
    00259                 }
    00260                 vec mean() const {
    00261                         const Vec<double> &w = eSmp._w();
    00262                         const Array<vec> &S = eSmp._samples();
    00263                         vec tmp = zeros (dim);
    00264                         for (int i = 0; i < Ns; i++) {
    00265                                 tmp += w (i) * S (i);
    00266                         }
    00267                         return tmp;
    00268                 }
    00269                 mat covariance() const {
    00270                         const vec &w = eSmp._w();
    00271                         const Array<vec> &S = eSmp._samples();
    00272
    00273                         vec mea = mean();
    00274
    00275                         cout << sum (w) << "," << w*w << endl;
    00276
    00277                         mat Tmp = zeros (dim, dim);
    00278                         for (int i = 0; i < Ns; i++) {
    00279                                 Tmp += w (i) * outer_product (S (i), S (i));
    00280                         }
    00281                         return Tmp -outer_product (mea, mea);
    00282                 }
    00283                 vec variance() const {
    00284                         const vec &w = eSmp._w();
    00285                         const Array<vec> &S = eSmp._samples();
    00286
    00287                         vec tmp = zeros (dim);
    00288                         for (int i = 0; i < Ns; i++) {
    00289                                 tmp += w (i) * pow (S (i), 2);
    00290                         }
    00291                         return tmp -pow (mean(), 2);
    00292                 }
    00293
    00295                 MixEF& _Mix() {return Mix;}
    00297                 emix* proposal() {emix* tmp = Mix.epredictor(); tmp->set_rv (rv); return tmp;}
    00299                 eEmp& _Smp() {return eSmp;}
    00300 };
    00301
    00302 }
    00303
    00304 #endif // MER_H