1 #ifndef STAN_MATH_PRIM_MAT_PROB_MULTI_NORMAL_CHOLESKY_LPDF_HPP 2 #define STAN_MATH_PRIM_MAT_PROB_MULTI_NORMAL_CHOLESKY_LPDF_HPP 40 template <
bool propto,
typename T_y,
typename T_loc,
typename T_covar>
42 const T_y& y,
const T_loc& mu,
const T_covar& L) {
43 static const char*
function =
"multi_normal_cholesky_lpdf";
48 typedef Eigen::Matrix<T_partials_return, Eigen::Dynamic, Eigen::Dynamic>
50 typedef Eigen::Matrix<T_partials_return, Eigen::Dynamic, 1> vector_partials_t;
51 typedef Eigen::Matrix<T_partials_return, 1, Eigen::Dynamic>
52 row_vector_partials_t;
57 if (number_of_y == 0 || number_of_mu == 0)
63 const int size_y = y_vec[0].size();
64 const int size_mu = mu_vec[0].size();
65 if (
likely(size_vec > 1)) {
67 int size_y_old = size_y;
68 for (
size_t i = 1, size_ =
length_mvt(y); i < size_; i++) {
69 int size_y_new = y_vec[i].size();
71 "Size of one of the vectors of " 72 "the random variable",
74 "Size of another vector of the " 77 size_y_old = size_y_new;
80 int size_mu_old = size_mu;
81 for (
size_t i = 1, size_ =
length_mvt(mu); i < size_; i++) {
82 int size_mu_new = mu_vec[i].size();
84 "Size of one of the vectors of " 85 "the location variable",
87 "Size of another vector of the " 90 size_mu_old = size_mu_new;
95 "size of location parameter", size_mu);
97 "rows of covariance parameter", L.rows());
99 "columns of covariance parameter", L.cols());
101 for (
size_t i = 0; i < size_vec; i++) {
102 check_finite(
function,
"Location parameter", mu_vec[i]);
109 T_partials_return logp(0);
115 const matrix_partials_t inv_L_dbl
116 = mdivide_left_tri<Eigen::Lower>(
value_of(L));
119 for (
size_t i = 0; i < size_vec; i++) {
120 vector_partials_t y_minus_mu_dbl(size_y);
121 for (
int j = 0; j < size_y; j++)
124 const row_vector_partials_t half
125 = (inv_L_dbl.template triangularView<Eigen::Lower>() * y_minus_mu_dbl)
127 const vector_partials_t scaled_diff
128 = (half * inv_L_dbl.template triangularView<Eigen::Lower>())
134 for (
int j = 0; j < size_y; j++)
135 ops_partials.
edge1_.partials_vec_[i](j) -= scaled_diff(j);
138 for (
int j = 0; j < size_y; j++)
139 ops_partials.
edge2_.partials_vec_[i](j) += scaled_diff(j);
142 ops_partials.
edge3_.partials_ += scaled_diff * half;
148 logp += inv_L_dbl.diagonal().array().log().sum() * size_vec;
150 ops_partials.
edge3_.partials_ -= size_vec * inv_L_dbl.transpose();
154 return ops_partials.
build(logp);
157 template <
typename T_y,
typename T_loc,
typename T_covar>
160 return multi_normal_cholesky_lpdf<false>(y, mu, L);
return_type< T_y, T_loc, T_covar >::type multi_normal_cholesky_lpdf(const T_y &y, const T_loc &mu, const T_covar &L)
The log of the multivariate normal density for the given y, mu, and a Cholesky factor L of the varian...
void check_finite(const char *function, const char *name, const T_y &y)
Check if y is finite.
boost::math::tools::promote_args< double, typename partials_type< typename scalar_type< T >::type >::type, typename partials_return_type< T_pack... >::type >::type type
size_t max_size_mvt(const T1 &x1, const T2 &x2)
T value_of(const fvar< T > &v)
Return the value of the specified variable.
Extends std::true_type when instantiated with zero or more template parameters, all of which extend t...
This template builds partial derivatives with respect to a set of operands.
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.
fvar< T > dot_self(const Eigen::Matrix< fvar< T >, R, C > &v)
Template metaprogram to calculate whether a summand needs to be included in a proportional (log) prob...
boost::math::tools::promote_args< double, typename scalar_type< T >::type, typename return_type< Types_pack... >::type >::type type
void check_not_nan(const char *function, const char *name, const T_y &y)
Check if y is not NaN.
T_return_type build(double value)
Build the node to be stored on the autodiff graph.
const double NEG_LOG_SQRT_TWO_PI
internal::ops_partials_edge< double, Op2 > edge2_
This class provides a low-cost wrapper for situations where you either need an Eigen Vector or RowVec...
void check_consistent_sizes_mvt(const char *function, const char *name1, const T1 &x1, const char *name2, const T2 &x2)
Check if the dimension of x1 is consistent with x2.
size_t length_mvt(const Eigen::Matrix< T, R, C > &)
internal::ops_partials_edge< double, Op3 > edge3_
internal::ops_partials_edge< double, Op1 > edge1_