Changeset 192

Timestamp:

10/22/08 10:46:40 (16 years ago)

Author:

smidl

Message:

modification of datalinks and switch mprod and merger to use them

Files:

• bdm/estim/merger.cpp (modified) (7 diffs)
• bdm/estim/merger.h (modified) (2 diffs)
• bdm/stat/emix.cpp (modified) (1 diff)
• bdm/stat/emix.h (modified) (4 diffs)
• bdm/stat/libBM.cpp (modified) (1 diff)
• bdm/stat/libBM.h (modified) (5 diffs)
• bdm/stat/libEF.h (modified) (2 diffs)
• tests/datalink_test.cpp (modified) (1 diff)

bdm/estim/merger.cpp

@@ -33 +33 @@
         vec &w = eSmp._w(); //aux
 
-        mat Smp_ex =ones ( rv.count() +1,Ns ); // Extended samples for the ARX model - the last row is ones 
+        mat Smp_ex =ones ( rv.count() +1,Ns ); // Extended samples for the ARX model - the last row is ones
         for ( int i=0;i<Ns;i++ ) {      set_col_part ( Smp_ex,i,Smp ( i ) );}
 
@@ -60 +60 @@
                 //Re-Initialize Mixture model
                 Mix.init ( &A0, Smp_ex, Nc );
-                Mix.bayesB ( Smp_ex, ones(Ns));//w*Ns );
-                Mpred = Mix.predictor(rv); // Allocation => must be deleted at the end!!
-        
+                Mix.bayesB ( Smp_ex, ones ( Ns ) );//w*Ns );
+                Mpred = Mix.predictor ( rv ); // Allocation => must be deleted at the end!!
+
                 // Generate new samples
                 eSmp.set_samples ( Mpred );
@@ -79 +79 @@
                         //======== Same RVs ===========
                         //Split according to dependency in rvs
