Changeset 1014

Timestamp:

05/27/10 23:07:57 (14 years ago)

Author:

smidl

Message:

Mixtures example

Files:

• applications/bdmtoolbox/tutorial/userguide/mixef_basic.m (modified) (4 diffs)
• library/bdm/estim/mixtures.cpp (modified) (3 diffs)
• library/bdm/estim/mixtures.h (modified) (7 diffs)

applications/bdmtoolbox/tutorial/userguide/mixef_basic.m

@@ -5 +5 @@
 com.rgr = RV({},[],[]);
 
-Data = [randn(2,50) randn(2,50)+10*ones(2,50)];
+Data = [randn(2,50) randn(2,50)+[5*ones(1,50); 10*ones(1,50)]];
 
-[Mix0,P0]=mixef_init(Data,com);
+[Mix0,P0]=mixef_init(Data,com,4);
 % show predictor
 
@@ -16 +16 @@
 hold on
 epdf_2dplot(Pred);
+hold on
+for i=1:length(Pred.pdfs)
+    mea=epdf_mean(Pred.pdfs{i});
+    plot(mea(1),mea(2),'+r');
+end
 
 % Do Quasi Bayes
@@ -32 +37 @@
     % normalize weights
     w = exp(log_w_nn-max(log_w_nn));
-    w = w/sum(w)
+    w = w/sum(w);
     
     for c=1:n
@@ -56 +61 @@
 epdf_2dplot(PredQB);
 
+PredQBc = bm_epredictor(MixQBc);
+figure(3);
+hold off
+plot(Data(1,:), Data(2,:),'.');
+hold on
+epdf_2dplot(PredQBc);
+
 %% relations between components
 Batta_dist=zeros(n);

library/bdm/estim/mixtures.cpp

@@ -18 +18 @@
 //      Coms(0)->set_evalll(false);
         Coms ( 0 )->bayes_batch ( Data );
-        // Flatten it to its original shape
-        shared_ptr<BMEF> SharpCom((BMEF*)Coms(0)->_copy());
+
         Coms ( 0 )->flatten ( Com0 ); 
 
@@ -34 +33 @@
                 } else { // pick at random
                         int ind = (int) floor ( ndat * UniRNG.sample() );
-                        Coms ( i )->bayes_weighted ( Data.get_col ( ind ), empty_vec, ndat/Coms.length() );
+                        Coms ( i )->bayes_weighted ( Data.get_col ( ind ), empty_vec, ndat );
+                        Coms (i)->flatten(Com0,ndat/Coms.length());
                 }
                 //sharpen to the sharp component
                 //Coms ( i )->flatten ( SharpCom.get(), 1.0/Coms.length() );
         }
+        MixEF_options old_opt =options;
+        MixEF_options ini_opt=options;
+        ini_opt.method = EM;
+        ini_opt.max_niter= 1;
+        bayes_batch(Data, empty_vec);
+        
+        for ( i = 0; i < Coms.length(); i++ ) {
+                Coms (i)->flatten(Com0);
+        }
+        
+        options = old_opt;
 }
 
@@ -83 +94 @@
 
                         maxll = max ( ll, maxi );
-                        switch ( method ) {
+                        switch ( options.method ) {
                         case QB:
                                 w = exp ( ll - maxll );

library/bdm/estim/mixtures.h

@@ -39 +39 @@
 This class uses EM-style type algorithms for estimation of its parameters. Under this simplification, the posterior density is a product of exponential family members, hence under EM-style approximate estimation this class itself belongs to the exponential family.
 
-TODO: Extend BM to use rvc.
+Two methods are provided:
+ - QB where the probability of being from one component is computed using predictors,
+ - EM where the data is assigned to component with the highest likelihood (winner takes all)
+
+These methods are stored in attribute options:
+ - method: ["EM","QB"] estimation method as mentioned above, QB is default
+ - max_niter: maximum of iterations in bayes_batch()
+
 */
 class MixEF: public BMEF {
@@ -60 +67 @@
         ////!indices of component rvc in common rvc
 
-        //! Flag for a method that is used in the inference
-        MixEF_METHOD method;
+        class MixEF_options: public root{
+                public:
+                //! Flag for a method that is used in the inference
+                MixEF_METHOD method;
 
-        //! maximum number of iterations
-        int max_niter;
+                //! maximum number of iterations
+                int max_niter;
+        
+        MixEF_options():method(QB),max_niter(10){};
+                
+                void from_setting(const Setting &set){
+                        string meth;
+                        UI::get(meth,set,"method",UI::optional);
+                        if (meth=="EM") {method=EM;}
+                        max_niter =10;
+                        UI::get(max_niter,set,"max_niter",UI::optional);
+                };
+                void to_setting(Setting &set)const{
+                        string meth=(method==EM ? "EM" : "QB");
+                        UI::save(meth,set,"method");
+                        UI::save(max_niter,set,"max_niter");
+                };
+        };
+        
+        MixEF_options options;
 public:
         //! Full constructor
         MixEF ( const Array<BMEF*> &Coms0, const vec &alpha0 ) :
                         BMEF ( ), Coms ( Coms0.length() ),
-                        weights (), est(*this), method ( QB ), max_niter(10) {
+                        weights (), est(*this), options() {
                 for ( int i = 0; i < Coms0.length(); i++ ) {
                         Coms ( i ) = ( BMEF* ) Coms0 ( i )->_copy();
@@ -80 +107 @@
         MixEF () :
                         BMEF ( ), Coms ( 0 ),
-                        weights (), est(*this), method ( QB ), max_niter(10) {
+                        weights (), est(*this), options() {
         }
         //! Copy constructor
         MixEF ( const MixEF &M2 ) : BMEF ( ),  Coms ( M2.Coms.length() ),
-                        weights ( M2.weights ), est(*this), method ( M2.method ), max_niter(M2.max_niter) {
+                        weights ( M2.weights ), est(*this), options(M2.options) {
                 for ( int i = 0; i < M2.Coms.length(); i++ ) {
                         Coms ( i ) = (BMEF*) M2.Coms ( i )->_copy();
@@ -94 +121 @@
         //! \param Data Data on which the initialization will be done
         //! \param c Initial number of components, default=5
+        //! when the number of data records (columns of Data) is equal to the number of requested components, each data is used, otherwise, they are picked randomly.
         void init ( BMEF* Com0, const mat &Data, const int c = 5 );
         //Destructor
@@ -119 +147 @@
         //!Set which method is to be used
         void set_method ( MixEF_METHOD M ) {
-                method = M;
+                options.method = M;
         }
 
@@ -126 +154 @@
         UI::save ( Coms, set, "Coms" );
         UI::save ( &weights, set, "weights" );
+        UI::save (options, set, "options");
 }
 /*! \brief reads data from setting
@@ -133 +162 @@
         UI::get ( Coms, set, "Coms" );
         UI::get ( weights, set, "weights" );
+        UI::get (options, set, "options",UI::optional);
 }
 };