1#ifndef STAN_MATH_PRIM_CONSTRAINT_STOCHASTIC_ROW_CONSTRAIN_HPP
2#define STAN_MATH_PRIM_CONSTRAINT_STOCHASTIC_ROW_CONSTRAIN_HPP
25template <
typename Mat, require_eigen_matrix_dynamic_t<Mat>* =
nullptr,
26 require_not_st_var<Mat>* =
nullptr>
29 const Eigen::Index N = y_ref.rows();
30 int Km1 = y_ref.cols();
32 using eigen_arr = Eigen::Array<scalar_type_t<Mat>, -1, 1>;
33 eigen_arr stick_len = eigen_arr::Constant(N, 1.0);
34 for (Eigen::Index k = 0; k < Km1; ++k) {
35 auto z_k =
inv_logit(y_ref.array().col(k) -
log(Km1 - k));
36 x.array().col(k) = stick_len * z_k;
37 stick_len -= x.array().col(k);
39 x.array().col(Km1) = stick_len;
55template <
typename Mat,
typename Lp,
61 const Eigen::Index N = y_ref.rows();
62 Eigen::Index Km1 = y_ref.cols();
64 Eigen::Array<scalar_type_t<Mat>, -1, 1> stick_len
65 = Eigen::Array<scalar_type_t<Mat>, -1, 1>::Constant(N, 1.0);
66 for (Eigen::Index k = 0; k < Km1; ++k) {
67 const auto eq_share = -
log(Km1 - k);
68 auto adj_y_k = (y_ref.array().col(k) + eq_share).
eval();
70 x.array().col(k) = stick_len * z_k;
73 stick_len -= x.array().col(k);
75 x.col(Km1).array() = stick_len;
89template <
typename T, require_std_vector_t<T>* =
nullptr>
108template <
typename T,
typename Lp, require_std_vector_t<T>* =
nullptr,
109 require_convertible_t<return_type_t<T>, Lp>* =
nullptr>
133template <
bool Jacobian,
typename Mat,
typename Lp,
136 if constexpr (Jacobian) {
require_t< std::is_convertible< std::decay_t< T >, std::decay_t< S > > > require_convertible_t
Require types T and S satisfies std::is_convertible.
require_t< is_eigen_matrix_dynamic< std::decay_t< T > > > require_eigen_matrix_dynamic_t
Require type satisfies is_eigen_matrix_dynamic.
require_not_t< is_var< scalar_type_t< std::decay_t< T > > > > require_not_st_var
Require scalar type does not satisfy is_var.
plain_type_t< Mat > stochastic_row_constrain(const Mat &y)
Return a row stochastic matrix.
T eval(T &&arg)
Inputs which have a plain_type equal to the own time are forwarded unmodified (for Eigen expressions ...
fvar< T > log(const fvar< T > &x)
fvar< T > log1p_exp(const fvar< T > &x)
ref_type_t< T && > to_ref(T &&a)
This evaluates expensive Eigen expressions.
auto sum(const std::vector< T > &m)
Return the sum of the entries of the specified standard vector.
fvar< T > inv_logit(const fvar< T > &x)
Returns the inverse logit function applied to the argument.
typename plain_type< T >::type plain_type_t
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...