- Main Page
- Related Pages
- Classes
- Files

- File List
- File Members

# bdmbase.h

Go to the documentation of this file.

    00001
    00013 #ifndef BDMBASE_H
    00014 #define BDMBASE_H
    00015
    00016 #include <map>
    00017
    00018 #include "../itpp_ext.h"
    00019 #include "../bdmroot.h"
    00020 #include "../shared_ptr.h"
    00021 #include "user_info.h"
    00022
    00023 using namespace libconfig;
    00024 using namespace itpp;
    00025 using namespace std;
    00026
    00027 namespace bdm {
    00029 class str {
    00030 public:
    00032         ivec ids;
    00034         ivec times;
    00036         str ( ivec ids0, ivec times0 ) : ids ( ids0 ), times ( times0 ) {
    00037                 bdm_assert ( times0.length() == ids0.length(), "Incompatible input" );
    00038         };
    00039 };
    00040
    00079 class RV : public root {
    00080 private:
    00081         typedef std::map<string, int> str2int_map;
    00082
    00084         static ivec SIZES;
    00086         static Array<string> NAMES;
    00087
    00089         const static int BUFFER_STEP;
    00090
    00092         static str2int_map MAP;
    00093
    00094 public:
    00095
    00096 protected:
    00098         int dsize;
    00100         int len;
    00102         ivec ids;
    00104         ivec times;
    00105
    00106 private:
    00107         enum enum_dummy {dummy};
    00108
    00110         void init ( const Array<std::string> &in_names, const ivec &in_sizes, const ivec &in_times );
    00111         int init ( const string &name, int size );
    00114         explicit RV ( const ivec &ids0, enum_dummy dum ) : dsize ( 0 ), len ( ids0.length() ), ids ( ids0 ), times ( zeros_i ( ids0.length() ) ) {
    00115                 dsize = countsize();
    00116         }
    00117 public:
    00120
    00122         RV ( const Array<std::string> &in_names, const ivec &in_sizes, const ivec &in_times ) {
    00123                 init ( in_names, in_sizes, in_times );
    00124         }
    00125
    00127         RV ( const Array<std::string> &in_names, const ivec &in_sizes ) {
    00128                 init ( in_names, in_sizes, zeros_i ( in_names.length() ) );
    00129         }
    00130
    00132         RV ( const Array<std::string> &in_names ) {
    00133                 init ( in_names, ones_i ( in_names.length() ), zeros_i ( in_names.length() ) );
    00134         }
    00135
    00137         RV() : dsize ( 0 ), len ( 0 ), ids ( 0 ), times ( 0 ) {}
    00138
    00140         RV ( string name, int sz, int tm = 0 );
    00141
    00142         // compiler-generated copy constructor is used
    00144
    00147
    00149         friend std::ostream &operator<< ( std::ostream &os, const RV &rv );
    00150
    00151         string to_string() {ostringstream o; o << this; return o.str();}
    00152
    00154         int _dsize() const {
    00155                 return dsize;
    00156         }
    00157
    00159         const ivec& _ids() const {
    00160                 return ids;
    00161         }
    00162
    00164         int countsize() const;
    00166         ivec cumsizes() const;
    00168         int length() const {
    00169                 return len;
    00170         }
    00171         int id ( int at ) const {
    00172                 return ids ( at );
    00173         }
    00174         int size ( int at ) const {
    00175                 return SIZES ( ids ( at ) );
    00176         }
    00177         int time ( int at ) const {
    00178                 return times ( at );
    00179         }
    00180         std::string name ( int at ) const {
    00181                 return NAMES ( ids ( at ) );
    00182         }
    00184         std::string scalarname ( int scalat ) const {
    00185                 bdm_assert(scalat<dsize,"Wrong input index");
    00186                 int id=0;
    00187                 int scalid=0;
    00188                 while (scalid+SIZES(ids(id))<=scalat)  { scalid+=SIZES(ids(id)); id++;};
    00189                 //now id is the id of variable of interest
    00190                 if (size(id)==1)
    00191                         return  NAMES ( ids ( id ) );
    00192                 else
    00193                         return  NAMES ( ids ( id ) )+ "_" + num2str(scalat-scalid);
    00194
    00195         }
    00196         void set_time ( int at, int time0 ) {
    00197                 times ( at ) = time0;
    00198         }
    00200
    00203
    00205         ivec findself ( const RV &rv2 ) const;
    00207         ivec findself_ids ( const RV &rv2 ) const;
    00209         bool equal ( const RV &rv2 ) const;
    00211         bool add ( const RV &rv2 );
    00213         RV subt ( const RV &rv2 ) const;
    00215         RV subselect ( const ivec &ind ) const;
    00216
    00218         RV operator() ( const ivec &ind ) const {
    00219                 return subselect ( ind );
    00220         }
    00221
    00223         RV operator() ( int di1, int di2 ) const;
    00224
    00226         void t_plus ( int delta );
    00228
    00231         RV remove_time() const {
    00232                 return RV ( unique ( ids ), dummy );
    00233         }
    00235         RV copy_t(int dt) const {
    00236                 RV tmp=*this;
    00237                 tmp.t_plus(dt);
    00238                 return tmp;
    00239         }
    00241         RV expand_delayes() const {
    00242                 RV rvt = this->remove_time(); //rv at t=0
    00243                 RV tmp = rvt;
    00244                 int td = mint();
    00245                 for ( int i = -1; i >= td; i-- ) {
    00246                         rvt.t_plus ( -1 );
    00247                         tmp.add ( rvt ); //shift u1
    00248                 }
    00249                 return tmp;
    00250         }
    00252
    00255
    00257         str tostr() const;
    00260         ivec dataind ( const RV &crv ) const;
    00262         ivec dataind_part ( const RV &crv ) const;
    00265         void dataind ( const RV &rv2, ivec &selfi, ivec &rv2i ) const;
    00267         int mint() const {
    00268                         return times.length()>0 ? min (times) : 0;
    00269         }
    00271
    00284         void from_setting ( const Setting &set );
    00285
    00286         // TODO dodelat void to_setting( Setting &set ) const;
    00287
    00289         static void clear_all();
    00291         string show_all();
    00292
    00293 };
    00294 UIREGISTER ( RV );
    00295 SHAREDPTR ( RV );
    00296
    00298 RV concat ( const RV &rv1, const RV &rv2 );
    00299
    00301 class fnc : public root {
    00302 protected:
    00304         int dimy;
    00305 public:
    00307         fnc() {};
    00309         virtual vec eval ( const vec &cond ) {
    00310                 return vec ( 0 );
    00311         };
    00312
    00314         virtual void condition ( const vec &val ) {};
    00315
    00317         int dimension() const {
    00318                 return dimy;
    00319         }
    00320 };
    00321
    00322 class mpdf;
    00323
    00325
    00326 class epdf : public root {
    00327 protected:
    00329         int dim;
    00331         RV rv;
    00332
    00333 public:
    00345         epdf() : dim ( 0 ), rv() {};
    00346         epdf ( const epdf &e ) : dim ( e.dim ), rv ( e.rv ) {};
    00347         epdf ( const RV &rv0 ) : dim ( rv0._dsize() ) {
    00348                 set_rv ( rv0 );
    00349         };
    00350         void set_parameters ( int dim0 ) {
    00351                 dim = dim0;
    00352         }
    00354
    00357
    00359         virtual vec sample() const {
    00360                 bdm_error ( "not implemented" );
    00361                 return vec();
    00362         }
    00363
    00365         virtual mat sample_m ( int N ) const;
    00366
    00369         virtual double evallog ( const vec &val ) const {
    00370                 bdm_error ( "not implemented" );
    00371                 return 0.0;
    00372         }
    00373
    00375         virtual vec evallog_m ( const mat &Val ) const;
    00376
    00378         virtual vec evallog_m ( const Array<vec> &Avec ) const;
    00379
    00381         virtual shared_ptr<mpdf> condition ( const RV &rv ) const;
    00382
    00384         virtual shared_ptr<epdf> marginal ( const RV &rv ) const;
    00385
    00387         virtual vec mean() const {
    00388                 bdm_error ( "not implemneted" );
    00389                 return vec();
    00390         }
    00391
    00393         virtual vec variance() const {
    00394                 bdm_error ( "not implemneted" );
    00395                 return vec();
    00396         }
    00397
    00399         virtual void qbounds ( vec &lb, vec &ub, double percentage = 0.95 ) const {
    00400                 vec mea = mean();
    00401                 vec std = sqrt ( variance() );
    00402                 lb = mea - 2 * std;
    00403                 ub = mea + 2 * std;
    00404         };
    00406
    00412
    00414         void set_rv ( const RV &rv0 ) {
    00415                 rv = rv0;
    00416         }
    00417
    00419         bool isnamed() const {
    00420                 bool b = ( dim == rv._dsize() );
    00421                 return b;
    00422         }
    00423
    00425         const RV& _rv() const {
    00426                 //bdm_assert_debug ( isnamed(), "" );
    00427                 return rv;
    00428         }
    00430
    00433
    00435         int dimension() const {
    00436                 return dim;
    00437         }
    00445         void from_setting ( const Setting &set ) {
    00446                 shared_ptr<RV> r = UI::build<RV> ( set, "rv", UI::optional );
    00447                 if ( r ) {
    00448                         set_rv ( *r );
    00449                 }
    00450         }
    00451
    00452 };
    00453 SHAREDPTR ( epdf );
    00454
    00456 class mpdf : public root {
    00457 protected:
    00459         int dimc;
    00461         RV rvc;
    00462 private:
    00464         epdf* ep;
    00465
    00466 protected:
    00468         void set_ep ( epdf &iepdf ) {
    00469                 ep = &iepdf;
    00470         }
    00472         void set_ep ( epdf *iepdfp ) {
    00473                 ep = iepdfp;
    00474         }
    00475
    00476 public:
    00479
    00480         mpdf() : dimc ( 0 ), rvc(), ep ( NULL ) { }
    00481
    00482         mpdf ( const mpdf &m ) : dimc ( m.dimc ), rvc ( m.rvc ), ep ( m.ep ) { }
    00484
    00487
    00489         virtual vec samplecond ( const vec &cond ) {
    00490                 bdm_error ( "Not implemented" );
    00491                 return vec();
    00492         }
    00493
    00495         virtual mat samplecond_m ( const vec &cond, int N );
    00496
    00498         virtual double evallogcond ( const vec &dt, const vec &cond ) {
    00499                 bdm_error ( "Not implemented" );
    00500                 return 0.0;
    00501         }
    00502
    00504         virtual vec evallogcond_m ( const mat &Dt, const vec &cond ) {
    00505                 vec v ( Dt.cols() );
    00506                 for ( int i = 0; i < Dt.cols(); i++ ) {
    00507                         v ( i ) = evallogcond ( Dt.get_col ( i ), cond );
    00508                 }
    00509                 return v;
    00510         }
    00511
    00513         virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond ) {
    00514                 bdm_error ( "Not implemented" );
    00515                 return vec();
    00516         }
    00517
    00520
    00521         const RV& _rv() const {
    00522                 return ep->_rv();
    00523         }
    00524         const RV& _rvc() const {
    00525                 return rvc;
    00526         }
    00527
    00528         int dimension() const {
    00529                 return ep->dimension();
    00530         }
    00531         int dimensionc() {
    00532                 return dimc;
    00533         }
    00534
    00544         void from_setting ( const Setting &set );
    00546
    00549         void set_rvc ( const RV &rvc0 ) {
    00550                 rvc = rvc0;
    00551         }
    00552         void set_rv ( const RV &rv0 ) {
    00553                 ep->set_rv ( rv0 );
    00554         }
    00555         bool isnamed() {
    00556                 return ( ep->isnamed() ) && ( dimc == rvc._dsize() );
    00557         }
    00559 };
    00560 SHAREDPTR ( mpdf );
    00561
    00563 template <class EPDF>
    00564 class mpdf_internal: public mpdf {
    00565 protected :
    00567         EPDF iepdf;
    00568 public:
    00570         mpdf_internal() : mpdf(), iepdf() {
    00571                 set_ep ( iepdf );
    00572         }
    00573
    00576         virtual void condition ( const vec &cond ) {
    00577                 bdm_error ( "Not implemented" );
    00578         }
    00579
    00581         EPDF& e() {
    00582                 return iepdf;
    00583         }
    00584
    00586         vec samplecond ( const vec &cond );
    00588         double evallogcond ( const vec &val, const vec &cond );
    00590         virtual vec evallogcond_m ( const mat &Dt, const vec &cond );
    00592         virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond );
    00594         virtual mat samplecond_m ( const vec &cond, int N );
    00595 };
    00596
    00622 class datalink {
    00623 protected:
    00625         int downsize;
    00626
    00628         int upsize;
    00629
    00631         ivec v2v_up;
    00632
    00633 public:
    00635         datalink() : downsize ( 0 ), upsize ( 0 ) { }
    00636
    00638         datalink ( const RV &rv, const RV &rv_up ) {
    00639                 set_connection ( rv, rv_up );
    00640         }
    00641
    00643         virtual void set_connection ( const RV &rv, const RV &rv_up );
    00644
    00646         virtual void set_connection ( int ds, int us, const ivec &upind );
    00647
    00649         vec pushdown ( const vec &val_up ) {
    00650                 vec tmp ( downsize );
    00651                 filldown ( val_up, tmp );
    00652                 return tmp;
    00653         }
    00655         virtual void filldown ( const vec &val_up, vec &val_down ) {
    00656                 bdm_assert_debug ( upsize == val_up.length(), "Wrong val_up" );
    00657                 val_down = val_up ( v2v_up );
    00658         }
    00660         virtual void pushup ( vec &val_up, const vec &val ) {
    00661                 bdm_assert_debug ( downsize == val.length(), "Wrong val" );
    00662                 bdm_assert_debug ( upsize == val_up.length(), "Wrong val_up" );
    00663                 set_subvector ( val_up, v2v_up, val );
    00664         }
    00666         int _upsize() {
    00667                 return upsize;
    00668         }
    00670         int _downsize() {
    00671                 return downsize;
    00672         }
    00673 };
    00674
    00677 class datalink_part : public datalink {
    00678 protected:
    00680         ivec v2v_down;
    00681 public:
    00682         void set_connection ( const RV &rv, const RV &rv_up );
    00684         void filldown ( const vec &val_up, vec &val_down ) {
    00685                 set_subvector ( val_down, v2v_down, val_up ( v2v_up ) );
    00686         }
    00687 };
    00688
    00693 class datalink_buffered: public datalink_part {
    00694 protected:
    00696         vec history;
    00698         RV Hrv;
    00700         ivec h2v_down;
    00702         ivec h2v_hist;
    00704         ivec v2h_up;
    00705 public:
    00706
    00707         datalink_buffered() : datalink_part(), history ( 0 ), h2v_down ( 0 ), h2v_hist ( 0 ) {};
    00709         void step ( const vec &val_up ) {
    00710                 if ( v2h_up.length() > 0 ) {
    00711                         history.shift_right ( 0, v2h_up.length() );
    00712                   history.set_subvector ( 0, val_up(v2h_up) );
    00713                 }
    00714         }
    00716         vec pushdown ( const vec &val_up ) {
    00717                 vec tmp ( downsize );
    00718                 filldown ( val_up, tmp );
    00719                 return tmp;
    00720         }
    00721
    00722         void filldown ( const vec &val_up, vec &val_down ) {
    00723                 bdm_assert_debug ( val_down.length() >= downsize, "short val_down" );
    00724
    00725                 set_subvector ( val_down, v2v_down, val_up ( v2v_up ) ); // copy direct values
    00726                 set_subvector ( val_down, h2v_down, history ( h2v_hist ) ); // copy delayed values
    00727         }
    00728
    00729         void set_connection ( const RV &rv, const RV &rv_up ) {
    00730                 // create link between up and down
    00731                 datalink_part::set_connection ( rv, rv_up );
    00732
    00733                 // create rvs of history
    00734                 // we can store only what we get in rv_up - everything else is removed
    00735                 ivec valid_ids = rv.findself_ids ( rv_up );
    00736                 RV rv_hist = rv.subselect ( find ( valid_ids >= 0 ) ); // select only rvs that are in rv_up
    00737                 RV rv_hist0 =rv_hist.remove_time();  // these RVs will form history at time =0
    00738                 // now we need to know what is needed from Up
    00739                 rv_hist = rv_hist.expand_delayes(); // full regressor - including time 0
    00740                 Hrv=rv_hist.subt(rv_hist0);   // remove time 0
    00741                 history = zeros ( Hrv._dsize() );
    00742
    00743                 // decide if we need to copy val to history
    00744                 if (Hrv._dsize() >0) {
    00745                         v2h_up = rv_hist0.dataind(rv_up); // indeces of elements of rv_up to be copied
    00746                 } // else v2h_up is empty
    00747
    00748                 Hrv.dataind ( rv, h2v_hist, h2v_down );
    00749
    00750                 downsize = v2v_down.length() + h2v_down.length();
    00751                 upsize = v2v_up.length();
    00752         }
    00754         void set_history(const RV& rv1, const vec &hist0){
    00755                 bdm_assert(rv1._dsize()==hist0.length(),"hist is not compatible with given rv1");
    00756                 ivec ind_H;
    00757                 ivec ind_h0;
    00758                 Hrv.dataind(rv1, ind_H, ind_h0); // find indeces of rv in
    00759                 set_subvector(history, ind_H, hist0(ind_h0)); // copy given hist to appropriate places
    00760         }
    00761 };
    00762
    00764 class datalink_2to1_buffered {
    00765 protected:
    00766         datalink_buffered dl1;
    00767         datalink_buffered dl2;
    00768 public:
    00770         void set_connection ( const RV &rv, const RV &rv_up1, const RV &rv_up2 ) {
    00771                 dl1.set_connection ( rv, rv_up1 );
    00772                 dl2.set_connection ( rv, rv_up2 );
    00773         }
    00774         void filldown ( const vec &val1, const vec &val2, vec &val_down ) {
    00775                 bdm_assert_debug ( val_down.length() >= dl1._downsize() + dl2._downsize(), "short val_down" );
    00776                 dl1.filldown ( val1, val_down );
    00777                 dl2.filldown ( val2, val_down );
    00778         }
    00779         void step ( const vec &dt, const vec &ut ) {
    00780                 dl1.step ( dt );
    00781                 dl2.step ( ut );
    00782         }
    00783 };
    00784
    00785
    00786
    00788 class datalink_m2e: public datalink {
    00789 protected:
    00791         int condsize;
    00792
    00794         ivec v2c_up;
    00795
    00797         ivec v2c_lo;
    00798
    00799 public:
    00801         datalink_m2e() : condsize ( 0 ) { }
    00802
    00804         void set_connection ( const RV &rv, const RV &rvc, const RV &rv_up );
    00805
    00807         vec get_cond ( const vec &val_up );
    00808
    00810         void pushup_cond ( vec &val_up, const vec &val, const vec &cond );
    00811 };
    00812
    00815 class datalink_m2m: public datalink_m2e {
    00816 protected:
    00818         ivec c2c_up;
    00820         ivec c2c_lo;
    00821
    00822 public:
    00824         datalink_m2m() {};
    00826         void set_connection ( const RV &rv, const RV &rvc, const RV &rv_up, const RV &rvc_up ) {
    00827                 datalink_m2e::set_connection ( rv, rvc, rv_up );
    00828                 //establish c2c connection
    00829                 rvc.dataind ( rvc_up, c2c_lo, c2c_up );
    00830                 bdm_assert_debug ( c2c_lo.length() + v2c_lo.length() == condsize, "cond is not fully given" );
    00831         }
    00832
    00834         vec get_cond ( const vec &val_up, const vec &cond_up ) {
    00835                 vec tmp ( condsize );
    00836                 set_subvector ( tmp, v2c_lo, val_up ( v2c_up ) );
    00837                 set_subvector ( tmp, c2c_lo, cond_up ( c2c_up ) );
    00838                 return tmp;
    00839         }
    00841
    00842 };
    00843
    00849 class logger : public root {
    00850 protected:
    00852         Array<RV> entries;
    00854         Array<string> names;
    00855 public:
    00857         logger() : entries ( 0 ), names ( 0 ) {}
    00858
    00861         virtual int add ( const RV &rv, string prefix = "" ) {
    00862                 int id;
    00863                 if ( rv._dsize() > 0 ) {
    00864                         id = entries.length();
    00865                         names = concat ( names, prefix ); // diff
    00866                         entries.set_length ( id + 1, true );
    00867                         entries ( id ) = rv;
    00868                 } else {
    00869                         id = -1;
    00870                 }
    00871                 return id; // identifier of the last entry
    00872         }
    00873
    00875         virtual void logit ( int id, const vec &v ) {
    00876                 bdm_error ( "Not implemented" );
    00877         };
    00879         virtual void logit ( int id, const double &d ) {
    00880                 bdm_error ( "Not implemented" );
    00881         };
    00882
    00884         virtual void step() {
    00885                 bdm_error ( "Not implemneted" );
    00886         };
    00887
    00889         virtual void finalize() {};
    00890
    00892         virtual void init() {};
    00893
    00894 };
    00895
    00899 class mepdf : public mpdf {
    00901         shared_ptr<epdf> iepdf;
    00902 public:
    00904         mepdf() { }
    00906         mepdf ( shared_ptr<epdf> em ) {
    00907                 iepdf = em;
    00908                 set_ep ( *iepdf.get() );
    00909                 dimc = 0;
    00910         }
    00911
    00913         vec samplecond ( const vec &cond ) {
    00914                 return iepdf->sample();
    00915         }
    00916         double evallogcond ( const vec &val, const vec &cond ) {
    00917                 return iepdf->evallog ( val );
    00918         }
    00919
    00927         void from_setting ( const Setting &set );
    00928 };
    00929 UIREGISTER ( mepdf );
    00930 SHAREDPTR ( mepdf );
    00931
    00933 RV get_composite_rv ( const Array<shared_ptr<mpdf> > &mpdfs, bool checkoverlap = false );
    00934
    00950 class DS : public root {
    00951 protected:
    00953         int dtsize;
    00955         int utsize;
    00957                 int ytsize;
    00959         RV Drv;
    00961         RV Urv; //
    00963                 RV Yrv; //
    00965         int L_dt, L_ut;
    00966 public:
    00968                 DS() : Drv(), Urv(),Yrv() {};
    00969
    00971         virtual int max_length() {return std::numeric_limits< int >::max();}
    00973         virtual void getdata ( vec &dt ) {
    00974                 bdm_error ( "abstract class" );
    00975         }
    00977         virtual void getdata ( vec &dt, const ivec &indeces ) {
    00978                 bdm_error ( "abstract class" );
    00979         }
    00980
    00982         virtual void write ( vec &ut ) {
    00983                 bdm_error ( "abstract class" );
    00984         }
    00985
    00987         virtual void write ( vec &ut, const ivec &indeces ) {
    00988                 bdm_error ( "abstract class" );
    00989         }
    00990
    00992         virtual void step()
    00993         {
    00994                 bdm_error ( "abstract class" );
    00995         }
    00996
    00998         virtual void log_add ( logger &L ) {
    00999                 bdm_assert ( dtsize == Drv._dsize(), "invalid DS: dtsize (" + num2str ( dtsize ) + ") different from Drv " + num2str ( Drv._dsize() ) );
    01000                 bdm_assert ( utsize == Urv._dsize(), "invalid DS: utsize (" + num2str ( utsize ) + ") different from Urv " + num2str ( Urv._dsize() ) );
    01001
    01002                 L_dt = L.add ( Drv, "" );
    01003                 L_ut = L.add ( Urv, "" );
    01004         }
    01006         virtual void logit ( logger &L ) {
    01007                 vec tmp ( Drv._dsize() + Urv._dsize() );
    01008                 getdata ( tmp );
    01009                 // d is first in getdata
    01010                 L.logit ( L_dt, tmp.left ( Drv._dsize() ) );
    01011                 // u follows after d in getdata
    01012                 L.logit ( L_ut, tmp.mid ( Drv._dsize(), Urv._dsize() ) );
    01013         }
    01015         virtual const RV& _drv() const {
    01016                 //              return concat (Drv, Urv);// why!!!
    01017                 return Drv;// why!!!
    01018         }
    01020         const RV& _urv() const {
    01021                 return Urv;
    01022         }
    01024                 const RV& _yrv() const {
    01025                         return Yrv;
    01026                 }
    01028                 virtual void set_drv (const  RV &yrv, const RV &urv) {
    01029                     Yrv = yrv;
    01030                         Drv = concat(yrv,urv);
    01031                 Urv = urv;
    01032         }
    01033 };
    01034
    01056 class BM : public root {
    01057 protected:
    01059         RV drv;
    01061         double ll;
    01063         bool evalll;
    01064 public:
    01067
    01068         BM() : ll ( 0 ), evalll ( true ), LIDs ( 4 ), LFlags ( 4 ) {
    01069                 LIDs = -1;/*empty IDs*/
    01070                 LFlags = 0;
    01071                 LFlags ( 0 ) = 1;  /*log only mean*/
    01072         };
    01073         BM ( const BM &B ) :  drv ( B.drv ), ll ( B.ll ), evalll ( B.evalll ) {}
    01076         virtual BM* _copy_() const {
    01077                 return NULL;
    01078         };
    01080
    01083
    01087         virtual void bayes ( const vec &dt ) = 0;
    01089         virtual void bayesB ( const mat &Dt );
    01092         virtual double logpred ( const vec &dt ) const {
    01093                 bdm_error ( "Not implemented" );
    01094                 return 0.0;
    01095         }
    01096
    01098         vec logpred_m ( const mat &dt ) const {
    01099                 vec tmp ( dt.cols() );
    01100                 for ( int i = 0; i < dt.cols(); i++ ) {
    01101                         tmp ( i ) = logpred ( dt.get_col ( i ) );
    01102                 }
    01103                 return tmp;
    01104         }
    01105
    01107         virtual epdf* epredictor() const {
    01108                 bdm_error ( "Not implemented" );
    01109                 return NULL;
    01110         };
    01112         virtual mpdf* predictor() const {
    01113                 bdm_error ( "Not implemented" );
    01114                 return NULL;
    01115         };
    01117
    01122
    01124         RV rvc;
    01126         const RV& _rvc() const {
    01127                 return rvc;
    01128         }
    01129
    01131         virtual void condition ( const vec &val ) {
    01132                 bdm_error ( "Not implemented!" );
    01133         }
    01134
    01136
    01137
    01140
    01141         const RV& _drv() const {
    01142                 return drv;
    01143         }
    01144         void set_drv ( const RV &rv ) {
    01145                 drv = rv;
    01146         }
    01147         void set_rv ( const RV &rv ) {
    01148                 const_cast<epdf&> ( posterior() ).set_rv ( rv );
    01149         }
    01150         double _ll() const {
    01151                 return ll;
    01152         }
    01153         void set_evalll ( bool evl0 ) {
    01154                 evalll = evl0;
    01155         }
    01156         virtual const epdf& posterior() const = 0;
    01158
    01161
    01163         virtual void set_options ( const string &opt ) {
    01164                 LFlags ( 0 ) = 1;
    01165                 if ( opt.find ( "logbounds" ) != string::npos ) {
    01166                         LFlags ( 1 ) = 1;
    01167                         LFlags ( 2 ) = 1;
    01168                 }
    01169                 if ( opt.find ( "logll" ) != string::npos ) {
    01170                         LFlags ( 3 ) = 1;
    01171                 }
    01172         }
    01174         ivec LIDs;
    01175
    01177         ivec LFlags;
    01179         virtual void log_add ( logger &L, const string &name = "" ) {
    01180                 // internal
    01181                 RV r;
    01182                 if ( posterior().isnamed() ) {
    01183                         r = posterior()._rv();
    01184                 } else {
    01185                         r = RV ( "est", posterior().dimension() );
    01186                 };
    01187
    01188                 // Add mean value
    01189                 if ( LFlags ( 0 ) ) LIDs ( 0 ) = L.add ( r, name + "mean_" );
    01190                 if ( LFlags ( 1 ) ) LIDs ( 1 ) = L.add ( r, name + "lb_" );
    01191                 if ( LFlags ( 2 ) ) LIDs ( 2 ) = L.add ( r, name + "ub_" );
    01192                 if ( LFlags ( 3 ) ) LIDs ( 3 ) = L.add ( RV ( "ll", 1 ), name );    //TODO: "local" RV
    01193         }
    01195         virtual void logit ( logger &L ) {
    01196                 L.logit ( LIDs ( 0 ), posterior().mean() );
    01197                 if ( LFlags ( 1 ) || LFlags ( 2 ) ) {  //if one of them is off, its LID==-1 and will not be stored
    01198                         vec ub, lb;
    01199                         posterior().qbounds ( lb, ub );
    01200                         L.logit ( LIDs ( 1 ), lb );
    01201                         L.logit ( LIDs ( 2 ), ub );
    01202                 }
    01203                 if ( LFlags ( 3 ) ) L.logit ( LIDs ( 3 ), ll );
    01204         }
    01206         void from_setting ( const Setting &set ) {
    01207                 shared_ptr<RV> r = UI::build<RV> ( set, "drv", UI::optional );
    01208                 if ( r ) {
    01209                         set_drv ( *r );
    01210                 }
    01211                 string opt;
    01212                 if ( !UI::get ( opt, set, "options", UI::optional ) ) {
    01213                         set_options ( opt );
    01214                 }
    01215         }
    01216
    01217 };
    01218
    01219 typedef Array<shared_ptr<epdf> > epdf_array;
    01220
    01221 typedef Array<shared_ptr<mpdf> > mpdf_array;
    01222
    01223 template<class EPDF>
    01224 vec mpdf_internal<EPDF>::samplecond ( const vec &cond ) {
    01225         condition ( cond );
    01226         vec temp = iepdf.sample();
    01227         return temp;
    01228 }
    01229
    01230 template<class EPDF>
    01231 mat mpdf_internal<EPDF>::samplecond_m ( const vec &cond, int N ) {
    01232         condition ( cond );
    01233         mat temp ( dimension(), N );
    01234         vec smp ( dimension() );
    01235         for ( int i = 0; i < N; i++ ) {
    01236                 smp = iepdf.sample();
    01237                 temp.set_col ( i, smp );
    01238         }
    01239
    01240         return temp;
    01241 }
    01242
    01243 template<class EPDF>
    01244 double mpdf_internal<EPDF>::evallogcond ( const vec &dt, const vec &cond ) {
    01245         double tmp;
    01246         condition ( cond );
    01247         tmp = iepdf.evallog ( dt );
    01248         return tmp;
    01249 }
    01250
    01251 template<class EPDF>
    01252 vec mpdf_internal<EPDF>::evallogcond_m ( const mat &Dt, const vec &cond ) {
    01253         condition ( cond );
    01254         return iepdf.evallog_m ( Dt );
    01255 }
    01256
    01257 template<class EPDF>
    01258 vec mpdf_internal<EPDF>::evallogcond_m ( const Array<vec> &Dt, const vec &cond ) {
    01259         condition ( cond );
    01260         return iepdf.evallog_m ( Dt );
    01261 }
    01262
    01263 }; //namespace
    01264 #endif // BDMBASE_H