Stan  1.3
probability, sampling & optimization
 All Classes Namespaces Files Functions Variables Typedefs Enumerator Friends Macros Pages
multi_student_t.hpp
Go to the documentation of this file.
1 #ifndef __STAN__PROB__DISTRIBUTIONS__MULTI_STUDENT_T_HPP__
2 #define __STAN__PROB__DISTRIBUTIONS__MULTI_STUDENT_T_HPP__
3 
4 #include <cstdlib>
5 
6 #include <boost/math/special_functions/gamma.hpp>
7 
8 #include <stan/math/error_handling.hpp>
9 #include <stan/math/matrix_error_handling.hpp>
10 #include <stan/prob/constants.hpp>
11 #include <stan/prob/traits.hpp>
12 #include <stan/prob/distributions/multivariate/continuous/multi_normal.hpp>
13 #include <stan/prob/distributions/univariate/continuous/inv_gamma.hpp>
14 
15 #include <boost/random/variate_generator.hpp>
16 
17 namespace stan {
18 
19  namespace prob {
20 
27  template <bool propto,
28  typename T_y, typename T_dof, typename T_loc, typename T_scale,
29  class Policy>
30  typename boost::math::tools::promote_args<T_y,T_dof,T_loc,T_scale>::type
31  multi_student_t_log(const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
32  const T_dof& nu,
33  const Eigen::Matrix<T_loc,Eigen::Dynamic,1>& mu,
34  const
35  Eigen::Matrix<T_scale,
36  Eigen::Dynamic,Eigen::Dynamic>& Sigma,
37  const Policy&) {
38  static const char* function = "stan::prob::multi_student_t(%1%)";
39 
40  using stan::math::check_size_match;
41  using stan::math::check_finite;
42  using stan::math::check_not_nan;
43  using stan::math::check_symmetric;
44  using stan::math::check_positive;
45  using stan::math::check_pos_definite;
46  using boost::math::tools::promote_args;
47  using boost::math::lgamma;
48 
49  typename promote_args<T_y,T_dof,T_loc,T_scale>::type lp(0.0);
50  if (!check_size_match(function,
51  y.size(), "Size of random variable",
52  mu.size(), "size of location parameter",
53  &lp, Policy()))
54  return lp;
55  if (!check_size_match(function,
56  y.size(), "Size of random variable",
57  Sigma.rows(), "rows of scale parameter",
58  &lp, Policy()))
59  return lp;
60  if (!check_size_match(function,
61  y.size(), "Size of random variable",
62  Sigma.cols(), "columns of scale parameter",
63  &lp, Policy()))
64  return lp;
65  if (!check_finite(function, mu, "Location parameter", &lp, Policy()))
66  return lp;
67  if (!check_not_nan(function, y, "Random variable", &lp, Policy()))
68  return lp;
69  if (!check_symmetric(function, Sigma, "Scale parameter", &lp, Policy()))
70  return lp;
71  if (!check_pos_definite(function, Sigma, "Scale parameter", &lp, Policy()))
72  return lp;
73 
74  // allows infinities
75  if (!check_not_nan(function, nu,
76  "Degrees of freedom parameter", &lp,
77  Policy()))
78  return lp;
79  if (!check_positive(function, nu,
80  "Degrees of freedom parameter", &lp,
81  Policy()))
82  return lp;
83 
84  using std::isinf;
85 
86  if (isinf(nu)) // already checked nu > 0
87  return multi_normal_log(y,mu,Sigma,Policy());
88 
89  double d = y.size();
90 
91  if (include_summand<propto,T_dof>::value) {
92  lp += lgamma(0.5 * (nu + d));
93  lp -= lgamma(0.5 * nu);
94  lp -= (0.5 * d) * log(nu);
95  }
96 
97  if (include_summand<propto>::value)
98  lp -= (0.5 * d) * LOG_PI;
99 
100  using stan::math::multiply;
101  using stan::math::dot_product;
102  using stan::math::subtract;
103  using Eigen::Array;
104  using stan::math::mdivide_left;
105  using stan::math::log_determinant;
106 
107 
108  if (include_summand<propto,T_scale>::value) {
109  lp -= 0.5*log_determinant(Sigma);
110  }
111 
112  if (include_summand<propto,T_y,T_dof,T_loc,T_scale>::value) {
113 
114  Eigen::Matrix<typename promote_args<T_y,T_loc>::type,
115  Eigen::Dynamic,
116  1> y_minus_mu = subtract(y,mu);
117  Eigen::Matrix<typename promote_args<T_scale,T_y,T_loc>::type,
118  Eigen::Dynamic,
119  1> invSigma_dy = mdivide_left(Sigma, y_minus_mu);
120  lp -= 0.5
121  * (nu + d)
122  * log(1.0 + dot_product(y_minus_mu,invSigma_dy) / nu);
123  }
124  return lp;
125  }
126 
127  template <bool propto,
128  typename T_y, typename T_dof, typename T_loc, typename T_scale>
129  inline
130  typename boost::math::tools::promote_args<T_y,T_dof,T_loc,T_scale>::type
131  multi_student_t_log(const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
132  const T_dof& nu,
133  const Eigen::Matrix<T_loc,Eigen::Dynamic,1>& mu,
134  const
135  Eigen::Matrix<T_scale,
136  Eigen::Dynamic,Eigen::Dynamic>& Sigma) {
137  return multi_student_t_log<propto>(y,nu,mu,Sigma,
138  stan::math::default_policy());
139  }
140 
141  template <typename T_y, typename T_dof, typename T_loc, typename T_scale,
142  class Policy>
143  inline
144  typename boost::math::tools::promote_args<T_y,T_dof,T_loc,T_scale>::type
145  multi_student_t_log(const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
146  const T_dof& nu,
147  const Eigen::Matrix<T_loc,Eigen::Dynamic,1>& mu,
148  const
149  Eigen::Matrix<T_scale,
150  Eigen::Dynamic,Eigen::Dynamic>& Sigma,
151  const Policy&) {
152  return multi_student_t_log<false>(y,nu,mu,Sigma,Policy());
153  }
154 
155 
156 
157  template <typename T_y, typename T_dof, typename T_loc, typename T_scale>
158  inline
159  typename boost::math::tools::promote_args<T_y,T_dof,T_loc,T_scale>::type
160  multi_student_t_log(const Eigen::Matrix<T_y,Eigen::Dynamic,1>& y,
161  const T_dof& nu,
162  const Eigen::Matrix<T_loc,Eigen::Dynamic,1>& mu,
163  const
164  Eigen::Matrix<T_scale,
165  Eigen::Dynamic,Eigen::Dynamic>& Sigma) {
166  return multi_student_t_log<false>(y,nu,mu,Sigma,
167  stan::math::default_policy());
168  }
169 
170 
171  template <class RNG>
172  inline Eigen::VectorXd
173  multi_student_t_rng(const double nu,
174  const Eigen::Matrix<double,Eigen::Dynamic,1>& mu,
175  const Eigen::Matrix<double,Eigen::Dynamic,Eigen::Dynamic>& s,
176  RNG& rng) {
177 
178  Eigen::VectorXd z(s.cols());
179  z.setZero();
180 
181  double w = stan::prob::inv_gamma_rng(nu / 2, nu / 2, rng);
182  return mu + std::sqrt(w) * stan::prob::multi_normal_rng(z, s, rng);
183  }
184  }
185 }
186 #endif
     [ Stan Home Page ] © 2011–2013, Stan Development Team.