Changeset 205 for bdm/estim/merger.h

Timestamp:

11/10/08 15:40:30 (16 years ago)

Author:

smidl

Message:

merger posledni verze

Files:

• bdm/estim/merger.h (modified) (4 diffs)

bdm/estim/merger.h

@@ -43 +43 @@
         //!Prior on the log-normal merging model
         double beta;
+        //! Projection to empirical density
+        eEmp eSmp;
+
 public:
 //!Default constructor
         merger ( const Array<mpdf*> &S ) :
                         compositepdf ( S ), epdf ( getrv ( false ) ),
-                        Mix ( Array<BMEF*> ( 0 ),vec ( 0 ) ), dls ( n ), rvzs ( n ), zdls ( n ) {
+                        Mix ( Array<BMEF*> ( 0 ),vec ( 0 ) ), dls ( n ), rvzs ( n ), zdls ( n ), eSmp(rv,0) {
                 RV ztmp;
                 // Extend rv by rvc!
@@ -68 +71 @@
         }
 //! Set internal parameters used in approximation
-        void set_parameters ( double beta0, int Ns0, int Nc0 ) {        beta=beta0;Ns=Ns0;Nc=Nc0;}
+        void set_parameters ( double beta0, int Ns0, int Nc0 ) {beta=beta0;Ns=Ns0;Nc=Nc0;eSmp.set_n(Ns0,false);}
 //!Initialize the proposal density. This function must be called before merge()!
         void init() { ////////////// NOT FINISHED
@@ -90 +93 @@
                 return Mix.logpred ( dtf );
         }
-        vec mean() const {return Mix._epdf().mean();}
+        vec mean() const {
+                const Vec<double> &w = eSmp._w();
+                const Array<vec> &S = eSmp._samples();
+                vec tmp=zeros ( rv.count() ); 
+                for ( int i=0; i<Ns; i++ ) {
+                        tmp+=w ( i ) *S ( i );
+                }
+                return tmp;
+        }
+        mat variance() const {
+                const vec &w = eSmp._w();
+                const Array<vec> &S = eSmp._samples();
+                
+                vec mea = mean(); 
+                
+                cout << sum(w) << "," << w*w <<endl;
+                
+                mat Tmp=zeros(rv.count(), rv.count());
+                for ( int i=0; i<Ns; i++ ) {
+                        Tmp+=w ( i ) *outer_product(S ( i ), S(i));
+                }
+                return Tmp-outer_product(mea,mea);
+        }
 //! for future use
         virtual ~merger() {
@@ -101 +126 @@
 //! Access function
         MixEF& _Mix() {return Mix;}
+//! Access function
+        eEmp& _Smp() {return eSmp;}
 };