IT++ Logo

mog_diag.cpp

Go to the documentation of this file.
00001 
00030 #include <itpp/base/math/log_exp.h>
00031 #include <itpp/stat/mog_diag.h>
00032 #include <cstdlib>
00033 
00034 
00035 namespace itpp {
00036 
00037   double MOG_diag::log_lhood_single_gaus_internal(const double * c_x_in, const int k) const {
00038 
00039     const double * c_mean = c_means[k];
00040     const double * c_diag_cov_inv_etc = c_diag_covs_inv_etc[k];
00041 
00042     double acc = 0.0;
00043 
00044     for(int d=0; d<D; d++) {
00045       double tmp_val = c_x_in[d] - c_mean[d];
00046       acc += (tmp_val*tmp_val) * c_diag_cov_inv_etc[d];
00047     }
00048     return(c_log_det_etc[k] - acc);
00049   }
00050 
00051 
00052   double MOG_diag::log_lhood_single_gaus_internal(const vec &x_in, const int k) const {
00053     return log_lhood_single_gaus_internal(x_in._data(), k);
00054   }
00055 
00056 
00057   double MOG_diag::log_lhood_single_gaus(const double * c_x_in, const int k) const {
00058     if(do_checks) {
00059       it_assert(valid, "MOG_diag::log_lhood_single_gaus(): model not valid");
00060       it_assert( ( (k>=0) && (k<K) ), "MOG::log_lhood_single_gaus(): k specifies a non-existant Gaussian");
00061     }
00062     return log_lhood_single_gaus_internal(c_x_in,k);
00063   }
00064 
00065 
00066   double MOG_diag::log_lhood_single_gaus(const vec &x_in, const int k) const {
00067     if(do_checks) {
00068       it_assert(valid, "MOG_diag::log_lhood_single_gaus(): model not valid");
00069       it_assert(check_size(x_in), "MOG_diag::log_lhood_single_gaus(): x has wrong dimensionality");
00070       it_assert( ( (k>=0) && (k<K) ), "MOG::log_lhood_single_gaus(): k specifies a non-existant Gaussian");
00071     }
00072     return log_lhood_single_gaus_internal(x_in._data(),k);
00073   }
00074 
00075 
00076   double MOG_diag::log_lhood_internal(const double * c_x_in) {
00077 
00078     bool danger = paranoid;
00079 
00080     for(int k=0;k<K;k++)  {
00081       double tmp = c_log_weights[k] + log_lhood_single_gaus_internal(c_x_in,k);
00082       c_tmpvecK[k] = tmp;
00083 
00084       if(tmp >= log_max_K)  danger = true;
00085     }
00086 
00087 
00088     if(danger) {
00089       double log_sum = c_tmpvecK[0];  for(int k=1; k<K; k++)  log_sum = log_add( log_sum, c_tmpvecK[k] );
00090       return(log_sum);
00091     }
00092     else {
00093       double sum = 0.0; for(int k=0;k<K;k++) sum += std::exp(c_tmpvecK[k]);
00094       return(std::log(sum));
00095     }
00096   }
00097 
00098 
00099   double MOG_diag::log_lhood_internal(const vec &x_in)  {
00100     return log_lhood_internal(x_in._data());
00101   }
00102 
00103 
00104   double MOG_diag::log_lhood(const vec &x_in) {
00105     if(do_checks) {
00106       it_assert(valid, "MOG_diag::log_lhood(): model not valid");
00107       it_assert(check_size(x_in), "MOG_diag::log_lhood(): x has wrong dimensionality");
00108     }
00109     return log_lhood_internal(x_in._data());
00110   }
00111 
00112 
00113   double MOG_diag::log_lhood(const double * c_x_in) {
00114     if(do_checks) {
00115       it_assert(valid, "MOG_diag::log_lhood(): model not valid");
00116       it_assert( (c_x_in != 0), "MOG_diag::log_lhood(): c_x_in is a null pointer");
00117     }
00118 
00119     return log_lhood_internal(c_x_in);
00120   }
00121 
00122 
00123   double MOG_diag::lhood_internal(const double * c_x_in) {
00124 
00125     bool danger = paranoid;
00126 
00127     for(int k=0;k<K;k++)  {
00128       double tmp = c_log_weights[k] + log_lhood_single_gaus_internal(c_x_in,k);
00129       c_tmpvecK[k] = tmp;
00130 
00131       if(tmp >= log_max_K)  danger = true;
00132     }
00133 
00134 
00135     if(danger) {
00136       double log_sum = c_tmpvecK[0];  for(int k=1; k<K; k++)  log_sum = log_add( log_sum, c_tmpvecK[k] );
00137       return(trunc_exp(log_sum));
00138     }
00139     else {
00140       double sum = 0.0; for(int k=0;k<K;k++) sum += std::exp(c_tmpvecK[k]);
00141       return(sum);
00142     }
00143   }
00144 
00145   double MOG_diag::lhood_internal(const vec &x_in) { return lhood_internal(x_in._data()); }
00146 
00147   double MOG_diag::lhood(const vec &x_in) {
00148     if(do_checks) {
00149       it_assert(valid, "MOG_diag::lhood(): model not valid");
00150       it_assert(check_size(x_in), "MOG_diag::lhood(): x has wrong dimensionality");
00151     }
00152     return lhood_internal(x_in._data());
00153   }
00154 
00155 
00156   double MOG_diag::lhood(const double * c_x_in) {
00157     if(do_checks) {
00158       it_assert(valid, "MOG_diag::lhood(): model not valid");
00159       it_assert( (c_x_in != 0), "MOG_diag::lhood(): c_x_in is a null pointer");
00160     }
00161 
00162     return lhood_internal(c_x_in);
00163   }
00164 
00165 
00166   double MOG_diag::avg_log_lhood(const double ** c_x_in, const int N) {
00167     if(do_checks) {
00168       it_assert(valid, "MOG_diag::avg_log_lhood(): model not valid");
00169       it_assert( (c_x_in != 0), "MOG_diag::avg_log_lhood(): c_x_in is a null pointer");
00170       it_assert( (N >= 0), "MOG_diag::avg_log_lhood(): N is zero or negative");
00171     }
00172 
00173     double acc = 0.0;  for(int n=0;n<N;n++) acc += log_lhood_internal(c_x_in[n]);
00174     return(acc/N);
00175   }
00176 
00177 
00178   double MOG_diag::avg_log_lhood(const Array<vec> &X_in) {
00179     if(do_checks) {
00180       it_assert(valid, "MOG_diag::avg_log_lhood(): model not valid");
00181       it_assert(check_size(X_in), "MOG_diag::avg_log_lhood(): X is empty or at least one vector has the wrong dimensionality");
00182     }
00183     const int N = X_in.size();
00184     double acc = 0.0;
00185     for(int n=0;n<N;n++)  acc += log_lhood_internal(X_in(n)._data());
00186     return(acc/N);
00187   }
00188 
00189   void MOG_diag::zero_all_ptrs() {
00190     c_means             = 0;
00191     c_diag_covs         = 0;
00192     c_diag_covs_inv_etc = 0;
00193     c_weights           = 0;
00194     c_log_weights       = 0;
00195     c_log_det_etc       = 0;
00196     c_tmpvecK           = 0;
00197   }
00198 
00199 
00200   void MOG_diag::free_all_ptrs() {
00201     c_means             = disable_c_access(c_means);
00202     c_diag_covs         = disable_c_access(c_diag_covs);
00203     c_diag_covs_inv_etc = disable_c_access(c_diag_covs_inv_etc);
00204     c_weights           = disable_c_access(c_weights);
00205     c_log_weights       = disable_c_access(c_log_weights);
00206     c_log_det_etc       = disable_c_access(c_log_det_etc);
00207     c_tmpvecK           = disable_c_access(c_tmpvecK);
00208   }
00209 
00210 
00211   void MOG_diag::setup_means() {
00212     MOG_generic::setup_means();
00213     disable_c_access(c_means);
00214     c_means = enable_c_access(means);
00215   }
00216 
00217 
00218   void MOG_diag::setup_covs() {
00219     MOG_generic::setup_covs();
00220     if(full) return;
00221 
00222     disable_c_access(c_diag_covs);
00223     disable_c_access(c_diag_covs_inv_etc);
00224     disable_c_access(c_log_det_etc);
00225 
00226     c_diag_covs         = enable_c_access(diag_covs);
00227     c_diag_covs_inv_etc = enable_c_access(diag_covs_inv_etc);
00228     c_log_det_etc       = enable_c_access(log_det_etc);
00229   }
00230 
00231 
00232   void MOG_diag::setup_weights() {
00233     MOG_generic::setup_weights();
00234 
00235     disable_c_access(c_weights);
00236     disable_c_access(c_log_weights);
00237 
00238     c_weights = enable_c_access(weights);
00239     c_log_weights = enable_c_access(log_weights);
00240   }
00241 
00242 
00243   void MOG_diag::setup_misc() {
00244     disable_c_access(c_tmpvecK);
00245     tmpvecK.set_size(K);
00246     c_tmpvecK = enable_c_access(tmpvecK);
00247 
00248     MOG_generic::setup_misc();
00249     if(full) convert_to_diag_internal();
00250   }
00251 
00252 
00253   void MOG_diag::load(const std::string &name_in) {
00254     MOG_generic::load(name_in);
00255     if(full) convert_to_diag();
00256   }
00257 
00258 
00259   double ** MOG_diag::enable_c_access(Array<vec> & A_in) {
00260     int rows = A_in.size();
00261     double ** A = (double **)std::malloc(rows*sizeof(double *));
00262     if(A)  for(int row=0;row<rows;row++)  A[row] = A_in(row)._data();
00263     return(A);
00264   }
00265 
00266   int ** MOG_diag::enable_c_access(Array<ivec> & A_in) {
00267     int rows = A_in.size();
00268     int ** A = (int **)std::malloc(rows*sizeof(int *));
00269     if(A)  for(int row=0;row<rows;row++)  A[row] = A_in(row)._data();
00270     return(A);
00271   }
00272 
00273   double ** MOG_diag::disable_c_access(double ** A_in) { if(A_in) std::free(A_in); return(0); }
00274   int ** MOG_diag::disable_c_access(int ** A_in) { if(A_in) std::free(A_in); return(0); }
00275 
00276   double * MOG_diag::enable_c_access(vec & v_in) { return v_in._data(); }
00277   int * MOG_diag::enable_c_access(ivec & v_in) { return v_in._data(); }
00278 
00279   double * MOG_diag::disable_c_access(double * v_in) { return(0); }
00280   int * MOG_diag::disable_c_access(int * v_in) { return(0); }
00281 
00282 }
SourceForge Logo

Generated on Sat May 3 16:10:44 2008 for IT++ by Doxygen 1.5.5