summaryrefslogtreecommitdiff
path: root/utils
diff options
context:
space:
mode:
authorChris Dyer <cdyer@cs.cmu.edu>2012-03-06 23:20:16 -0500
committerChris Dyer <cdyer@cs.cmu.edu>2012-03-06 23:20:16 -0500
commit27e0de58bf49a4fc74bbf58718d1b89525a154a6 (patch)
treea31c4f4f62f65eb68d0e19cfe686b09867666d80 /utils
parentde34b1493df93169c991a1828f951ca5abc00cae (diff)
a few statistical helpers i'm using to figure some algorithms out
Diffstat (limited to 'utils')
-rw-r--r--utils/m.h45
-rw-r--r--utils/m_test.cc16
2 files changed, 61 insertions, 0 deletions
diff --git a/utils/m.h b/utils/m.h
index 5e45efee..dc881b36 100644
--- a/utils/m.h
+++ b/utils/m.h
@@ -4,6 +4,10 @@
#include <cassert>
#include <cmath>
#include <boost/math/special_functions/digamma.hpp>
+#include <boost/math/constants/constants.hpp>
+
+// TODO right now I sometimes assert that x is in the support of the distributions
+// should be configurable to return -inf instead
template <typename F>
struct M {
@@ -59,6 +63,47 @@ struct M {
return (alpha-1)*std::log(x)+(beta-1)*std::log(1-x) - log_beta_fn(alpha, beta);
}
+ // support x \in R
+ static inline F log_laplace_density(const F& x, const F& mu, const F& b) {
+ assert(b > 0.0);
+ return -std::log(2*b) - std::fabs(x - mu) / b;
+ }
+
+ // support x \in R
+ // this is NOT the "log normal" density, it is the log of the "normal density at x"
+ static inline F log_gaussian_density(const F& x, const F& mu, const F& var) {
+ assert(var > 0.0);
+ return -0.5 * std::log(var * 2 * boost::math::constants::pi<F>()) - (x - mu)*(x - mu) / (2 * var);
+ }
+
+ // (x1,x2) \in R^2
+ // parameterized in terms of two means, a two "variances", a correlation < 1
+ static inline F log_bivariate_gaussian_density(const F& x1, const F& x2,
+ const F& mu1, const F& mu2,
+ const F& var1, const F& var2,
+ const F& cor) {
+ assert(var1 > 0);
+ assert(var2 > 0);
+ assert(std::fabs(cor) < 1.0);
+ const F cor2 = cor*cor;
+ const F var1var22 = var1 * var2;
+ const F Z = 0.5 * std::log(var1var22 * (1 - cor2)) + std::log(2 * boost::math::constants::pi<F>());
+ return -Z -1.0 / (2 * (1 - cor2)) * ((x1 - mu1)*(x1-mu1) / var1 + (x2-mu2)*(x2-mu2) / var2 - 2*cor*(x1 - mu1)*(x2-mu2) / std::sqrt(var1var22));
+ }
+
+ // support x \in [a,b]
+ static inline F log_triangle_density(const F& x, const F& a, const F& b, const F& c) {
+ assert(a < b);
+ assert(a <= c);
+ assert(c <= b);
+ assert(x >= a);
+ assert(x <= b);
+ if (x <= c)
+ return std::log(2) + std::log(x - a) - std::log(b - a) - std::log(c - a);
+ else
+ return std::log(2) + std::log(b - x) - std::log(b - a) - std::log(b - c);
+ }
+
// note: this has been adapted so that 0 is in the support of the distribution
// support [0, 1, 2 ...)
static inline F log_yule_simon(unsigned x, const F& rho) {
diff --git a/utils/m_test.cc b/utils/m_test.cc
index fca8f895..c4d6a166 100644
--- a/utils/m_test.cc
+++ b/utils/m_test.cc
@@ -14,6 +14,22 @@ class MTest : public testing::Test {
virtual void TearDown() { }
};
+TEST_F(MTest, Densities) {
+ double px1 = Md::log_gaussian_density(1.0, 0.0, 1.0);
+ double px2 = Md::log_gaussian_density(-1.0, 0.0, 1.0);
+ double py1 = Md::log_laplace_density(1.0, 0.0, 1.0);
+ double py2 = Md::log_laplace_density(1.0, 0.0, 1.0);
+ double pz1 = Md::log_triangle_density(1.0, -2.0, 2.0, 0.0);
+ double pz2 = Md::log_triangle_density(1.0, -2.0, 2.0, 0.0);
+ cerr << px1 << " " << py1 << " " << pz2 << endl;
+ EXPECT_FLOAT_EQ(px1, px2);
+ EXPECT_FLOAT_EQ(py1, py2);
+ EXPECT_FLOAT_EQ(pz1, pz2);
+ double b1 = Md::log_bivariate_gaussian_density(1.0, -1.0, 0.0, 0.0, 1.0, 1.0, -0.8);
+ double b2 = Md::log_bivariate_gaussian_density(-1.0, 1.0, 0.0, 0.0, 1.0, 1.0, -0.8);
+ cerr << b1 << " " << b2 << endl;
+}
+
TEST_F(MTest, Poisson) {
double prev = 1.0;
double tot = 0;