Changeset 700

Show
Ignore:
Timestamp:
11/04/09 22:54:58 (15 years ago)
Author:
smidl
Message:

Making tutorial/userguide example work again (changes of mpdf and bayes)

Files:
2 removed
20 modified
1 moved

Legend:

Unmodified
Added
Removed

  • applications/bdmtoolbox/mex/estimator.cpp

    r685 r700  
    4848 */ 
    4949 
    50 #include "estim/arx.h" 
    51 #include "stat/emix.h" 
    52 #include "base/datasources.h" 
    53 #include "base/loggers.h" 
    54 #include "estim/particles.h" 
    55 #include "estim/kalman.h" 
     50#include <estim/arx.h> 
     51#include <stat/emix.h> 
     52#include <base/datasources.h> 
     53#include <base/loggers.h> 
     54#include <estim/particles.h> 
     55#include <estim/kalman.h> 
    5656 
    5757//#include "mex_datasource.h" 
    … …  
    6161#ifdef MEX 
    6262#include <itpp/itmex.h> 
    63 #include "mex/mex_BM.h" 
    64 #include "mex/mex_logger.h" 
    65 #include "mex/mex_datasource.h" 
     63#include <mex/mex_BM.h> 
     64#include <mex/mex_logger.h> 
     65#include <mex/mex_datasource.h> 
    6666 
    6767void mexFunction ( int n_output, mxArray *output[], int n_input, const mxArray *input[] ) { 
    … …  
    120120        Array<shared_ptr<BM> > Es; 
    121121        UI::get ( Es,Cfg, "estimators" ); 
    122         long Ndat=10; 
     122        int Ndat=10; 
    123123        if ( Cfg.exists ( "experiment" ) ) { 
    124                 if ( Cfg.lookupValue ( "experiment.ndat",Ndat ) ) { 
     124                if ( Cfg.getRoot().lookupValue ( "experiment.ndat",Ndat ) ) { 
    125125                        bdm_assert ( Ndat<=Ds->max_length(), "Data source has less data then required" ); 
    126126                }; 

  • applications/bdmtoolbox/mex/merger.cpp

    r685 r700  
    5858#endif 
    5959        // Sources 
    60         Array<shared_ptr<mpdf> >  Sources; 
     60        Array<shared_ptr<pdf> >  Sources; 
    6161        UI::get(Sources, Cfg, "Sources", UI::compulsory); 
    6262        shared_ptr<merger_base> Merger = UI::build<merger_base> ( Cfg, "Merger" ); 

  • applications/bdmtoolbox/mex/simulator.cpp

    r676 r700  
    3939 */ 
    4040 
    41 #include "estim/arx.h" 
    42 #include "stat/emix.h" 
    43 #include "base/datasources.h" 
    44 #include "base/loggers.h" 
     41#include <estim/arx.h> 
     42#include <stat/emix.h> 
     43#include <base/datasources.h> 
     44#include <base/loggers.h> 
    4545 
    4646//#include "mex_datasource.h" 
    … …  
    104104 
    105105        shared_ptr<DS> Ds = UI::build<DS> ( Cfg, "system" ); 
    106         long Ndat=10; 
     106        int Ndat=10; 
    107107        if ( Cfg.exists ( "experiment" ) ) { 
    108108                if ( Cfg.lookupValue ( "experiment.ndat",Ndat ) ) { 

  • applications/bdmtoolbox/tutorial/merging/pdfs.m

    r568 r700  
    3838    'rv', ab); 
    3939Gb_a=struct('class','mgamma',  'beta',2,  'k',1,  'rv',b,  'rvc',a); % f(b|a) 
    40 Ga_ =struct('class','mepdf','epdf',Ga); % convert f(a) to f(a|) 
    41 Gba = struct('class','mprod',  'mpdfs',{{Gb_a,Ga_}}); 
     40Gba = struct('class','mprod',  'pdfs',{{Gb_a,Ga}}); 
    4241 
    4342pd.Ga=Ga; 

  • applications/bdmtoolbox/tutorial/userguide/arx_basic_example.m

    r631 r700  
    11% load data created by the MpdfDS_example 
    2 load mpdfds_results 
     2load pdfds_results 
    33 
    44DS.class   = 'MemDS'; 
    … …  
    1717%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
    1818% plot results 
    19 ndat = size(M.u,1); 
     19ndat = size(M.DS_u,1); 
    2020 
    2121subplot(1,2,1); 
    2222hold off 
    23 plot((true_theta' *ones(1,ndat))','-.'); 
    24 title(' Regression parameters \theta') 
    25 hold on 
    26 plot(M.mean_theta(1:ndat,:)); 
    27 co = get(gca,'ColorOrder'); 
    28 for i=1:size(M.mean_theta,2) 
    29     ind =1:10:ndat; 
    30     h=errorbar(ind,M.mean_theta(ind,i),... 
    31     M.mean_theta(ind,i)-M.lb_theta(ind,i),M.mean_theta(ind,i)-M.ub_theta(ind,i),'.'); 
    32     set(h,'color',co(i,:)); 
    33 end 
    34  
     23plotestimates(true_theta, ... 
     24    M.Est0_apost_mean_theta, ... 
     25    M.Est0_apost_lb_theta, ... 
     26    M.Est0_apost_ub_theta); 
    3527set(gca,'YLim',[-1.5,1]); 
    3628 
    3729subplot(1,2,2); 
    3830hold off 
    39 plot(true_R*ones(1,ndat),'-.'); 
     31plotestimates(true_R, ... 
     32    M.Est0_apost_mean_r, ... 
     33    M.Est0_apost_lb_r, ... 
     34    M.Est0_apost_ub_r); 
     35 
    4036title('Variance parameters r') 
    41 hold on 
    42 plot(M.mean_r(1:ndat,:)); 
    43 co = get(gca,'ColorOrder'); 
    44 ind =1:10:ndat; 
    45 for i=1:size(M.mean_r, 2) 
    46     h=errorbar(ind,M.mean_r(ind,i),... 
    47     M.mean_r(ind,i)-M.lb_r(ind,i),M.mean_r(ind,i)-M.ub_r(ind,i),'.'); 
    48     set(h,'color',co(i,:)); 
    49 end 

  • applications/bdmtoolbox/tutorial/userguide/arx_selection_example.m

    r631 r700  
    11% load data created by the MpdfDS_example 
    2 load mpdfds_results 
     2load pdfds_results 
    33 
    44DS.class   = 'MemDS'; 
    … …  
    3030 
    3131%%%% Process results 
    32 lls = [sum(M.A1ll) sum(M.A2ll) sum(M.A3ll)] 
     32lls = [sum(M.A1_ll_ll) sum(M.A2_ll_ll) sum(M.A3_ll_ll)] 
    3333 
    3434ells=exp(lls-max(lls)); 

  • applications/bdmtoolbox/tutorial/userguide/epdfds_example.m

    r618 r700  
    1212M=simulator(DS,experiment); 
    1313 
    14 M.a 
     14M.DS_a 

  • applications/bdmtoolbox/tutorial/userguide/frg_estim.m

    r661 r700  
    1717%%%%%% Random walk on frg - Dirichlet  
    1818phi_pdf.class = 'mDirich';         % random walk on coefficient phi 
    19 phi_pdf.rv    = RV('phi',2);       % 2D random walk - frg is the first element 
     19phi_pdf.rv    = RV({'phi','1_phi'});       % 2D random walk - frg is the first element 
    2020phi_pdf.k     = 0.01;              % width of the random walk 
    2121phi_pdf.betac = [0.01 0.01];       % stabilizing elememnt of random walk 
    … …  
    2626E.parameter_pdf = phi_pdf;         % Random walk is the parameter evolution model 
    2727E.res_threshold = 1.0;             % resampling parameter 
    28 E.n = 20;                          % number of particles 
     28E.n = 10;                          % number of particles 
    2929E.prior.class = 'eDirich';         % prior on non-linear part 
    3030E.prior.beta  = [2 1]; %  
    … …  
    3636%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
    3737% plot results 
    38 ndat = size(M.u,1); 
     38ndat = size(M.DS_u,1); 
    3939 
    4040figure(1); 
    4141subplot(2,2,1); 
    42 plotestimates(true_theta, M.MPFmean_theta, M.MPFlb_theta, M.MPFub_theta); 
     42plotestimates(true_theta, M.MPF_apost_mean_theta, M.MPF_apost_lb_theta, M.MPF_apost_ub_theta); 
    4343title(' Regression parameters \theta') 
    4444set(gca,'YLim',[-1.5,1]); 
    4545 
    4646subplot(2,2,2); 
    47 plotestimates(true_R, M.MPFmean_r,M.MPFlb_r,M.MPFub_r); 
     47plotestimates(true_R, M.MPF_apost_mean_r,M.MPF_apost_lb_r,M.MPF_apost_ub_r); 
    4848title('Variance parameters r') 
    4949 
    5050subplot(2,2,3); 
    51 plotestimates(1, M.MPFmean_phi(:,1),M.MPFlb_phi(:,1),M.MPFub_phi(:,1)); 
     51plotestimates(1, M.MPF_apost_mean_phi(:,1),M.MPF_apost_lb_phi(:,1),M.MPF_apost_ub_phi(:,1)); 
    5252title('Forgetting factor') 
    5353 

  • applications/bdmtoolbox/tutorial/userguide/pdfds_example.m

    r631 r700  
    2020DS.class = 'MpdfDS'; 
    2121DS.mpdf.class  = 'mprod'; 
    22 DS.mpdf.mpdfs  = {fy, epdf2mpdf(fu)}; 
     22DS.mpdf.mpdfs  = {fy, fu}; 
    2323DS.init_rv = RVtimes([y,y,y], [-1,-2,-3]); 
    2424DS.init_values = [0.1, 0.2, 0.3]; 
    … …  
    3030 
    3131%%% store results 
    32 Data=[M.y'; M.u']; 
     32Data=[M.DS_y'; M.DS_u']; 
    3333drv = RVjoin([y,u]); 
    3434true_theta=[fy.A fy.const];  

  • applications/pmsm/filters.h

    r686 r700  
    1919 
    2020//! modelling theta as normal random walk with 2pi correction 
    21         class rwtheta: public mpdf_internal<enorm<fsqmat> > { 
     21        class rwtheta: public pdf_internal<enorm<fsqmat> > { 
    2222                double om_hat_dt; 
    2323                double om_var; 
    … …  
    120120                pf->prior_from_set ( set ); 
    121121                pf->resmethod_from_set ( set ); 
    122                 pf->set_model ( rwt,new mpdf() ); 
     122                pf->set_model ( rwt,new pdf() ); 
    123123 
    124124                shared_ptr<RV> rv = UI::build<RV> ( set, "rv", UI::optional ); 

  • library/bdm/base/bdmbase.cpp

    r693 r700  
    5353                id = iter->second; 
    5454                if (size>0){ 
    55                         bdm_assert ( SIZES ( id ) == size, "RV " + name + " of different size already exists" ); 
     55                        bdm_assert ( SIZES ( id ) == size, "RV " + name + " of size " + num2str(SIZES(id)) + " exists, requested size " + num2str(size) + "can not be assigned" ); 
    5656                } 
    5757        } 

  • library/bdm/base/bdmbase.h

    r693 r700  
    266266        //! Minimum time-offset 
    267267        int mint() const { 
    268                         return times.length()>0 ? min (times) : 0; 
     268                return times.length()>0 ? min (times) : 0; 
     269        } 
     270        //! Minimum time-offset of ids of given RVs 
     271        int mint(const RV &rv) const { 
     272                bvec belong=zeros_b(len); 
     273                for (int r=0; r<rv.length(); r++){ 
     274                        belong = belong | (ids == rv.id(r)); 
     275                } 
     276                return times.length()>0 ? min (times(belong)) : 0; 
    269277        } 
    270278        //!@} 
    … …  
    598606                        logrec->L.logit( logrec->ids(0), mean() ); 
    599607                } 
    600                 if (log_level>2) { 
     608                if (log_level>1) { 
    601609                        vec lb; vec ub; 
    602610                        qbounds(lb,ub); 
    … …  
    913921                //establish c2c connection 
    914922                rvc.dataind ( rvc_up, c2c_lo, c2c_up ); 
    915                 bdm_assert_debug ( c2c_lo.length() + v2c_lo.length() == condsize, "cond is not fully given" ); 
     923//              bdm_assert_debug ( c2c_lo.length() + v2c_lo.length() == condsize, "cond is not fully given" ); 
    916924        } 
    917925 
    … …  
    919927        vec get_cond ( const vec &val_up, const vec &cond_up ) { 
    920928                vec tmp ( condsize ); 
    921                 set_subvector ( tmp, v2c_lo, val_up ( v2c_up ) ); 
    922                 set_subvector ( tmp, c2c_lo, cond_up ( c2c_up ) ); 
     929                fill_cond (val_up, cond_up, tmp ); 
    923930                return tmp; 
     931        } 
     932        //! fill condition  
     933        void fill_cond ( const vec &val_up, const vec &cond_up, vec& cond_out){ 
     934                bdm_assert_debug(cond_out.length()>=condsize,"dl.fill_cond: cond_out is too small"); 
     935                set_subvector ( cond_out, v2c_lo, val_up ( v2c_up ) ); 
     936                set_subvector ( cond_out, c2c_lo, cond_up ( c2c_up ) ); 
    924937        } 
    925938        //! Fill 
    … …  
    10771090 
    10781091        BM() : yrv(),dimy(0),rvc(),dimc(0), ll ( 0 ), evalll ( true ) { }; 
    1079         BM ( const BM &B ) :  yrv ( B.yrv ), dimy(B.dimy), rvc ( B.rvc ), ll ( B.ll ), evalll ( B.evalll ) {} 
     1092//      BM ( const BM &B ) :  yrv ( B.yrv ), dimy(B.dimy), rvc ( B.rvc ),dimc(B.dimc), ll ( B.ll ), evalll ( B.evalll ) {} 
    10801093        //! \brief Copy function required in vectors, Arrays of BM etc. Have to be DELETED manually! 
    10811094        //! Prototype: \code BM* _copy_() const {return new BM(*this);} \endcode 

  • library/bdm/base/datasources.h

    r697 r700  
    132132                        dt =  zeros(iepdf->dimension()); 
    133133                        dtsize=dt.length(); 
     134                        ytsize=dt.length(); 
    134135                        set_drv(iepdf->_rv(),RV()); 
    135136                        utsize =0; 
    … …  
    365366                class = "stateDS"; 
    366367                //Internal model 
    367                 IM = { type = "mpdf-offspring"; }; 
     368                IM = { type = "pdf-offspring"; }; 
    368369                //Observation model 
    369                 OM = { type = "mpdf-offspring"; } 
     370                OM = { type = "pdf-offspring"; } 
    370371                //initial state 
    371372                x0 = [...]; //vector of data 

  • library/bdm/estim/arx.cpp

    r679 r700  
    22namespace bdm { 
    33 
    4         void ARX::bayes_weighted ( const vec &yt, const vec &cond, const double w ) { 
     4void ARX::bayes_weighted ( const vec &yt, const vec &cond, const double w ) { 
     5                 
     6        bdm_assert_debug(yt.length()>=dimy,"ARX::bayes yt is smaller then dimc");        
     7        bdm_assert_debug(cond.length()>=dimc,"ARX::bayes cond is smaller then dimc");    
    58        double lnc; 
    69        //cache 

  • library/bdm/estim/arx.h

    r679 r700  
    5454        //! \name Constructors 
    5555        //!@{ 
    56         ARX ( const double frg0 = 1.0 ) : BMEF ( frg0 ),  have_constant(true), dyad(), est() {}; 
    57         ARX ( const ARX &A0 ) : BMEF (A0.frg),  have_constant(A0.have_constant), dyad(A0.dyad),est(est) { }; 
     56        ARX ( const double frg0 = 1.0 ) : BMEF ( frg0 ),  have_constant(true), dyad(), est(), alter_est() {}; 
     57        ARX ( const ARX &A0 ) : BMEF (A0),  have_constant(A0.have_constant), dyad(A0.dyad),est(A0.est),alter_est(A0.alter_est) { }; 
    5858        ARX* _copy_() const; 
    5959        void set_parameters ( double frg0 ) { 
    … …  
    167167/*! ARX model conditined by knowledge of the forgetting factor 
    168168\f[ f(\theta| d_1 \ldots d_t , \phi_t) \f] 
     169 
     170The symbol \f$ \phi \f$ is assumed to be the last of the conditioning variables. 
    169171*/ 
    170172class ARXfrg : public ARX{ 
    … …  
    174176                ARXfrg(const ARXfrg &A0):ARX(A0){}; 
    175177                ARXfrg* _copy_() const {ARXfrg *A = new ARXfrg(*this); return A;} 
    176         void condition(const vec &val){ 
    177                 frg = val(0); 
     178 
     179        void bayes(const vec &val, const vec &cond){ 
     180                frg = cond(dimc-1); //  last in cond is phi 
     181                ARX::bayes(val,cond); 
     182        } 
     183        void validate() { 
     184                ARX::validate(); 
     185                rvc.add(RV("{phi }",vec_1(1))); 
     186                dimc +=1; 
    178187        } 
    179188}; 

  • library/bdm/estim/particles.cpp

    r679 r700  
    7878        int n=pf->__w().length(); 
    7979        vec &lls = pf->_lls(); 
     80        Array<vec> &samples=pf->__samples(); 
    8081         
    8182        // generate samples - time step 
    … …  
    8485        #pragma parallel for 
    8586        for ( i = 0; i < n; i++ ) { 
    86                 BMs(i) -> bayes(this2bm.pushdown(yt), this2bm.get_cond(yt,cond)); 
     87                vec bm_cond(BMs(i)->dimensionc()); 
     88                this2bm.fill_cond(yt,cond, bm_cond); 
     89                pf2bm.filldown(samples(i), bm_cond); 
     90                BMs(i) -> bayes(this2bm.pushdown(yt), bm_cond); 
    8791                lls ( i ) += BMs(i)->_ll(); 
    8892        } 

  • library/bdm/estim/particles.h

    r693 r700  
    157157                        res_threshold=0.5; 
    158158                } 
    159                 validate(); 
     159                //validate(); 
    160160        } 
    161161        //! load prior information from set and set internal structures accordingly 
    … …  
    177177         
    178178        void validate(){ 
     179                bdm_assert(par,"PF::parameter_pdf not specified."); 
    179180                n=_w.length(); 
    180181                lls=zeros(n); 
     182                 
    181183                if (par->_rv()._dsize()>0) { 
    182184                        bdm_assert(par->_rv()._dsize()==est.dimension(),"Mismatch of RV and dimension of posterior" ); 
    … …  
    301303        //! datalink from global yt and cond (Up) to PFs yt and cond (Down) 
    302304        datalink_m2m this2pf; 
     305        //!datalink from PF part to BM 
     306        datalink_part pf2bm; 
    303307         
    304308public: 
    … …  
    367371                } 
    368372                //set drv 
     373                //??set_yrv(concat(BM0->_yrv(),u) ); 
     374                set_yrv(BM0->_yrv() ); 
     375                rvc = BM0->_rvc().subt(par->_rv()); 
    369376                //find potential input - what remains in rvc when we subtract rv 
    370                 RV u = par->_rvc().remove_time().subt( par->_rv() );             
    371                 set_yrv(concat(BM0->_yrv(),u) ); 
     377                RV u = par->_rvc().subt( par->_rv().copy_t(-1) );                
     378                rvc.add(u); 
     379                dimc = rvc._dsize(); 
    372380                validate(); 
    373381        } 
    374382        void validate(){ 
    375383                try{ 
    376                 pf->validate(); 
     384                        pf->validate(); 
    377385                } catch (std::exception &e){ 
    378386                        throw UIException("Error in PF part of MPF:"); 
    … …  
    380388                jest.read_parameters(); 
    381389                this2bm.set_connection(BMs(0)->_yrv(), BMs(0)->_rvc(), yrv, rvc); 
    382                 this2bm.set_connection(pf->_yrv(), pf->_rvc(), yrv, rvc); 
     390                this2pf.set_connection(pf->_yrv(), pf->_rvc(), yrv, rvc); 
     391                pf2bm.set_connection(BMs(0)->_rvc(), pf->posterior()._rv()); 
    383392        } 
    384393 

  • library/bdm/math/square_mat.h

    r619 r700  
    221221        //! print full matrix 
    222222        friend std::ostream &operator<< ( std::ostream &os, const fsqmat &sq ); 
     223        //!access function 
     224        mat & _M (  ) {return M;}; 
    223225 
    224226}; 

  • library/bdm/mex/config2mxstruct.h

    r665 r700  
    22#define CFGSTR_H 
    33 
    4 #include "../base/libconfig/libconfig.h++" 
     4#include "../base/libconfig/lib/libconfig.h++" 
    55#include <itpp/itbase.h> 
    66#include <itpp/itmex.h> 

  • library/bdm/mex/mex_parser.h

    r691 r700  
    33 
    44 
    5 #include "../base/libconfig/lib/lib/libconfig.h++" 
     5#include "base/libconfig/lib/libconfig.h++" 
    66#include <itpp/itbase.h> 
    77#include <itpp/itmex.h> 

  • library/tests/CMakeLists.txt

    r693 r700  
    4141# using UnitTest++ 
    4242 
    43 SET(unit_test_configurations  mepdf.cfg egiw.cfg mlnorm.cfg edirich.cfg mprod.cfg generator.cfg epdfds.cfg pdfds.cfg test_user_info_input.cfg merger.cfg erroneous.cfg egamma.cfg test_user_info_external.cfg merger_error.cfg mgamma.cfg enorm.cfg) 
     43SET(unit_test_configurations  egiw.cfg mlnorm.cfg edirich.cfg mprod.cfg generator.cfg epdfds.cfg pdfds.cfg test_user_info_input.cfg merger.cfg erroneous.cfg egamma.cfg test_user_info_external.cfg merger_error.cfg mgamma.cfg enorm.cfg) 
    4444 
    4545add_executable(testsuite arx_straux_test.cpp datalink_test.cpp  datasource_test.cpp egiw_test.cpp emix_test.cpp epdf_test.cpp logger_test.cpp LQG_test.cpp merger_test.cpp