Automatic Differentiation
 
Loading...
Searching...
No Matches
multi_normal_cholesky_rng.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_PROB_MULTI_NORMAL_CHOLESKY_RNG_HPP
2#define STAN_MATH_PRIM_PROB_MULTI_NORMAL_CHOLESKY_RNG_HPP
3
10#include <boost/random/normal_distribution.hpp>
11#include <boost/random/variate_generator.hpp>
12
13namespace stan {
14namespace math {
15
32template <typename T_loc, class RNG>
35 const T_loc& mu,
36 const Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>& L, RNG& rng) {
37 using boost::normal_distribution;
38 using boost::variate_generator;
39
40 static constexpr const char* function = "multi_normal_cholesky_rng";
41 vector_seq_view<T_loc> mu_vec(mu);
42 size_t size_mu = mu_vec[0].size();
43
44 size_t N = size_mvt(mu);
45 for (size_t i = 1; i < N; i++) {
46 check_size_match(function,
47 "Size of one of the vectors of "
48 "the location variable",
49 mu_vec[i].size(),
50 "Size of the first vector of the "
51 "location variable",
52 size_mu);
53 }
54
55 for (size_t i = 0; i < N; i++) {
56 check_finite(function, "Location parameter", mu_vec[i]);
57 }
58 check_cholesky_factor(function, "Cholesky factor of covariance matrix", L);
59
60 const auto& L_ref = to_ref(L);
61
63
64 variate_generator<RNG&, normal_distribution<> > std_normal_rng(
65 rng, normal_distribution<>(0, 1));
66
67 for (size_t n = 0; n < N; ++n) {
68 Eigen::VectorXd z(L.cols());
69 for (int i = 0; i < L.cols(); i++) {
70 z(i) = std_normal_rng();
71 }
72
73 output[n] = as_column_vector_or_scalar(mu_vec[n]) + L_ref * z;
74 }
75
76 return output.data();
77}
78} // namespace math
79} // namespace stan
80#endif
typename helper::type type
StdVectorBuilder allocates type T1 values to be used as intermediate values.
This class provides a low-cost wrapper for situations where you either need an Eigen Vector or RowVec...
StdVectorBuilder< true, Eigen::VectorXd, T_loc >::type multi_normal_cholesky_rng(const T_loc &mu, const Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > &L, RNG &rng)
Return a multivariate normal random variate with the given location and Cholesky factorization of the...
auto as_column_vector_or_scalar(T &&a)
as_column_vector_or_scalar of a kernel generator expression.
double std_normal_rng(RNG &rng)
Return a standard Normal random variate using the specified random number generator.
int64_t size_mvt(const ScalarT &)
Provides the size of a multivariate argument.
Definition size_mvt.hpp:25
int64_t size(const T &m)
Returns the size (number of the elements) of a matrix_cl or var_value<matrix_cl<T>>.
Definition size.hpp:19
ref_type_t< T && > to_ref(T &&a)
This evaluates expensive Eigen expressions.
Definition to_ref.hpp:17
void check_finite(const char *function, const char *name, const T_y &y)
Return true if all values in y are finite.
void check_cholesky_factor(const char *function, const char *name, const Mat &y)
Throw an exception if the specified matrix is not a valid Cholesky factor.
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.
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...