Changeset 229 for bdm/stat

Timestamp:

01/15/09 10:53:54 (15 years ago)

Author:

smidl

Message:

epdf has a new function: variance()

Location:

bdm/stat

Files:

• emix.h (modified) (2 diffs)
• libBM.h (modified) (1 diff)
• libEF.h (modified) (9 diffs)

bdm/stat/emix.h

@@ -104 +104 @@
                 for ( i = 0;i < w.length();i++ ) {mu += w ( i ) * Coms ( i )->mean(); }
                 return mu;
+        }
+        vec variance() const {
+                //non-central moment
+                vec mom2 = zeros(rv.count());
+                for ( int i = 0;i < w.length();i++ ) {mom2 += w ( i ) * pow(Coms ( i )->mean(),2); }
+                //central moment
+                return mom2-pow(mean(),2);
         }
         double evallog ( const vec &val ) const {
@@ -243 +250 @@
                 return tmp;
         }
+        vec variance() const {
+                vec tmp ( rv.count() ); //second moment
+                for ( int i=0;i<epdfs.length();i++ ) {
+                        vec pom = epdfs ( i )->mean();
+                        dls ( i )->fill_val ( tmp, pow(pom,2) );
+                }
+                return tmp-pow(mean(),2);
+        }
         vec sample() const {
                 vec tmp ( rv.count() );

bdm/stat/libBM.h

@@ -197 +197 @@
         virtual vec mean() const =0;
 
+        //! return expected variance (not covariance!)
+        virtual vec variance() const = 0;
+        
         //! Destructor for future use;
         virtual ~epdf() {};

bdm/stat/libEF.h

@@ -132 +132 @@
         double lognc () const;
         vec mean() const {return mu;}
+        vec variance() const {return diag(R.to_mat());}
 //      mlnorm<sq_T>* condition ( const RV &rvn ) const ;
         mpdf* condition ( const RV &rvn ) const ;
@@ -192 +193 @@
         vec sample() const;
         vec mean() const;
+        vec variance() const{it_error("Not implemented"); return vec(0);};
         void mean_mat ( mat &M, mat&R ) const;
         //! In this instance, val= [theta, r]. For multivariate instances, it is stored columnwise val = [theta_1 theta_2 ... r_1 r_2 ]
@@ -220 +222 @@
         //!sufficient statistics
         vec beta;
+        //!speedup variable
+        double gamma;
 public:
         //!Default constructor
@@ -226 +230 @@
         eDirich ( const eDirich &D0 ) : eEF ( D0.rv ),beta ( D0.beta ) {};
         vec sample() const {it_error ( "Not implemented" );return vec_1 ( 0.0 );};
-        vec mean() const {return beta/sum ( beta );};
+        vec mean() const {return beta/gamma;};
+        vec variance() const {return elem_mult(beta,(beta+1))/ (gamma*(gamma+1));}
         //! In this instance, val is ...
         double evallog_nn ( const vec &val ) const {double tmp; tmp=( beta-1 ) *log ( val );            it_assert_debug(std::isfinite(tmp),"Infinite value");
@@ -248 +253 @@
                 }
                 beta= beta0;
+                gamma = sum(beta);
         }
 };
@@ -327 +333 @@
         //! Returns poiter to alpha and beta. Potentially dengerous: use with care!
         void _param ( vec* &a, vec* &b ) {a=&alpha;b=&beta;};
-        vec mean() const {vec pom ( alpha ); pom/=beta; return pom;}
+        vec mean() const {return elem_div(alpha,beta);}
+        vec variance() const {return elem_div(alpha,elem_mult(beta,beta)); }
 };
 
@@ -365 +372 @@
                 void _param ( vec* &a, vec* &b ) {a=alpha;b=beta;};
                 vec mean() const {return elem_div(*beta,*alpha-1);}
+                vec variance() const {vec mea=mean(); return elem_div(elem_mult(mea,mea),*alpha-2);}
 };
 /*
@@ -416 +424 @@
                 lnk = log ( nk );
         }
-        vec mean() const {vec pom=high; pom-=low; pom/=2.0; return pom;}
+        vec mean() const {return (high-low)/2.0;}
+        vec variance() const {return (pow(high,2)+pow(low,2)+elem_mult(high,low))/3.0;}
 };
 
@@ -646 +655 @@
                 return pom;
         }
+        vec variance() const {
+                vec pom=zeros ( rv.count() );
+                for ( int i=0;i<n;i++ ) {pom+=pow(samples ( i ),2) *w ( i );}
+                return pom-pow(mean(),2);
+        }
 };