Stan Math Library  2.20.0
reverse mode automatic differentiation
lkj_corr_cholesky_lpdf.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_PRIM_MAT_PROB_LKJ_CORR_CHOLESKY_LPDF_HPP
2 #define STAN_MATH_PRIM_MAT_PROB_LKJ_CORR_CHOLESKY_LPDF_HPP
3 
9 
10 namespace stan {
11 namespace math {
12 
13 // LKJ_Corr(L|eta) [ L Cholesky factor of correlation matrix
14 // eta > 0; eta == 1 <-> uniform]
15 template <bool propto, typename T_covar, typename T_shape>
16 typename boost::math::tools::promote_args<T_covar, T_shape>::type
18  const Eigen::Matrix<T_covar, Eigen::Dynamic, Eigen::Dynamic>& L,
19  const T_shape& eta) {
20  static const char* function = "lkj_corr_cholesky_lpdf";
21 
22  using boost::math::tools::promote_args;
23 
24  typedef typename promote_args<T_covar, T_shape>::type lp_ret;
25  lp_ret lp(0.0);
26  check_positive(function, "Shape parameter", eta);
27  check_lower_triangular(function, "Random variable", L);
28 
29  const unsigned int K = L.rows();
30  if (K == 0)
31  return 0.0;
32 
34  lp += do_lkj_constant(eta, K);
36  const int Km1 = K - 1;
37  Eigen::Matrix<T_covar, Eigen::Dynamic, 1> log_diagonals
38  = L.diagonal().tail(Km1).array().log();
39  Eigen::Matrix<lp_ret, Eigen::Dynamic, 1> values(Km1);
40  for (int k = 0; k < Km1; k++)
41  values(k) = (Km1 - k - 1) * log_diagonals(k);
42  if ((eta == 1.0)
44  lp += sum(values);
45  return (lp);
46  }
47  values += multiply(2.0 * eta - 2.0, log_diagonals);
48  lp += sum(values);
49  }
50 
51  return lp;
52 }
53 
54 template <typename T_covar, typename T_shape>
55 inline typename boost::math::tools::promote_args<T_covar, T_shape>::type
57  const Eigen::Matrix<T_covar, Eigen::Dynamic, Eigen::Dynamic>& L,
58  const T_shape& eta) {
59  return lkj_corr_cholesky_lpdf<false>(L, eta);
60 }
61 
62 } // namespace math
63 } // namespace stan
64 #endif
fvar< T > sum(const std::vector< fvar< T > > &m)
Return the sum of the entries of the specified standard vector.
Definition: sum.hpp:20
void check_lower_triangular(const char *function, const char *name, const Eigen::Matrix< T_y, Eigen::Dynamic, Eigen::Dynamic > &y)
Check if the specified matrix is lower triangular.
Extends std::true_type when instantiated with zero or more template parameters, all of which extend t...
Definition: conjunction.hpp:14
Eigen::Matrix< fvar< T >, R1, C1 > multiply(const Eigen::Matrix< fvar< T >, R1, C1 > &m, const fvar< T > &c)
Definition: multiply.hpp:14
Template metaprogram to calculate whether a summand needs to be included in a proportional (log) prob...
boost::math::tools::promote_args< T_covar, T_shape >::type lkj_corr_cholesky_lpdf(const Eigen::Matrix< T_covar, Eigen::Dynamic, Eigen::Dynamic > &L, const T_shape &eta)
void check_positive(const char *function, const char *name, const T_y &y)
Check if y is positive.
boost::math::tools::promote_args< double, T_shape >::type do_lkj_constant(const T_shape &eta, const unsigned int &K)

     [ Stan Home Page ] © 2011–2018, Stan Development Team.