-                        if ( rvsinrv ( i ).length() ==rv.count() ) {
+                        if ( mpdfs ( i )->_rv().count() ==rv.count() ) {
                                 // no need for conditioning or marginalization
                                 for ( int j=0;j<Ns; j++ ) { // Smp is Array<> => for cycle
@@ -86 +86 @@
                         }
                         else {
-                                vec smpk;
                                 // compute likelihood of marginal on the conditional variable
                                 if ( mpdfs ( i )->_rvc().count() >0 ) {
                                         // Make marginal on rvc_i
                                         epdf* tmp_marg = Mpred->marginal ( mpdfs ( i )->_rvc() );
-                                        for(int k=0;k<Ns;k++){
-                                        lw_src(k) += tmp_marg->evalpdflog ( get_vec(Smp(i), irv_rvcs(i) ) );}
+                                        //compute vector of lw_src
+                                        for ( int k=0;k<Ns;k++ ) {
+                                                lw_src ( k ) += tmp_marg->evalpdflog ( dls ( i )->get_val ( Smp ( i ) ) );
+                                        }
                                         delete tmp_marg;
                                 }
                                 // Compute likelihood of the missing variable
-                                if ( rv.count() > mpdfs ( i )->_rv().count() + mpdfs ( i )->_rvc().count() ) {
-                                        // There are variales unknown to mpdfs
-                                        RV z = ( rv.subt ( mpdfs ( i )->_rv() ) ).subt ( mpdfs ( i )->_rvc() );
-                                        mpdf* tmp_cond = Mpred->condition ( z ); 
-                                        // Indeces of rest rv in Smp
-                                        ivec zinrv=z.dataind ( rv ) ;
-                                        // Indeces of rest rvc in Smp
-                                        ivec zinrvc=tmp_cond->_rvc().dataind ( rv );
+                                if ( rv.count() > ( mpdfs ( i )->_rv().count() + mpdfs ( i )->_rvc().count() ) ) {
+                                        // There are variales unknown to mpdfs(i) : rvzs
+                                        mpdf* tmp_cond = Mpred->condition ( rvzs ( i ) );
                                         // Compute likelihood
                                         for ( int k= 0; k<Ns; k++ ) {
-                                                smpk=Smp( k );
-                                                lw_src ( k ) += log( tmp_cond->evalcond ( get_vec ( smpk,zinrv ), get_vec ( smpk,zinrvc ) ));
+                                                lw_src ( k ) += log (
+                                                                    tmp_cond->evalcond (
+                                                                        zdls ( i )->get_val ( Smp ( k ) ),
+                                                                        zdls ( i )->get_cond ( Smp ( k ) ) ) );
                                         }
                                         delete tmp_cond;
@@ -113 +111 @@
                                 // Compute likelihood of the partial source
                                 for ( int k= 0; k<Ns; k++ ) {
-                                        smpk=Smp( k );
-                                        mpdfs ( i )->condition ( get_vec ( smpk,irv_rvcs(i) ) );
-                                        lw_src ( k ) += mpdfs ( i )->_epdf().evalpdflog ( get_vec ( smpk, rvsinrv ( i ) ) );
+                                        mpdfs ( i )->condition ( dls ( i )->get_cond ( Smp(k) ) );
+                                        lw_src ( k ) += mpdfs ( i )->_epdf().evalpdflog ( dls ( i )->get_val ( Smp(k) ) );
                                 }
                         }
@@ -128 +125 @@
 
                 w = lognorm_merge ( lW ); //merge
-                
+
                 sprintf ( str,"w%d",niter );
                 dbg << Name ( str ) << w;
@@ -150 +147 @@
                 // ==== stopping rule ===
                 niter++;
-                converged = ( niter>20);
+                converged = ( niter>20 );
         }
 

bdm/estim/merger.h

@@ -30 +30 @@
         //!Internal mixture of EF models
         MixEF Mix;
+        //! Data link for each mpdf in mpdfs
+        Array<datalink_m2e*> dls;
+        //! Array of rvs that are not modelled by mpdfs at all (aux)
+        Array<RV> rvzs;
+        //! Data Links of rv0 mpdfs - these will be conditioned the (rv,rvc) of mpdfs
+        Array<datalink_m2e*> zdls; 
+        
         //!Number of samples used in approximation
         int Ns;
@@ -40 +47 @@
         merger ( const Array<mpdf*> &S ) :
                         compositepdf ( S ), epdf ( getrv ( false ) ),
-                        Mix ( Array<BMEF*> ( 0 ),vec ( 0 ) )
-        {setindices(rv); beta=2.0; Ns=100; Nc=10; Mix.set_method(EM);}
-        //! Set internal parameters used in approximation
+                        Mix ( Array<BMEF*> ( 0 ),vec ( 0 ) ), dls(n), rvzs(n), zdls(n) {
+                                RV ztmp;
+                for ( int i=0;i<n;i++ ) {
+                        //Establich connection between mpdfs and merger
+                        dls ( i ) = new datalink_m2e ( mpdfs ( i )->_rv(), mpdfs(i)->_rvc(), rv );
+                        // find out what is missing in each mpdf
+                        ztmp= mpdfs(i)->_rv();
+                        ztmp.add(mpdfs(i)->_rvc());
+                        rvzs(i)=rv.subt(ztmp);
+                        zdls(i) = new datalink_m2e (rvzs(i), ztmp, rv ) ;
+                };
+                //Set Default values of parameters
+                beta=2.0;
+                Ns=100;
+                Nc=10;
+                Mix.set_method ( EM );
+        }
+//! Set internal parameters used in approximation
         void set_parameters ( double beta0, int Ns0, int Nc0 ) {        beta=beta0;Ns=Ns0;Nc=Nc0;}
-        //!Initialize the proposal density. This function must be called before merge()!
-        void init() {
+//!Initialize the proposal density. This function must be called before merge()!
+        void init() { ////////////// NOT FINISHED
                 Array<vec> Smps ( n );
                 //Gibbs sampling
                 for ( int i=0;i<n;i++ ) {Smps ( i ) =zeros ( 0 );}
         }
-        //!Create a mixture density using known proposal
+//!Create a mixture density using known proposal
         void merge ( const epdf* g0 );
-        //!Create a mixture density, make sure to call init() before the first call
+//!Create a mixture density, make sure to call init() before the first call
         void merge () {merge ( & ( Mix._epdf() ) );};
 
-        //! Merge log-likelihood values
+//! Merge log-likelihood values
         vec lognorm_merge ( mat &lW );
-        //! sample from merged density
-        //! weight w is a
-        vec sample ( )const { return Mix._epdf().sample();}
-        double evalpdflog ( const vec &dt ) const{ 
-                vec dtf=ones(dt.length()+1);    
-                dtf.set_subvector(0,dt); 
-                return Mix.logpred ( dtf );}
-        vec mean()const {return Mix._epdf().mean();}
-        //! for future use
-        virtual ~merger() {};
-        
-        //! Access function
+//! sample from merged density
+//! weight w is a
+        vec sample ( ) const { return Mix._epdf().sample();}
+        double evalpdflog ( const vec &dt ) const {
+                vec dtf=ones ( dt.length() +1 );
+                dtf.set_subvector ( 0,dt );
+                return Mix.logpred ( dtf );
+        }
+        vec mean() const {return Mix._epdf().mean();}
+//! for future use
+        virtual ~merger() {
+                for (int i=0; i<n; i++){
+                        delete dls(i);
+                        delete zdls(i);                 
+                }
+        };
+
+//! Access function
         MixEF& _Mix() {return Mix;}
 };

bdm/stat/emix.cpp

@@ -11 +11 @@
         if ( copy ) {
                 Coms.set_length(Coms0.length());
-                for ( i=0;i<w.length();i++ ) {*Coms ( i ) =*Coms0 ( i );}
+                for ( i=0;i<w.length();i++ ) {it_error("Not imp...");
+                        *Coms ( i ) =*Coms0 ( i );}
                 destroyComs=true;
         }

bdm/stat/emix.h

@@ -35 +35 @@
  */
 class mratio: public mpdf {
-        protected:
+protected:
         //! Nominator in the form of mpdf
-                const epdf* nom;
+        const epdf* nom;
         //!Denominator in the form of epdf
-                epdf* den;
-        //!flag for destructor 
-                bool destroynom;
-        public:
-        //!Default constructor. By default, the given epdf is not copied! 
+        epdf* den;
+        //!flag for destructor
+        bool destroynom;
+public:
+        //!Default constructor. By default, the given epdf is not copied!
         //! It is assumed that this function will be used only temporarily.
-                mratio(const epdf* nom0, const RV &rv, bool copy=false):mpdf(rv,RV()){
-                        if (copy){
+        mratio ( const epdf* nom0, const RV &rv, bool copy=false ) :mpdf ( rv,RV() ) {
+                if ( copy ) {
 //                      nom = nom0->_copy_();
-                                destroynom=true;
-                        } else {
-                                nom = nom0;
-                                destroynom = false;
-                        }
-                        rvc = nom->_rv().subt(rv);
-                        it_assert_debug(rvc.length()>0,"Makes no sense to use this object!");
-                        den = nom->marginal(rvc);
-                };
-                double evalcond(const vec &val, const vec &cond) {
-                        return exp(nom->evalpdflog(concat(val,cond)) - den->evalpdflog(cond));
-                }
+                        destroynom=true;
+                }
+                else {
+                        nom = nom0;
+                        destroynom = false;
+                }
+                rvc = nom->_rv().subt ( rv );
+                it_assert_debug ( rvc.length() >0,"Makes no sense to use this object!" );
+                den = nom->marginal ( rvc );
+        };
+        double evalcond ( const vec &val, const vec &cond ) {
+                return exp ( nom->evalpdflog ( concat ( val,cond ) ) - den->evalpdflog ( cond ) );
+        }
         //! Object takes ownership of nom and will destroy it
-                void ownnom(){destroynom=true;}
+        void ownnom() {destroynom=true;}
         //! Default destructor
-                ~mratio(){delete den; if (destroynom) {delete nom;}}
+        ~mratio() {delete den; if ( destroynom ) {delete nom;}}
 };
 
@@ -100 +101 @@
                 int i;
                 double sum = 0.0;
-                for ( i = 0;i < w.length();i++ ) {sum += w ( i ) * exp(Coms ( i )->evalpdflog ( val ));}
+                for ( i = 0;i < w.length();i++ ) {sum += w ( i ) * exp ( Coms ( i )->evalpdflog ( val ) );}
                 return log ( sum );
         };
         vec evalpdflog_m ( const mat &Val ) const {
-                vec x=ones(Val.cols());
-                vec y(Val.cols());
-                for (int i = 0; i < w.length(); i++ ) {
-                        y = w ( i )*exp(Coms ( i )->evalpdflog_m ( Val ) );
-                        elem_mult_inplace(y,x); //result is in x
-                }
-                return log(x);
+                vec x=ones ( Val.cols() );
+                vec y ( Val.cols() );
+                for ( int i = 0; i < w.length(); i++ ) {
+                        y = w ( i ) *exp ( Coms ( i )->evalpdflog_m ( Val ) );
+                        elem_mult_inplace ( y,x ); //result is in x
+                }
+                return log ( x );
         };
         mat evalpdflog_M ( const mat &Val ) const {
-                mat X(w.length(), Val.cols());
-                for (int i = 0; i < w.length(); i++ ) {
-                        X.set_row(i, w ( i )*exp(Coms ( i )->evalpdflog_m ( Val ) ));
+                mat X ( w.length(), Val.cols() );
+                for ( int i = 0; i < w.length(); i++ ) {
+                        X.set_row ( i, w ( i ) *exp ( Coms ( i )->evalpdflog_m ( Val ) ) );
                 }
                 return X;
@@ -143 +144 @@
         //! pointers to epdfs - shortcut to mpdfs()._epdf()
         Array<epdf*> epdfs;
-        //! Indeces of rvc in common rvc
-        Array<ivec> irvcs_rvc;
-        //! Indeces of common rvc in rvcs
-        Array<ivec> irvc_rvcs;
+        //! Data link for each mpdfs
+        Array<datalink_m2m*> dls;
 public:
         /*!\brief Constructor from list of mFacs,
         */
-        mprod ( Array<mpdf*> mFacs ) : compositepdf ( mFacs ), mpdf ( getrv ( true ),RV() ), epdfs ( n ), irvcs_rvc ( n ), irvc_rvcs ( n ) {
+        mprod ( Array<mpdf*> mFacs ) : compositepdf ( mFacs ), mpdf ( getrv ( true ),RV() ), epdfs ( n ), dls ( n ) {
                 setrvc ( rv,rvc );
-                setindices ( rv );
+                // rv and rvc established = > we can link them with mpdfs
                 for ( int i = 0;i < n;i++ ) {
-                        // establish link between rvc and rvcs
-                        rvc.dataind ( mpdfs ( i )->_rvc(), irvc_rvcs ( i ), irvcs_rvc ( i ) );
+                        dls ( i ) = new datalink_m2m ( mpdfs ( i )->_rv(), mpdfs(i)->_rvc(), rv, rvc );
                 }
 
@@ -166 +164 @@
                 int i;
                 double res = 0.0;
-                vec condi;
                 for ( i = n - 1;i >= 0;i-- ) {
-                        if ( ( irvcs_rvc ( i ).length() > 0 ) || ( irvcs_rv ( i ).length() > 0 ) ) {
-                                condi = zeros ( irvcs_rvc ( i ).length() +irvcs_rv ( i ).length() );
-                                //copy part of the rvc into condi
-                                set_subvector ( condi, irvcs_rvc ( i ), get_vec ( cond,irvc_rvcs ( i ) ) );
-                                //copy part of the rv into condi
-                                set_subvector ( condi, irvcs_rv ( i ), get_vec ( val,irv_rvcs ( i ) ) );
-                                mpdfs ( i )->condition ( condi );
+                        if ( mpdfs(i)->_rvc().count() >0) {
+                                mpdfs ( i )->condition ( dls ( i )->get_cond ( val,cond ) );
                         }
                         // add logarithms
-                        res += epdfs ( i )->evalpdflog ( val ( rvsinrv ( i ) ) );
+                        res += epdfs ( i )->evalpdflog ( dls ( i )->get_val ( val ) );
                 }
                 return exp ( res );
         }
         vec samplecond ( const vec &cond, double &ll ) {
-                vec smp=zeros ( rv.count() );
-                vec condi;
+                //! Ugly hack to help to discover if mpfs are not in proper order. Correct solution = check that explicitely.
+                vec smp= std::numeric_limits<double>::infinity() * ones ( rv.count() );
                 vec smpi;
                 ll = 0;
+                // Hard assumption here!!! We are going backwards, to assure that samples that are needed from smp are already generated!
                 for ( int i = ( n - 1 );i >= 0;i-- ) {
-                        if ( ( irvcs_rvc ( i ).length() > 0 ) || ( irvcs_rv ( i ).length() > 0 ) ) {
-                                condi = zeros ( irvcs_rv ( i ).length() + irvcs_rvc ( i ).length() );
-                                // copy data in condition
-                                set_subvector ( condi,irvcs_rvc ( i ), get_vec ( cond, irvc_rvcs ( i ) ) );
-                                // copy data in already generated sample
-                                set_subvector ( condi,irvcs_rv ( i ), get_vec ( smp,irv_rvcs ( i ) ) );
-
-                                mpdfs ( i )->condition ( condi );
+                        if ( mpdfs(i)->_rvc().count() ) {
+                                mpdfs ( i )->condition ( dls ( i )->get_cond ( smp ,cond ) ); // smp is val here!!
                         }
                         smpi = epdfs ( i )->sample();
                         // copy contribution of this pdf into smp
-                        set_subvector ( smp,rvsinrv ( i ), smpi );
+                        dls(i)->fill_val(smp, smpi);
                         // add ith likelihood contribution
                         ll+=epdfs ( i )->evalpdflog ( smpi );

bdm/stat/libBM.cpp

@@ -234 +234 @@
 }
 
-void compositepdf::setindices(const RV &rv){
-        for (int i = 0;i < n;i++ ) {
-                // find ith rv in common rv
-                rvsinrv ( i ) = mpdfs ( i )->_rv().dataind ( rv );
-                // establish link between rv and rvc
-                rv.dataind(mpdfs(i)->_rvc(), irv_rvcs(i), irvcs_rv(i));
-        }
-}
-
 void BM::bayesB(const mat &Data){
         for(int t=0;t<Data.cols();t++){bayes(Data.get_col(t));}

bdm/stat/libBM.h

@@ -247 +247 @@
 };
 
-//!DataLink is a connection between epdf and its superordinate (Up)
+//!DataLink is a connection between an epdf and its superordinate epdf (Up)
 //! It is assumed that my val is fully present in "Up"
-class datalink {
+class datalink_e2e {
 protected:
         //! Remember how long val should be
@@ -257 +257 @@
         //! val-to-val link, indeces of the upper val
         ivec v2v_up;
-        public:
+public:
         //! Constructor
-        datalink ( const RV &rv, const RV &rv_up ) :
-                valsize(rv.count()), valupsize(rv_up.count()), v2v_up ( rv.dataind ( rv_up ) )  {it_assert_debug(v2v_up.length()==valsize,"rv is not fully in rv_up");}
+        datalink_e2e ( const RV &rv, const RV &rv_up ) :
+                        valsize ( rv.count() ), valupsize ( rv_up.count() ), v2v_up ( rv.dataind ( rv_up ) )  {
+                it_assert_debug ( v2v_up.length() ==valsize,"rv is not fully in rv_up" );
+        }
         //! Get val for myself from val of "Up"
-        vec get_val ( const vec &val_up ) {it_assert_debug(valupsize==val_up.length(),"Wrong val_up"); return get_vec ( val_up,v2v_up );}
+        vec get_val ( const vec &val_up ) {it_assert_debug ( valupsize==val_up.length(),"Wrong val_up" ); return get_vec ( val_up,v2v_up );}
         //! Fill val of "Up" by my pieces
-        void fill_val ( vec &val_up, const vec &val ) {it_assert_debug(valsize==val.length(),"Wrong val");it_assert_debug(valupsize==val_up.length(),"Wrong val_up");set_subvector ( val_up, v2v_up, val );}
-};
-
-//!DataLink is a connection between mpdf and its superordinate (Up)
-//! This class links 
-class datalink_mpdf: public datalink {
-protected:
+        void fill_val ( vec &val_up, const vec &val ) {
+                it_assert_debug ( valsize==val.length(),"Wrong val" );
+                it_assert_debug ( valupsize==val_up.length(),"Wrong val_up" );
+                set_subvector ( val_up, v2v_up, val );
+        }
+};
+
+//! data link between
+class datalink_m2e: public datalink_e2e {
+protected:
+        //! Remember how long cond should be
+        int condsize;
         //!upper_val-to-local_cond link, indeces of the upper val
         ivec v2c_up;
         //!upper_val-to-local_cond link, ideces of the local cond
         ivec v2c_lo;
+
+public:
+        //! Constructor
+        datalink_m2e ( const RV &rv,  const RV &rvc, const RV &rv_up ) :
+                        datalink_e2e ( rv,rv_up ), condsize ( rvc.count() ) {
+                //establish v2c connection
+                rvc.dataind ( rv_up, v2c_lo, v2c_up );
+        }
+        //!Construct condition
+        vec get_cond ( const vec &val_up ) {
+                vec tmp ( condsize );
+                set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) );
+                return tmp;
+        }
+};
+//!DataLink is a connection between mpdf and its superordinate (Up)
+//! This class links
+class datalink_m2m: public datalink_m2e {
+protected:
         //!cond-to-cond link, indeces of the upper cond
         ivec c2c_up;
         //!cond-to-cond link, indeces of the local cond
         ivec c2c_lo;
-        //! size of cond
-        int csize;
 public:
         //! Constructor
-        datalink_mpdf ( const mpdf &Me, const mpdf &Up ) :
-                        datalink ( Me._rv(),Up._rv() ) {
-                //establish v2c connection
-                Me._rvc().dataind ( Up._rv(), v2c_lo, v2c_up );
+        datalink_m2m ( const RV &rv, const RV &rvc, const RV &rv_up, const RV &rvc_up ) :
+                        datalink_m2e ( rv, rvc, rv_up) {
                 //establish c2c connection
-                Me._rvc().dataind ( Up._rvc(), c2c_lo, c2c_up );
-                csize = Me._rvc().count();
+                rvc.dataind ( rvc_up, c2c_lo, c2c_up );
+                it_assert_debug(c2c_lo.length()+v2c_lo.length()==condsize, "cond is not fully given");
         }
         //! Get cond for myself from val and cond of "Up"
         vec get_cond ( const vec &val_up, const vec &cond_up ) {
-                vec tmp(csize);
-                set_subvector(tmp,v2c_lo,val_up(v2c_up));
-                set_subvector(tmp,c2c_lo,cond_up(c2c_up));
+                vec tmp ( condsize );
+                set_subvector ( tmp,v2c_lo,val_up ( v2c_up ) );
+                set_subvector ( tmp,c2c_lo,cond_up ( c2c_up ) );
                 return tmp;
         }
-        //! Fill 
+        //! Fill
 
 };
@@ -304 +326 @@
 /*! \brief Unconditional mpdf, allows using epdf in the role of mpdf.
 
-WARNING: the class does not check validity of the \c ep pointer nor its existence.
 */
 class mepdf : public mpdf {
@@ -314 +335 @@
 
 //!\brief Abstract composition of pdfs, a base for specific classes
+//!this abstract class is common to epdf and mpdf
 class compositepdf {
 protected:
@@ -320 +342 @@
         //! Elements of composition
         Array<mpdf*> mpdfs;
-        //! Indeces of rvs in common rv
-        Array<ivec> rvsinrv;
-        //! Indeces of rvc in common rv
-        Array<ivec> irvcs_rv;
-        //! Indeces of common rv in rvc
-        Array<ivec> irv_rvcs;
-public:
-        compositepdf ( Array<mpdf*> A0 ) : n ( A0.length() ), mpdfs ( A0 ), rvsinrv ( n ), irvcs_rv ( n ),irv_rvcs ( n ) {};
+public:
+        compositepdf ( Array<mpdf*> A0 ) : n ( A0.length() ), mpdfs ( A0 ) {};
         //! find common rv, flag \param checkoverlap modifies whether overlaps are acceptable
         RV getrv ( bool checkoverlap=false );
         //! common rvc of all mpdfs is written to rvc
         void setrvc ( const RV &rv, RV &rvc );
-        //! fill all rv*inrv* according to
-        void setindices ( const RV &rv );
 };
 

bdm/stat/libEF.h

@@ -380 +380 @@
         //! Set \c A and \c R
         void set_parameters ( const  mat &A, const vec &mu0, const sq_T &R );
-        //!Generate one sample of the posterior
-        vec samplecond (const vec &cond, double &lik );
-        //!Generate matrix of samples of the posterior
-        mat samplecond (const vec &cond, vec &lik, int n );
+//      //!Generate one sample of the posterior
+//      vec samplecond (const vec &cond, double &lik );
+//      //!Generate matrix of samples of the posterior
+//      mat samplecond (const vec &cond, vec &lik, int n );
         //! Set value of \c rvc . Result of this operation is stored in \c epdf use function \c _ep to access it.
         void condition (const vec &cond );
@@ -567 +567 @@
 }
 
-template<class sq_T>
-vec mlnorm<sq_T>::samplecond (const  vec &cond, double &lik ) {
-        this->condition ( cond );
-        vec smp = epdf.sample();
-        lik = epdf.eval ( smp );
-        return smp;
-}
-
-template<class sq_T>
-mat mlnorm<sq_T>::samplecond (const vec &cond, vec &lik, int n ) {
-        int i;
-        int dim = rv.count();
-        mat Smp ( dim,n );
-        vec smp ( dim );
-        this->condition ( cond );
-
-        for ( i=0; i<n; i++ ) {
-                smp = epdf.sample();
-                lik ( i ) = epdf.eval ( smp );
-                Smp.set_col ( i ,smp );
-        }
-
-        return Smp;
-}
+// template<class sq_T>
+// vec mlnorm<sq_T>::samplecond (const  vec &cond, double &lik ) {
+//      this->condition ( cond );
+//      vec smp = epdf.sample();
+//      lik = epdf.eval ( smp );
+//      return smp;
+// }
+
+// template<class sq_T>
+// mat mlnorm<sq_T>::samplecond (const vec &cond, vec &lik, int n ) {
+//      int i;
+//      int dim = rv.count();
+//      mat Smp ( dim,n );
+//      vec smp ( dim );
+//      this->condition ( cond );
+// 
+//      for ( i=0; i<n; i++ ) {
+//              smp = epdf.sample();
+//              lik ( i ) = epdf.eval ( smp );
+//              Smp.set_col ( i ,smp );
+//      }
+// 
+//      return Smp;
+// }
 
 template<class sq_T>

tests/datalink_test.cpp

@@ -10 +10 @@
 using std::endl;
 
-int main()
-{
-        RV a = RV ( "{a }","2");
-        RV b = RV ( "{b }");
-        RV c = RV ( "{c }");
-        
-        mpdf M1(concat(a,b),c);
-        mpdf M2(a, concat(b,c));
-        
-        datalink_mpdf dl(M2,M1);
-        
-        vec val("1 1.5 2");
-        vec cond("3");
-        
+int main() {
+        RV a = RV ( "{a }","2" );
+        RV b = RV ( "{b }" );
+        RV c = RV ( "{c }" );
+
+        mpdf M1 ( concat ( a,b ),c );
+        mpdf M2 ( a, concat ( b,c ) );
+
+        datalink_m2m dl ( a,M2._rvc(),M1._rv(), M1._rvc() );
+
+        vec val ( "1 1.5 2" );
+        vec cond ( "3" );
+
         cout << "val: " << val  << endl;
         cout << "cond: " << cond  << endl;
-        
-        cout << "lo val: " << dl.get_val(val) <<endl;
-        cout << "lo cond: " << dl.get_cond(val,cond) <<endl;
-        
+
+        cout << "lo val: " << dl.get_val ( val ) <<endl;
+        cout << "lo cond: " << dl.get_cond ( val,cond ) <<endl;
+
         //getchar();
         return 0;