Skip to content

Instantly share code, notes, and snippets.

@SteveBronder
Created March 12, 2020 18:21
Show Gist options
  • Save SteveBronder/191bc7af9da1384b42da923539031ab9 to your computer and use it in GitHub Desktop.
Save SteveBronder/191bc7af9da1384b42da923539031ab9 to your computer and use it in GitHub Desktop.
data {
int<lower = 1> nSubjects;
int nIIV;
}
parameters {
cholesky_factor_corr[nIIV] L;
matrix[nIIV, nSubjects] etaStd;
}
transformed parameters {
matrix[nIIV, nSubjects] thetaM;
thetaM = L * etaStd;
}
model { to_vector(etaStd) ~normal(0, 1); }
// Code generated by Stan version 2.22.0
#include <stan/model/model_header.hpp>
namespace debug2_model_namespace {
using stan::io::dump;
using stan::math::lgamma;
using stan::model::prob_grad;
using std::istream;
using std::string;
using std::stringstream;
using std::vector;
using namespace stan::math;
static int current_statement_begin__;
stan::io::program_reader prog_reader__() {
stan::io::program_reader reader;
reader.add_event(0, 0, "start", "examples/debug2.stan");
reader.add_event(20, 18, "end", "examples/debug2.stan");
return reader;
}
class debug2_model : public stan::model::model_base_crtp<debug2_model> {
private:
int nSubjects;
int nIIV;
public:
debug2_model(stan::io::var_context &context__, std::ostream *pstream__ = 0)
: model_base_crtp(0) {
ctor_body(context__, 0, pstream__);
}
debug2_model(stan::io::var_context &context__, unsigned int random_seed__,
std::ostream *pstream__ = 0)
: model_base_crtp(0) {
ctor_body(context__, random_seed__, pstream__);
}
void ctor_body(stan::io::var_context &context__, unsigned int random_seed__,
std::ostream *pstream__) {
typedef double local_scalar_t__;
boost::ecuyer1988 base_rng__ =
stan::services::util::create_rng(random_seed__, 0);
(void)base_rng__; // suppress unused var warning
current_statement_begin__ = -1;
static const char *function__ = "debug2_model_namespace::debug2_model";
(void)function__; // dummy to suppress unused var warning
size_t pos__;
(void)pos__; // dummy to suppress unused var warning
std::vector<int> vals_i__;
std::vector<double> vals_r__;
local_scalar_t__ DUMMY_VAR__(std::numeric_limits<double>::quiet_NaN());
(void)DUMMY_VAR__; // suppress unused var warning
try {
// initialize data block variables from context__
current_statement_begin__ = 2;
context__.validate_dims("data initialization", "nSubjects", "int",
context__.to_vec());
nSubjects = int(0);
vals_i__ = context__.vals_i("nSubjects");
pos__ = 0;
nSubjects = vals_i__[pos__++];
check_greater_or_equal(function__, "nSubjects", nSubjects, 1);
current_statement_begin__ = 3;
context__.validate_dims("data initialization", "nIIV", "int",
context__.to_vec());
nIIV = int(0);
vals_i__ = context__.vals_i("nIIV");
pos__ = 0;
nIIV = vals_i__[pos__++];
// initialize transformed data variables
// execute transformed data statements
// validate transformed data
// validate, set parameter ranges
num_params_r__ = 0U;
param_ranges_i__.clear();
current_statement_begin__ = 7;
validate_non_negative_index("L", "nIIV", nIIV);
validate_non_negative_index("L", "nIIV", nIIV);
num_params_r__ += ((nIIV * (nIIV - 1)) / 2);
current_statement_begin__ = 8;
validate_non_negative_index("etaStd", "nIIV", nIIV);
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
num_params_r__ += (nIIV * nSubjects);
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, current_statement_begin__,
prog_reader__());
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
}
~debug2_model() {}
void transform_inits(const stan::io::var_context &context__,
std::vector<int> &params_i__,
std::vector<double> &params_r__,
std::ostream *pstream__) const {
typedef double local_scalar_t__;
stan::io::writer<double> writer__(params_r__, params_i__);
size_t pos__;
(void)pos__; // dummy call to supress warning
std::vector<double> vals_r__;
std::vector<int> vals_i__;
current_statement_begin__ = 7;
if (!(context__.contains_r("L")))
stan::lang::rethrow_located(
std::runtime_error(std::string("Variable L missing")),
current_statement_begin__, prog_reader__());
vals_r__ = context__.vals_r("L");
pos__ = 0U;
validate_non_negative_index("L", "nIIV", nIIV);
validate_non_negative_index("L", "nIIV", nIIV);
context__.validate_dims("parameter initialization", "L", "matrix_d",
context__.to_vec(nIIV, nIIV));
Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L(nIIV, nIIV);
size_t L_j_2_max__ = nIIV;
size_t L_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < L_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < L_j_1_max__; ++j_1__) {
L(j_1__, j_2__) = vals_r__[pos__++];
}
}
try {
writer__.cholesky_factor_corr_unconstrain(L);
} catch (const std::exception &e) {
stan::lang::rethrow_located(
std::runtime_error(std::string("Error transforming variable L: ") +
e.what()),
current_statement_begin__, prog_reader__());
}
current_statement_begin__ = 8;
if (!(context__.contains_r("etaStd")))
stan::lang::rethrow_located(
std::runtime_error(std::string("Variable etaStd missing")),
current_statement_begin__, prog_reader__());
vals_r__ = context__.vals_r("etaStd");
pos__ = 0U;
validate_non_negative_index("etaStd", "nIIV", nIIV);
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
context__.validate_dims("parameter initialization", "etaStd", "matrix_d",
context__.to_vec(nIIV, nSubjects));
Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> etaStd(nIIV,
nSubjects);
size_t etaStd_j_2_max__ = nSubjects;
size_t etaStd_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < etaStd_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < etaStd_j_1_max__; ++j_1__) {
etaStd(j_1__, j_2__) = vals_r__[pos__++];
}
}
try {
writer__.matrix_unconstrain(etaStd);
} catch (const std::exception &e) {
stan::lang::rethrow_located(
std::runtime_error(
std::string("Error transforming variable etaStd: ") + e.what()),
current_statement_begin__, prog_reader__());
}
params_r__ = writer__.data_r();
params_i__ = writer__.data_i();
}
void transform_inits(const stan::io::var_context &context,
Eigen::Matrix<double, Eigen::Dynamic, 1> &params_r,
std::ostream *pstream__) const {
std::vector<double> params_r_vec;
std::vector<int> params_i_vec;
transform_inits(context, params_i_vec, params_r_vec, pstream__);
params_r.resize(params_r_vec.size());
for (int i = 0; i < params_r.size(); ++i)
params_r(i) = params_r_vec[i];
}
template <bool propto__, bool jacobian__, typename T__>
T__ log_prob(std::vector<T__> &params_r__, std::vector<int> &params_i__,
std::ostream *pstream__ = 0) const {
typedef T__ local_scalar_t__;
local_scalar_t__ DUMMY_VAR__(std::numeric_limits<double>::quiet_NaN());
(void)DUMMY_VAR__; // dummy to suppress unused var warning
T__ lp__(0.0);
stan::math::accumulator<T__> lp_accum__;
try {
stan::io::reader<local_scalar_t__> in__(params_r__, params_i__);
// model parameters
current_statement_begin__ = 7;
Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> L;
(void)L; // dummy to suppress unused var warning
if (jacobian__)
L = in__.cholesky_factor_corr_constrain(nIIV, lp__);
else
L = in__.cholesky_factor_corr_constrain(nIIV);
current_statement_begin__ = 8;
Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> etaStd;
(void)etaStd; // dummy to suppress unused var warning
if (jacobian__)
etaStd = in__.matrix_constrain(nIIV, nSubjects, lp__);
else
etaStd = in__.matrix_constrain(nIIV, nSubjects);
// transformed parameters
current_statement_begin__ = 12;
validate_non_negative_index("thetaM", "nIIV", nIIV);
validate_non_negative_index("thetaM", "nSubjects", nSubjects);
Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> thetaM(
nIIV, nSubjects);
stan::math::initialize(thetaM, DUMMY_VAR__);
stan::math::fill(thetaM, DUMMY_VAR__);
// transformed parameters block statements
current_statement_begin__ = 13;
stan::math::assign(thetaM, multiply(L, etaStd));
// validate transformed parameters
const char *function__ = "validate transformed params";
(void)function__; // dummy to suppress unused var warning
current_statement_begin__ = 12;
size_t thetaM_j_1_max__ = nIIV;
size_t thetaM_j_2_max__ = nSubjects;
for (size_t j_1__ = 0; j_1__ < thetaM_j_1_max__; ++j_1__) {
for (size_t j_2__ = 0; j_2__ < thetaM_j_2_max__; ++j_2__) {
if (stan::math::is_uninitialized(thetaM(j_1__, j_2__))) {
std::stringstream msg__;
msg__ << "Undefined transformed parameter: thetaM"
<< "(" << j_1__ << ", " << j_2__ << ")";
stan::lang::rethrow_located(
std::runtime_error(
std::string("Error initializing variable thetaM: ") +
msg__.str()),
current_statement_begin__, prog_reader__());
}
}
}
// model body
current_statement_begin__ = 17;
lp_accum__.add(normal_log<propto__>(to_vector(etaStd), 0, 1));
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, current_statement_begin__,
prog_reader__());
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
lp_accum__.add(lp__);
return lp_accum__.sum();
} // log_prob()
template <bool propto, bool jacobian, typename T_>
T_ log_prob(Eigen::Matrix<T_, Eigen::Dynamic, 1> &params_r,
std::ostream *pstream = 0) const {
std::vector<T_> vec_params_r;
vec_params_r.reserve(params_r.size());
for (int i = 0; i < params_r.size(); ++i)
vec_params_r.push_back(params_r(i));
std::vector<int> vec_params_i;
return log_prob<propto, jacobian, T_>(vec_params_r, vec_params_i, pstream);
}
void get_param_names(std::vector<std::string> &names__) const {
names__.resize(0);
names__.push_back("L");
names__.push_back("etaStd");
names__.push_back("thetaM");
}
void get_dims(std::vector<std::vector<size_t>> &dimss__) const {
dimss__.resize(0);
std::vector<size_t> dims__;
dims__.resize(0);
dims__.push_back(nIIV);
dims__.push_back(nIIV);
dimss__.push_back(dims__);
dims__.resize(0);
dims__.push_back(nIIV);
dims__.push_back(nSubjects);
dimss__.push_back(dims__);
dims__.resize(0);
dims__.push_back(nIIV);
dims__.push_back(nSubjects);
dimss__.push_back(dims__);
}
template <typename RNG>
void write_array(RNG &base_rng__, std::vector<double> &params_r__,
std::vector<int> &params_i__, std::vector<double> &vars__,
bool include_tparams__ = true, bool include_gqs__ = true,
std::ostream *pstream__ = 0) const {
typedef double local_scalar_t__;
vars__.resize(0);
stan::io::reader<local_scalar_t__> in__(params_r__, params_i__);
static const char *function__ = "debug2_model_namespace::write_array";
(void)function__; // dummy to suppress unused var warning
// read-transform, write parameters
Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L =
in__.cholesky_factor_corr_constrain(nIIV);
size_t L_j_2_max__ = nIIV;
size_t L_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < L_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < L_j_1_max__; ++j_1__) {
vars__.push_back(L(j_1__, j_2__));
}
}
Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> etaStd =
in__.matrix_constrain(nIIV, nSubjects);
size_t etaStd_j_2_max__ = nSubjects;
size_t etaStd_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < etaStd_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < etaStd_j_1_max__; ++j_1__) {
vars__.push_back(etaStd(j_1__, j_2__));
}
}
double lp__ = 0.0;
(void)lp__; // dummy to suppress unused var warning
stan::math::accumulator<double> lp_accum__;
local_scalar_t__ DUMMY_VAR__(std::numeric_limits<double>::quiet_NaN());
(void)DUMMY_VAR__; // suppress unused var warning
if (!include_tparams__ && !include_gqs__)
return;
try {
// declare and define transformed parameters
current_statement_begin__ = 12;
validate_non_negative_index("thetaM", "nIIV", nIIV);
validate_non_negative_index("thetaM", "nSubjects", nSubjects);
Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> thetaM(nIIV,
nSubjects);
stan::math::initialize(thetaM, DUMMY_VAR__);
stan::math::fill(thetaM, DUMMY_VAR__);
// do transformed parameters statements
current_statement_begin__ = 13;
stan::math::assign(thetaM, multiply(L, etaStd));
if (!include_gqs__ && !include_tparams__)
return;
// validate transformed parameters
const char *function__ = "validate transformed params";
(void)function__; // dummy to suppress unused var warning
// write transformed parameters
if (include_tparams__) {
size_t thetaM_j_2_max__ = nSubjects;
size_t thetaM_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < thetaM_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < thetaM_j_1_max__; ++j_1__) {
vars__.push_back(thetaM(j_1__, j_2__));
}
}
}
if (!include_gqs__)
return;
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, current_statement_begin__,
prog_reader__());
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
}
template <typename RNG>
void write_array(RNG &base_rng,
Eigen::Matrix<double, Eigen::Dynamic, 1> &params_r,
Eigen::Matrix<double, Eigen::Dynamic, 1> &vars,
bool include_tparams = true, bool include_gqs = true,
std::ostream *pstream = 0) const {
std::vector<double> params_r_vec(params_r.size());
for (int i = 0; i < params_r.size(); ++i)
params_r_vec[i] = params_r(i);
std::vector<double> vars_vec;
std::vector<int> params_i_vec;
write_array(base_rng, params_r_vec, params_i_vec, vars_vec, include_tparams,
include_gqs, pstream);
vars.resize(vars_vec.size());
for (int i = 0; i < vars.size(); ++i)
vars(i) = vars_vec[i];
}
std::string model_name() const { return "debug2_model"; }
void constrained_param_names(std::vector<std::string> &param_names__,
bool include_tparams__ = true,
bool include_gqs__ = true) const {
std::stringstream param_name_stream__;
size_t L_j_2_max__ = nIIV;
size_t L_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < L_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < L_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "L" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
param_names__.push_back(param_name_stream__.str());
}
}
size_t etaStd_j_2_max__ = nSubjects;
size_t etaStd_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < etaStd_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < etaStd_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "etaStd" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
param_names__.push_back(param_name_stream__.str());
}
}
if (!include_gqs__ && !include_tparams__)
return;
if (include_tparams__) {
size_t thetaM_j_2_max__ = nSubjects;
size_t thetaM_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < thetaM_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < thetaM_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "thetaM" << '.' << j_1__ + 1 << '.'
<< j_2__ + 1;
param_names__.push_back(param_name_stream__.str());
}
}
}
if (!include_gqs__)
return;
}
void unconstrained_param_names(std::vector<std::string> &param_names__,
bool include_tparams__ = true,
bool include_gqs__ = true) const {
std::stringstream param_name_stream__;
size_t L_j_1_max__ = ((nIIV * (nIIV - 1)) / 2);
for (size_t j_1__ = 0; j_1__ < L_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "L" << '.' << j_1__ + 1;
param_names__.push_back(param_name_stream__.str());
}
size_t etaStd_j_2_max__ = nSubjects;
size_t etaStd_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < etaStd_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < etaStd_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "etaStd" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
param_names__.push_back(param_name_stream__.str());
}
}
if (!include_gqs__ && !include_tparams__)
return;
if (include_tparams__) {
size_t thetaM_j_2_max__ = nSubjects;
size_t thetaM_j_1_max__ = nIIV;
for (size_t j_2__ = 0; j_2__ < thetaM_j_2_max__; ++j_2__) {
for (size_t j_1__ = 0; j_1__ < thetaM_j_1_max__; ++j_1__) {
param_name_stream__.str(std::string());
param_name_stream__ << "thetaM" << '.' << j_1__ + 1 << '.'
<< j_2__ + 1;
param_names__.push_back(param_name_stream__.str());
}
}
}
if (!include_gqs__)
return;
}
}; // model
} // namespace debug2_model_namespace
typedef debug2_model_namespace::debug2_model stan_model;
#ifndef USING_R
stan::model::model_base &new_model(stan::io::var_context &data_context,
unsigned int seed,
std::ostream *msg_stream) {
stan_model *m = new stan_model(data_context, seed, msg_stream);
return *m;
}
#endif
// Code generated by %%NAME%% %%VERSION%%
#include <stan/model/model_header.hpp>
namespace debug_model_namespace {
template <typename T, typename S>
std::vector<T> resize_to_match__(std::vector<T> &dst,
const std::vector<S> &src) {
dst.resize(src.size());
return dst;
}
template <typename T>
Eigen::Matrix<T, -1, -1>
resize_to_match__(Eigen::Matrix<T, -1, -1> &dst,
const Eigen::Matrix<T, -1, -1> &src) {
dst.resize(src.rows(), src.cols());
return dst;
}
template <typename T>
Eigen::Matrix<T, 1, -1> resize_to_match__(Eigen::Matrix<T, 1, -1> &dst,
const Eigen::Matrix<T, 1, -1> &src) {
dst.resize(src.size());
return dst;
}
template <typename T>
Eigen::Matrix<T, -1, 1> resize_to_match__(Eigen::Matrix<T, -1, 1> &dst,
const Eigen::Matrix<T, -1, 1> &src) {
dst.resize(src.size());
return dst;
}
std::vector<double> to_doubles__(std::initializer_list<double> x) { return x; }
std::vector<stan::math::var>
to_vars__(std::initializer_list<stan::math::var> x) {
return x;
}
inline void validate_positive_index(const char *var_name, const char *expr,
int val) {
if (val < 1) {
std::stringstream msg;
msg << "Found dimension size less than one in simplex declaration"
<< "; variable=" << var_name << "; dimension size expression=" << expr
<< "; expression value=" << val;
std::string msg_str(msg.str());
throw std::invalid_argument(msg_str.c_str());
}
}
inline void validate_unit_vector_index(const char *var_name, const char *expr,
int val) {
if (val <= 1) {
std::stringstream msg;
if (val == 1) {
msg << "Found dimension size one in unit vector declaration."
<< " One-dimensional unit vector is discrete"
<< " but the target distribution must be continuous."
<< " variable=" << var_name << "; dimension size expression=" << expr;
} else {
msg << "Found dimension size less than one in unit vector declaration"
<< "; variable=" << var_name << "; dimension size expression=" << expr
<< "; expression value=" << val;
}
std::string msg_str(msg.str());
throw std::invalid_argument(msg_str.c_str());
}
}
using stan::io::dump;
using stan::math::lgamma;
using stan::model::cons_list;
using stan::model::index_max;
using stan::model::index_min;
using stan::model::index_min_max;
using stan::model::index_multi;
using stan::model::index_omni;
using stan::model::index_uni;
using stan::model::model_base_crtp;
using stan::model::nil_index_list;
using stan::model::rvalue;
using std::istream;
using std::string;
using std::stringstream;
using std::vector;
using namespace stan::math;
static int current_statement__ = 0;
static const std::vector<string> locations_array__ = {
" (found before start of program)",
" (in 'examples/debug.stan', line 7, column 2 to column 31)",
" (in 'examples/debug.stan', line 8, column 2 to column 33)",
" (in 'examples/debug.stan', line 12, column 2 to column 33)",
" (in 'examples/debug.stan', line 13, column 2 to column 22)",
" (in 'examples/debug.stan', line 17, column 2 to column 35)",
" (in 'examples/debug.stan', line 2, column 2 to column 27)",
" (in 'examples/debug.stan', line 3, column 2 to column 11)"};
class debug_model : public model_base_crtp<debug_model> {
private:
int pos__;
int nSubjects;
int nIIV;
public:
~debug_model() {}
std::string model_name() const { return "debug_model"; }
debug_model(stan::io::var_context &context__, unsigned int random_seed__ = 0,
std::ostream *pstream__ = nullptr)
: model_base_crtp(0) {
typedef double local_scalar_t__;
boost::ecuyer1988 base_rng__ =
stan::services::util::create_rng(random_seed__, 0);
(void)base_rng__; // suppress unused var warning
static const char *function__ = "debug_model_namespace::debug_model";
(void)function__; // suppress unused var warning
try {
pos__ = 1;
context__.validate_dims("data initialization", "nSubjects", "int",
context__.to_vec());
current_statement__ = 6;
nSubjects = context__.vals_i("nSubjects")[(1 - 1)];
context__.validate_dims("data initialization", "nIIV", "int",
context__.to_vec());
current_statement__ = 7;
nIIV = context__.vals_i("nIIV")[(1 - 1)];
current_statement__ = 6;
current_statement__ = 6;
check_greater_or_equal(function__, "nSubjects", nSubjects, 1);
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, locations_array__[current_statement__]);
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
num_params_r__ = 0U;
try {
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
num_params_r__ += ((nIIV * (nIIV - 1)) / 2);
current_statement__ = 2;
validate_non_negative_index("etaStd", "nIIV", nIIV);
current_statement__ = 2;
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
num_params_r__ += nIIV * nSubjects;
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, locations_array__[current_statement__]);
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
}
template <bool propto__, bool jacobian__, typename T__>
T__ log_prob(std::vector<T__> &params_r__, std::vector<int> &params_i__,
std::ostream *pstream__ = 0) const {
typedef T__ local_scalar_t__;
T__ lp__(0.0);
stan::math::accumulator<T__> lp_accum__;
static const char *function__ = "debug_model_namespace::log_prob";
(void)function__; // suppress unused var warning
stan::io::reader<local_scalar_t__> in__(params_r__, params_i__);
try {
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
Eigen::Matrix<local_scalar_t__, -1, -1> L;
L = Eigen::Matrix<local_scalar_t__, -1, -1>(nIIV, nIIV);
Eigen::Matrix<local_scalar_t__, -1, 1> L_in__;
L_in__ =
Eigen::Matrix<local_scalar_t__, -1, 1>(((nIIV * (nIIV - 1)) / 2));
current_statement__ = 1;
L_in__ = in__.vector(((nIIV * (nIIV - 1)) / 2));
current_statement__ = 1;
if (jacobian__) {
current_statement__ = 1;
assign(L, nil_index_list(),
stan::math::cholesky_corr_constrain(L_in__, nIIV, lp__),
"assigning variable L");
} else {
current_statement__ = 1;
assign(L, nil_index_list(),
stan::math::cholesky_corr_constrain(L_in__, nIIV),
"assigning variable L");
}
current_statement__ = 2;
validate_non_negative_index("etaStd", "nIIV", nIIV);
current_statement__ = 2;
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
Eigen::Matrix<local_scalar_t__, -1, -1> etaStd;
etaStd = Eigen::Matrix<local_scalar_t__, -1, -1>(nIIV, nSubjects);
current_statement__ = 2;
etaStd = in__.matrix(nIIV, nSubjects);
current_statement__ = 3;
validate_non_negative_index("thetaM", "nIIV", nIIV);
current_statement__ = 3;
validate_non_negative_index("thetaM", "nSubjects", nSubjects);
Eigen::Matrix<local_scalar_t__, -1, -1> thetaM;
thetaM = Eigen::Matrix<local_scalar_t__, -1, -1>(nIIV, nSubjects);
current_statement__ = 3;
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
current_statement__ = 3;
for (size_t sym2__ = 1; sym2__ <= nSubjects; ++sym2__) {
current_statement__ = 3;
assign(thetaM,
cons_list(index_uni(sym1__),
cons_list(index_uni(sym2__), nil_index_list())),
std::numeric_limits<double>::quiet_NaN(),
"assigning variable thetaM");
}
}
current_statement__ = 4;
assign(thetaM, nil_index_list(), multiply(L, etaStd),
"assigning variable thetaM");
{
current_statement__ = 5;
lp_accum__.add(normal_log<propto__>(to_vector(etaStd), 0, 1));
}
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, locations_array__[current_statement__]);
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
lp_accum__.add(lp__);
return lp_accum__.sum();
} // log_prob()
template <typename RNG>
void write_array(RNG &base_rng__, std::vector<double> &params_r__,
std::vector<int> &params_i__, std::vector<double> &vars__,
bool emit_transformed_parameters__ = true,
bool emit_generated_quantities__ = true,
std::ostream *pstream__ = 0) const {
typedef double local_scalar_t__;
vars__.resize(0);
stan::io::reader<local_scalar_t__> in__(params_r__, params_i__);
static const char *function__ = "debug_model_namespace::write_array";
(void)function__; // suppress unused var warning
(void)function__; // suppress unused var warning
double lp__ = 0.0;
(void)lp__; // dummy to suppress unused var warning
stan::math::accumulator<double> lp_accum__;
try {
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
Eigen::Matrix<double, -1, -1> L;
L = Eigen::Matrix<double, -1, -1>(nIIV, nIIV);
Eigen::Matrix<local_scalar_t__, -1, 1> L_in__;
L_in__ =
Eigen::Matrix<local_scalar_t__, -1, 1>(((nIIV * (nIIV - 1)) / 2));
current_statement__ = 1;
L_in__ = in__.vector(((nIIV * (nIIV - 1)) / 2));
current_statement__ = 1;
assign(L, nil_index_list(),
stan::math::cholesky_corr_constrain(L_in__, nIIV),
"assigning variable L");
current_statement__ = 2;
validate_non_negative_index("etaStd", "nIIV", nIIV);
current_statement__ = 2;
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
Eigen::Matrix<double, -1, -1> etaStd;
etaStd = Eigen::Matrix<double, -1, -1>(nIIV, nSubjects);
current_statement__ = 2;
etaStd = in__.matrix(nIIV, nSubjects);
current_statement__ = 3;
validate_non_negative_index("thetaM", "nIIV", nIIV);
current_statement__ = 3;
validate_non_negative_index("thetaM", "nSubjects", nSubjects);
Eigen::Matrix<double, -1, -1> thetaM;
thetaM = Eigen::Matrix<double, -1, -1>(nIIV, nSubjects);
current_statement__ = 3;
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
current_statement__ = 3;
for (size_t sym2__ = 1; sym2__ <= nSubjects; ++sym2__) {
current_statement__ = 3;
assign(thetaM,
cons_list(index_uni(sym1__),
cons_list(index_uni(sym2__), nil_index_list())),
std::numeric_limits<double>::quiet_NaN(),
"assigning variable thetaM");
}
}
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
vars__.push_back(
rvalue(L,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
"L"));
}
}
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
vars__.push_back(
rvalue(etaStd,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
"etaStd"));
}
}
if (logical_negation((primitive_value(emit_transformed_parameters__) ||
primitive_value(emit_generated_quantities__)))) {
return;
}
current_statement__ = 4;
assign(thetaM, nil_index_list(), multiply(L, etaStd),
"assigning variable thetaM");
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
vars__.push_back(
rvalue(thetaM,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
"thetaM"));
}
}
if (logical_negation(emit_generated_quantities__)) {
return;
}
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, locations_array__[current_statement__]);
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
} // write_array()
void transform_inits(const stan::io::var_context &context__,
std::vector<int> &params_i__,
std::vector<double> &vars__,
std::ostream *pstream__) const {
typedef double local_scalar_t__;
vars__.resize(0);
vars__.reserve(num_params_r__);
try {
int pos__;
pos__ = 1;
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
current_statement__ = 1;
validate_non_negative_index("L", "nIIV", nIIV);
Eigen::Matrix<double, -1, -1> L;
L = Eigen::Matrix<double, -1, -1>(nIIV, nIIV);
{
std::vector<local_scalar_t__> L_flat__;
current_statement__ = 1;
assign(L_flat__, nil_index_list(), context__.vals_r("L"),
"assigning variable L_flat__");
current_statement__ = 1;
pos__ = 1;
current_statement__ = 1;
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
current_statement__ = 1;
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
current_statement__ = 1;
assign(L,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
L_flat__[(pos__ - 1)], "assigning variable L");
current_statement__ = 1;
pos__ = (pos__ + 1);
}
}
}
current_statement__ = 1;
assign(L, nil_index_list(), stan::math::cholesky_corr_free(L),
"assigning variable L");
current_statement__ = 2;
validate_non_negative_index("etaStd", "nIIV", nIIV);
current_statement__ = 2;
validate_non_negative_index("etaStd", "nSubjects", nSubjects);
Eigen::Matrix<double, -1, -1> etaStd;
etaStd = Eigen::Matrix<double, -1, -1>(nIIV, nSubjects);
{
std::vector<local_scalar_t__> etaStd_flat__;
current_statement__ = 2;
assign(etaStd_flat__, nil_index_list(), context__.vals_r("etaStd"),
"assigning variable etaStd_flat__");
current_statement__ = 2;
pos__ = 1;
current_statement__ = 2;
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
current_statement__ = 2;
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
current_statement__ = 2;
assign(etaStd,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
etaStd_flat__[(pos__ - 1)], "assigning variable etaStd");
current_statement__ = 2;
pos__ = (pos__ + 1);
}
}
}
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
vars__.push_back(
rvalue(L,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
"L"));
}
}
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
vars__.push_back(
rvalue(etaStd,
cons_list(index_uni(sym2__),
cons_list(index_uni(sym1__), nil_index_list())),
"etaStd"));
}
}
} catch (const std::exception &e) {
stan::lang::rethrow_located(e, locations_array__[current_statement__]);
// Next line prevents compiler griping about no return
throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***");
}
} // transform_inits()
void get_param_names(std::vector<std::string> &names__) const {
names__.resize(0);
names__.push_back("L");
names__.push_back("etaStd");
names__.push_back("thetaM");
} // get_param_names()
void get_dims(std::vector<std::vector<size_t>> &dimss__) const {
dimss__.resize(0);
std::vector<size_t> dims__;
dims__.push_back(nIIV);
dims__.push_back(nIIV);
dimss__.push_back(dims__);
dims__.resize(0);
dims__.push_back(nIIV);
dims__.push_back(nSubjects);
dimss__.push_back(dims__);
dims__.resize(0);
dims__.push_back(nIIV);
dims__.push_back(nSubjects);
dimss__.push_back(dims__);
dims__.resize(0);
} // get_dims()
void constrained_param_names(std::vector<std::string> &param_names__,
bool emit_transformed_parameters__ = true,
bool emit_generated_quantities__ = true) const {
for (size_t sym1__ = 1; sym1__ <= nIIV; ++sym1__) {
{
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
{
param_names__.push_back(std::string() + "L" + '.' +
std::to_string(sym2__) + '.' +
std::to_string(sym1__));
}
}
}
}
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
{
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
{
param_names__.push_back(std::string() + "etaStd" + '.' +
std::to_string(sym2__) + '.' +
std::to_string(sym1__));
}
}
}
}
if (emit_transformed_parameters__) {
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
{
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
{
param_names__.push_back(std::string() + "thetaM" + '.' +
std::to_string(sym2__) + '.' +
std::to_string(sym1__));
}
}
}
}
}
if (emit_generated_quantities__) {
}
} // constrained_param_names()
void
unconstrained_param_names(std::vector<std::string> &param_names__,
bool emit_transformed_parameters__ = true,
bool emit_generated_quantities__ = true) const {
for (size_t sym1__ = 1; sym1__ <= ((nIIV * (nIIV - 1)) / 2); ++sym1__) {
{
param_names__.push_back(std::string() + "L" + '.' +
std::to_string(sym1__));
}
}
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
{
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
{
param_names__.push_back(std::string() + "etaStd" + '.' +
std::to_string(sym2__) + '.' +
std::to_string(sym1__));
}
}
}
}
if (emit_transformed_parameters__) {
for (size_t sym1__ = 1; sym1__ <= nSubjects; ++sym1__) {
{
for (size_t sym2__ = 1; sym2__ <= nIIV; ++sym2__) {
{
param_names__.push_back(std::string() + "thetaM" + '.' +
std::to_string(sym2__) + '.' +
std::to_string(sym1__));
}
}
}
}
}
if (emit_generated_quantities__) {
}
} // unconstrained_param_names()
std::string get_constrained_sizedtypes() const {
stringstream s__;
s__ << "[{\"name\":\"L\",\"type\":{\"name\":\"matrix\",\"rows\":" << nIIV
<< ",\"cols\":" << nIIV
<< "},\"block\":\"parameters\"},{\"name\":\"etaStd\",\"type\":{"
"\"name\":\"matrix\",\"rows\":"
<< nIIV << ",\"cols\":" << nSubjects
<< "},\"block\":\"parameters\"},{\"name\":\"thetaM\",\"type\":{"
"\"name\":\"matrix\",\"rows\":"
<< nIIV << ",\"cols\":" << nSubjects
<< "},\"block\":\"transformed_parameters\"}]";
return s__.str();
} // get_constrained_sizedtypes()
std::string get_unconstrained_sizedtypes() const {
stringstream s__;
s__ << "[{\"name\":\"L\",\"type\":{\"name\":\"vector\",\"length\":"
<< ((nIIV * (nIIV - 1)) / 2)
<< "},\"block\":\"parameters\"},{\"name\":\"etaStd\",\"type\":{"
"\"name\":\"matrix\",\"rows\":"
<< nIIV << ",\"cols\":" << nSubjects
<< "},\"block\":\"parameters\"},{\"name\":\"thetaM\",\"type\":{"
"\"name\":\"matrix\",\"rows\":"
<< nIIV << ",\"cols\":" << nSubjects
<< "},\"block\":\"transformed_parameters\"}]";
return s__.str();
} // get_unconstrained_sizedtypes()
// Begin method overload boilerplate
template <typename RNG>
void write_array(RNG &base_rng__,
Eigen::Matrix<double, Eigen::Dynamic, 1> &params_r,
Eigen::Matrix<double, Eigen::Dynamic, 1> &vars,
bool emit_transformed_parameters__ = true,
bool emit_generated_quantities__ = true,
std::ostream *pstream = 0) const {
std::vector<double> params_r_vec(params_r.size());
for (int i = 0; i < params_r.size(); ++i)
params_r_vec[i] = params_r(i);
std::vector<double> vars_vec;
std::vector<int> params_i_vec;
write_array(base_rng__, params_r_vec, params_i_vec, vars_vec,
emit_transformed_parameters__, emit_generated_quantities__,
pstream);
vars.resize(vars_vec.size());
for (int i = 0; i < vars.size(); ++i)
vars(i) = vars_vec[i];
}
template <bool propto__, bool jacobian__, typename T_>
T_ log_prob(Eigen::Matrix<T_, Eigen::Dynamic, 1> &params_r,
std::ostream *pstream = 0) const {
std::vector<T_> vec_params_r;
vec_params_r.reserve(params_r.size());
for (int i = 0; i < params_r.size(); ++i)
vec_params_r.push_back(params_r(i));
std::vector<int> vec_params_i;
return log_prob<propto__, jacobian__, T_>(vec_params_r, vec_params_i,
pstream);
}
void transform_inits(const stan::io::var_context &context,
Eigen::Matrix<double, Eigen::Dynamic, 1> &params_r,
std::ostream *pstream__) const {
std::vector<double> params_r_vec;
std::vector<int> params_i_vec;
transform_inits(context, params_i_vec, params_r_vec, pstream__);
params_r.resize(params_r_vec.size());
for (int i = 0; i < params_r.size(); ++i)
params_r(i) = params_r_vec[i];
}
};
} // namespace debug_model_namespace
typedef debug_model_namespace::debug_model stan_model;
#ifndef USING_R
// Boilerplate
stan::model::model_base &new_model(stan::io::var_context &data_context,
unsigned int seed,
std::ostream *msg_stream) {
stan_model *m = new stan_model(data_context, seed, msg_stream);
return *m;
}
#endif
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment