MLPACK  1.0.4
phi.hpp
Go to the documentation of this file.
00001 
00023 #ifndef __MLPACK_METHODS_MOG_PHI_HPP
00024 #define __MLPACK_METHODS_MOG_PHI_HPP
00025 
00026 #include <mlpack/core.hpp>
00027 
00028 namespace mlpack {
00029 namespace gmm {
00030 
00046 inline double phi(const double x, const double mean, const double var)
00047 {
00048   return exp(-1.0 * ((x - mean) * (x - mean) / (2 * var)))
00049       / sqrt(2 * M_PI * var);
00050 }
00051 
00068 inline double phi(const arma::vec& x,
00069                   const arma::vec& mean,
00070                   const arma::mat& cov)
00071 {
00072   arma::vec diff = mean - x;
00073 
00074   // Parentheses required for Armadillo 3.0.0 bug.
00075   arma::vec exponent = -0.5 * (trans(diff) * inv(cov) * diff);
00076 
00077   // TODO: What if det(cov) < 0?
00078   return pow(2 * M_PI, (double) x.n_elem / -2.0) * pow(det(cov), -0.5) *
00079       exp(exponent[0]);
00080 }
00081 
00094 inline double phi(const arma::vec& x,
00095                   const arma::vec& mean,
00096                   const arma::mat& cov,
00097                   const std::vector<arma::mat>& d_cov,
00098                   arma::vec& g_mean,
00099                   arma::vec& g_cov)
00100 {
00101   // We don't call out to another version of the function to avoid inverting the
00102   // covariance matrix more than once.
00103   arma::mat cinv = inv(cov);
00104 
00105   arma::vec diff = mean - x;
00106   // Parentheses required for Armadillo 3.0.0 bug.
00107   arma::vec exponent = -0.5 * (trans(diff) * inv(cov) * diff);
00108 
00109   long double f = pow(2 * M_PI, (double) x.n_elem / 2) * pow(det(cov), -0.5)
00110       * exp(exponent[0]);
00111 
00112   // Calculate the g_mean values; this is a (1 x dim) vector.
00113   arma::vec invDiff = cinv * diff;
00114   g_mean = f * invDiff;
00115 
00116   // Calculate the g_cov values; this is a (1 x (dim * (dim + 1) / 2)) vector.
00117   for (size_t i = 0; i < d_cov.size(); i++)
00118   {
00119     arma::mat inv_d = cinv * d_cov[i];
00120 
00121     g_cov[i] = f * dot(d_cov[i] * invDiff, invDiff) +
00122         accu(inv_d.diag()) / 2;
00123   }
00124 
00125   return f;
00126 }
00127 
00138 inline void phi(const arma::mat& x,
00139                 const arma::vec& mean,
00140                 const arma::mat& cov,
00141                 arma::vec& probabilities)
00142 {
00143   // Column i of 'diffs' is the difference between x.col(i) and the mean.
00144   arma::mat diffs = x - (mean * arma::ones<arma::rowvec>(x.n_cols));
00145 
00146   // Now, we only want to calculate the diagonal elements of (diffs' * cov^-1 *
00147   // diffs).  We just don't need any of the other elements.  We can calculate
00148   // the right hand part of the equation (instead of the left side) so that
00149   // later we are referencing columns, not rows -- that is faster.
00150   arma::mat rhs = -0.5 * inv(cov) * diffs;
00151   arma::vec exponents(diffs.n_cols); // We will now fill this.
00152   for (size_t i = 0; i < diffs.n_cols; i++)
00153     exponents(i) = exp(accu(diffs.unsafe_col(i) % rhs.unsafe_col(i)));
00154 
00155   probabilities = pow(2 * M_PI, (double) mean.n_elem / -2.0) *
00156       pow(det(cov), -0.5) * exponents;
00157 }
00158 
00159 }; // namespace gmm
00160 }; // namespace mlpack
00161 
00162 #endif