Automatic Differentiation
 
Loading...
Searching...
No Matches
wishart_lpdf.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_PROB_WISHART_LPDF_HPP
2#define STAN_MATH_PRIM_PROB_WISHART_LPDF_HPP
3
12
13namespace stan {
14namespace math {
15
46template <bool propto, typename T_y, typename T_dof, typename T_scale,
47 require_stan_scalar_t<T_dof>* = nullptr,
48 require_all_matrix_t<T_y, T_scale>* = nullptr>
50 const T_scale& S) {
51 using Eigen::Dynamic;
52 using Eigen::Lower;
53 using Eigen::Matrix;
54 using T_W_ref = ref_type_t<T_y>;
55 using T_nu_ref = ref_type_t<T_dof>;
56 using T_S_ref = ref_type_t<T_scale>;
57 static constexpr const char* function = "wishart_lpdf";
58 Eigen::Index k = W.rows();
59 check_size_match(function, "Rows of random variable", W.rows(),
60 "columns of scale parameter", S.rows());
61
62 T_W_ref W_ref = W;
63 T_nu_ref nu_ref = nu;
64 T_S_ref S_ref = S;
65
66 check_greater(function, "Degrees of freedom parameter", nu_ref, k - 1);
67 check_square(function, "random variable", W_ref);
68 check_square(function, "scale parameter", S_ref);
69 check_symmetric(function, "random variable", W_ref);
70 check_symmetric(function, "scale parameter", S_ref);
71
72 auto ldlt_W = make_ldlt_factor(W_ref);
73 check_ldlt_factor(function, "LDLT_Factor of random variable", ldlt_W);
74 auto ldlt_S = make_ldlt_factor(S_ref);
75 check_ldlt_factor(function, "LDLT_Factor of scale parameter", ldlt_S);
76
78
80 lp -= nu_ref * k * HALF_LOG_TWO;
81 }
82
84 lp -= lmgamma(k, 0.5 * nu_ref);
85 }
86
88 lp -= 0.5 * nu_ref * log_determinant_ldlt(ldlt_S);
89 }
90
92 lp -= 0.5 * trace(mdivide_left_ldlt(ldlt_S, W_ref));
93 }
94
95 if (include_summand<propto, T_y, T_dof>::value && nu_ref != (k + 1)) {
96 lp += 0.5 * (nu_ref - k - 1.0) * log_determinant_ldlt(ldlt_W);
97 }
98 return lp;
99}
100
101template <typename T_y, typename T_dof, typename T_scale>
103 const T_dof& nu,
104 const T_scale& S) {
105 return wishart_lpdf<false>(W, nu, S);
106}
107
108} // namespace math
109} // namespace stan
110#endif
void check_symmetric(const char *function, const char *name, const matrix_cl< T > &y)
Check if the matrix_cl is symmetric.
return_type_t< T_y, T_dof, T_scale > wishart_lpdf(const T_y &W, const T_dof &nu, const T_scale &S)
The log of the Wishart density for the given W, degrees of freedom, and scale matrix.
typename return_type< Ts... >::type return_type_t
Convenience type for the return type of the specified template parameters.
static constexpr double HALF_LOG_TWO
The value of half the natural logarithm 2, .
value_type_t< T > log_determinant_ldlt(LDLT_factor< T > &A)
Returns log(abs(det(A))) given a LDLT_factor of A.
auto make_ldlt_factor(const T &A)
Make an LDLT_factor with matrix type plain_type_t<T>
void check_square(const char *function, const char *name, const T_y &y)
Check if the specified matrix is square.
value_type_t< T > trace(const T &m)
Calculates trace (sum of diagonal) of given kernel generator expression.
Definition trace.hpp:22
Eigen::Matrix< value_type_t< EigMat >, Eigen::Dynamic, EigMat::ColsAtCompileTime > mdivide_left_ldlt(LDLT_factor< T > &A, const EigMat &b)
Returns the solution of the system Ax=b given an LDLT_factor of A.
fvar< return_type_t< T, int > > lmgamma(int x1, const fvar< T > &x2)
Definition lmgamma.hpp:14
void check_ldlt_factor(const char *function, const char *name, LDLT_factor< T > &A)
Raise domain error if the specified LDLT factor is invalid.
void check_size_match(const char *function, const char *name_i, T_size1 i, const char *name_j, T_size2 j)
Check if the provided sizes match.
void check_greater(const char *function, const char *name, const T_y &y, const T_low &low, Idxs... idxs)
Throw an exception if y is not strictly greater than low.
typename ref_type_if< true, T >::type ref_type_t
Definition ref_type.hpp:55
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...
Template metaprogram to calculate whether a summand needs to be included in a proportional (log) prob...