mlpack  1.0.12
phi.hpp
Go to the documentation of this file.
1 
15 #ifndef __MLPACK_METHODS_MOG_PHI_HPP
16 #define __MLPACK_METHODS_MOG_PHI_HPP
17 
18 #include <mlpack/core.hpp>
19 
20 namespace mlpack {
21 namespace gmm {
22 
38 inline double phi(const double x, const double mean, const double var)
39 {
40  return exp(-1.0 * ((x - mean) * (x - mean) / (2 * var)))
41  / sqrt(2 * M_PI * var);
42 }
43 
60 inline double phi(const arma::vec& x,
61  const arma::vec& mean,
62  const arma::mat& cov)
63 {
64  arma::vec diff = mean - x;
65 
66  // Parentheses required for Armadillo 3.0.0 bug.
67  arma::vec exponent = -0.5 * (trans(diff) * inv(cov) * diff);
68 
69  // TODO: What if det(cov) < 0?
70  return pow(2 * M_PI, (double) x.n_elem / -2.0) * pow(det(cov), -0.5) *
71  exp(exponent[0]);
72 }
73 
86 inline double phi(const arma::vec& x,
87  const arma::vec& mean,
88  const arma::mat& cov,
89  const std::vector<arma::mat>& d_cov,
90  arma::vec& g_mean,
91  arma::vec& g_cov)
92 {
93  // We don't call out to another version of the function to avoid inverting the
94  // covariance matrix more than once.
95  arma::mat cinv = inv(cov);
96 
97  arma::vec diff = mean - x;
98  // Parentheses required for Armadillo 3.0.0 bug.
99  arma::vec exponent = -0.5 * (trans(diff) * inv(cov) * diff);
100 
101  long double f = pow(2 * M_PI, (double) x.n_elem / 2) * pow(det(cov), -0.5)
102  * exp(exponent[0]);
103 
104  // Calculate the g_mean values; this is a (1 x dim) vector.
105  arma::vec invDiff = cinv * diff;
106  g_mean = f * invDiff;
107 
108  // Calculate the g_cov values; this is a (1 x (dim * (dim + 1) / 2)) vector.
109  for (size_t i = 0; i < d_cov.size(); i++)
110  {
111  arma::mat inv_d = cinv * d_cov[i];
112 
113  g_cov[i] = f * dot(d_cov[i] * invDiff, invDiff) +
114  accu(inv_d.diag()) / 2;
115  }
116 
117  return f;
118 }
119 
130 inline void phi(const arma::mat& x,
131  const arma::vec& mean,
132  const arma::mat& cov,
133  arma::vec& probabilities)
134 {
135  // Column i of 'diffs' is the difference between x.col(i) and the mean.
136  arma::mat diffs = x - (mean * arma::ones<arma::rowvec>(x.n_cols));
137 
138  // Now, we only want to calculate the diagonal elements of (diffs' * cov^-1 *
139  // diffs). We just don't need any of the other elements. We can calculate
140  // the right hand part of the equation (instead of the left side) so that
141  // later we are referencing columns, not rows -- that is faster.
142  arma::mat rhs = -0.5 * inv(cov) * diffs;
143  arma::vec exponents(diffs.n_cols); // We will now fill this.
144  for (size_t i = 0; i < diffs.n_cols; i++)
145  exponents(i) = exp(accu(diffs.unsafe_col(i) % rhs.unsafe_col(i)));
146 
147  probabilities = pow(2 * M_PI, (double) mean.n_elem / -2.0) *
148  pow(det(cov), -0.5) * exponents;
149 }
150 
151 }; // namespace gmm
152 }; // namespace mlpack
153 
154 #endif
Linear algebra utility functions, generally performed on matrices or vectors.
Definition: load.hpp:23
#define M_PI
Definition: prereqs.hpp:42
double phi(const double x, const double mean, const double var)
Calculates the univariate Gaussian probability density function.
Definition: phi.hpp:38