// Code generated by Stan version 2.21.0

#include <stan/model/model_header.hpp>

namespace model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb_namespace {

using std::istream;
using std::string;
using std::stringstream;
using std::vector;
using stan::io::dump;
using stan::math::lgamma;
using stan::model::prob_grad;
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", "model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb");
    reader.add_event(300, 298, "end", "model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb");
    return reader;
}

template <typename T0__, typename T1__, typename T2__, typename T3__>
Eigen::Matrix<typename boost::math::tools::promote_args<T0__, T1__, T2__, T3__>::type, Eigen::Dynamic, Eigen::Dynamic>
scale_r_cor_cov(const Eigen::Matrix<T0__, Eigen::Dynamic, Eigen::Dynamic>& z,
                    const Eigen::Matrix<T1__, Eigen::Dynamic, 1>& SD,
                    const Eigen::Matrix<T2__, Eigen::Dynamic, Eigen::Dynamic>& L,
                    const Eigen::Matrix<T3__, Eigen::Dynamic, Eigen::Dynamic>& Lcov, std::ostream* pstream__) {
    typedef typename boost::math::tools::promote_args<T0__, T1__, T2__, T3__>::type local_scalar_t__;
    typedef local_scalar_t__ fun_return_scalar_t__;
    const static bool propto__ = true;
    (void) propto__;
        local_scalar_t__ DUMMY_VAR__(std::numeric_limits<double>::quiet_NaN());
        (void) DUMMY_VAR__;  // suppress unused var warning

    int current_statement_begin__ = -1;
    try {
        {
        current_statement_begin__ = 14;
        validate_non_negative_index("z_flat", "num_elements(z)", num_elements(z));
        Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> z_flat(num_elements(z));
        stan::math::initialize(z_flat, DUMMY_VAR__);
        stan::math::fill(z_flat, DUMMY_VAR__);
        stan::math::assign(z_flat,to_vector(z));

        current_statement_begin__ = 15;
        validate_non_negative_index("r", "num_elements(z)", num_elements(z));
        Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r(num_elements(z));
        stan::math::initialize(r, DUMMY_VAR__);
        stan::math::fill(r, DUMMY_VAR__);
        stan::math::assign(r,rep_vector(0, num_elements(z)));

        current_statement_begin__ = 16;
        validate_non_negative_index("LC", "rows(L)", rows(L));
        validate_non_negative_index("LC", "cols(L)", cols(L));
        Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> LC(rows(L), cols(L));
        stan::math::initialize(LC, DUMMY_VAR__);
        stan::math::fill(LC, DUMMY_VAR__);
        stan::math::assign(LC,diag_pre_multiply(SD, L));

        current_statement_begin__ = 17;
        int rows_z(0);
        (void) rows_z;  // dummy to suppress unused var warning
        stan::math::fill(rows_z, std::numeric_limits<int>::min());
        stan::math::assign(rows_z,rows(z));

        current_statement_begin__ = 18;
        int rows_L(0);
        (void) rows_L;  // dummy to suppress unused var warning
        stan::math::fill(rows_L, std::numeric_limits<int>::min());
        stan::math::assign(rows_L,rows(L));


        current_statement_begin__ = 20;
        for (int icov = 1; icov <= rows(Lcov); ++icov) {

            current_statement_begin__ = 21;
            for (int jcov = 1; jcov <= icov; ++jcov) {

                current_statement_begin__ = 22;
                if (as_bool(logical_gt(get_base1(Lcov, icov, jcov, "Lcov", 1), 1e-10))) {

                    current_statement_begin__ = 24;
                    for (int i = 1; i <= rows_L; ++i) {

                        current_statement_begin__ = 25;
                        for (int j = 1; j <= i; ++j) {
                            {
                            current_statement_begin__ = 27;
                            int k(0);
                            (void) k;  // dummy to suppress unused var warning
                            stan::math::fill(k, std::numeric_limits<int>::min());
                            stan::math::assign(k,((rows_L * (icov - 1)) + i));

                            current_statement_begin__ = 29;
                            int l(0);
                            (void) l;  // dummy to suppress unused var warning
                            stan::math::fill(l, std::numeric_limits<int>::min());
                            stan::math::assign(l,((rows_L * (jcov - 1)) + j));


                            current_statement_begin__ = 30;
                            stan::model::assign(r, 
                                        stan::model::cons_list(stan::model::index_uni(k), stan::model::nil_index_list()), 
                                        (get_base1(r, k, "r", 1) + ((get_base1(Lcov, icov, jcov, "Lcov", 1) * get_base1(LC, i, j, "LC", 1)) * get_base1(z_flat, l, "z_flat", 1))), 
                                        "assigning variable r");
                            }
                        }
                    }
                }
            }
        }
        current_statement_begin__ = 37;
        return stan::math::promote_scalar<fun_return_scalar_t__>(to_matrix(r, cols(z), rows(z), 0));
        }
    } 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 ***");
    }
}


struct scale_r_cor_cov_functor__ {
    template <typename T0__, typename T1__, typename T2__, typename T3__>
        Eigen::Matrix<typename boost::math::tools::promote_args<T0__, T1__, T2__, T3__>::type, Eigen::Dynamic, Eigen::Dynamic>
    operator()(const Eigen::Matrix<T0__, Eigen::Dynamic, Eigen::Dynamic>& z,
                    const Eigen::Matrix<T1__, Eigen::Dynamic, 1>& SD,
                    const Eigen::Matrix<T2__, Eigen::Dynamic, Eigen::Dynamic>& L,
                    const Eigen::Matrix<T3__, Eigen::Dynamic, Eigen::Dynamic>& Lcov, std::ostream* pstream__) const {
        return scale_r_cor_cov(z, SD, L, Lcov, pstream__);
    }
};

template <typename T0__, typename T1__, typename T2__>
Eigen::Matrix<typename boost::math::tools::promote_args<T0__, T1__, T2__>::type, Eigen::Dynamic, Eigen::Dynamic>
scale_r_cor(const Eigen::Matrix<T0__, Eigen::Dynamic, Eigen::Dynamic>& z,
                const Eigen::Matrix<T1__, Eigen::Dynamic, 1>& SD,
                const Eigen::Matrix<T2__, Eigen::Dynamic, Eigen::Dynamic>& L, std::ostream* pstream__) {
    typedef typename boost::math::tools::promote_args<T0__, T1__, T2__>::type local_scalar_t__;
    typedef local_scalar_t__ fun_return_scalar_t__;
    const static bool propto__ = true;
    (void) propto__;
        local_scalar_t__ DUMMY_VAR__(std::numeric_limits<double>::quiet_NaN());
        (void) DUMMY_VAR__;  // suppress unused var warning

    int current_statement_begin__ = -1;
    try {

        current_statement_begin__ = 49;
        return stan::math::promote_scalar<fun_return_scalar_t__>(transpose(multiply(diag_pre_multiply(SD, L), z)));
    } 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 ***");
    }
}


struct scale_r_cor_functor__ {
    template <typename T0__, typename T1__, typename T2__>
        Eigen::Matrix<typename boost::math::tools::promote_args<T0__, T1__, T2__>::type, Eigen::Dynamic, Eigen::Dynamic>
    operator()(const Eigen::Matrix<T0__, Eigen::Dynamic, Eigen::Dynamic>& z,
                const Eigen::Matrix<T1__, Eigen::Dynamic, 1>& SD,
                const Eigen::Matrix<T2__, Eigen::Dynamic, Eigen::Dynamic>& L, std::ostream* pstream__) const {
        return scale_r_cor(z, SD, L, pstream__);
    }
};

class model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb
  : public stan::model::model_base_crtp<model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb> {
private:
        int N;
        int N_bwt;
        vector_d Y_bwt;
        int K_bwt;
        matrix_d X_bwt;
        int K_sigma_bwt;
        matrix_d X_sigma_bwt;
        int N_tarsus;
        vector_d Y_tarsus;
        int K_tarsus;
        matrix_d X_tarsus;
        int K_sigma_tarsus;
        matrix_d X_sigma_tarsus;
        int nresp;
        int nrescor;
        int N_1;
        int M_1;
        std::vector<int> J_1_bwt;
        std::vector<int> J_1_tarsus;
        matrix_d Lcov_1;
        vector_d Z_1_bwt_1;
        vector_d Z_1_bwt_2;
        vector_d Z_1_tarsus_3;
        vector_d Z_1_tarsus_4;
        int NC_1;
        int N_2;
        int M_2;
        std::vector<int> J_2_bwt;
        std::vector<int> J_2_tarsus;
        vector_d Z_2_bwt_1;
        vector_d Z_2_bwt_2;
        vector_d Z_2_tarsus_3;
        vector_d Z_2_tarsus_4;
        int NC_2;
        int N_3;
        int M_3;
        std::vector<int> J_3_bwt;
        std::vector<int> J_3_tarsus;
        vector_d Z_3_bwt_1;
        vector_d Z_3_bwt_2;
        vector_d Z_3_tarsus_3;
        vector_d Z_3_tarsus_4;
        int NC_3;
        int prior_only;
        int Kc_bwt;
        matrix_d Xc_bwt;
        vector_d means_X_bwt;
        int Kc_tarsus;
        matrix_d Xc_tarsus;
        vector_d means_X_tarsus;
        std::vector<vector_d> Y;
public:
    model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb(rstan::io::rlist_ref_var_context& context__,
        std::ostream* pstream__ = 0)
        : model_base_crtp(0) {
        ctor_body(context__, 0, pstream__);
    }

    model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb(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__ = "model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb_namespace::model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb";
        (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__ = 53;
            context__.validate_dims("data initialization", "N", "int", context__.to_vec());
            N = int(0);
            vals_i__ = context__.vals_i("N");
            pos__ = 0;
            N = vals_i__[pos__++];
            check_greater_or_equal(function__, "N", N, 1);

            current_statement_begin__ = 54;
            context__.validate_dims("data initialization", "N_bwt", "int", context__.to_vec());
            N_bwt = int(0);
            vals_i__ = context__.vals_i("N_bwt");
            pos__ = 0;
            N_bwt = vals_i__[pos__++];
            check_greater_or_equal(function__, "N_bwt", N_bwt, 1);

            current_statement_begin__ = 55;
            validate_non_negative_index("Y_bwt", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Y_bwt", "vector_d", context__.to_vec(N_bwt));
            Y_bwt = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Y_bwt");
            pos__ = 0;
            size_t Y_bwt_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Y_bwt_j_1_max__; ++j_1__) {
                Y_bwt(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 56;
            context__.validate_dims("data initialization", "K_bwt", "int", context__.to_vec());
            K_bwt = int(0);
            vals_i__ = context__.vals_i("K_bwt");
            pos__ = 0;
            K_bwt = vals_i__[pos__++];
            check_greater_or_equal(function__, "K_bwt", K_bwt, 1);

            current_statement_begin__ = 57;
            validate_non_negative_index("X_bwt", "N_bwt", N_bwt);
            validate_non_negative_index("X_bwt", "K_bwt", K_bwt);
            context__.validate_dims("data initialization", "X_bwt", "matrix_d", context__.to_vec(N_bwt,K_bwt));
            X_bwt = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_bwt, K_bwt);
            vals_r__ = context__.vals_r("X_bwt");
            pos__ = 0;
            size_t X_bwt_j_2_max__ = K_bwt;
            size_t X_bwt_j_1_max__ = N_bwt;
            for (size_t j_2__ = 0; j_2__ < X_bwt_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < X_bwt_j_1_max__; ++j_1__) {
                    X_bwt(j_1__, j_2__) = vals_r__[pos__++];
                }
            }

            current_statement_begin__ = 58;
            context__.validate_dims("data initialization", "K_sigma_bwt", "int", context__.to_vec());
            K_sigma_bwt = int(0);
            vals_i__ = context__.vals_i("K_sigma_bwt");
            pos__ = 0;
            K_sigma_bwt = vals_i__[pos__++];
            check_greater_or_equal(function__, "K_sigma_bwt", K_sigma_bwt, 1);

            current_statement_begin__ = 59;
            validate_non_negative_index("X_sigma_bwt", "N_bwt", N_bwt);
            validate_non_negative_index("X_sigma_bwt", "K_sigma_bwt", K_sigma_bwt);
            context__.validate_dims("data initialization", "X_sigma_bwt", "matrix_d", context__.to_vec(N_bwt,K_sigma_bwt));
            X_sigma_bwt = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_bwt, K_sigma_bwt);
            vals_r__ = context__.vals_r("X_sigma_bwt");
            pos__ = 0;
            size_t X_sigma_bwt_j_2_max__ = K_sigma_bwt;
            size_t X_sigma_bwt_j_1_max__ = N_bwt;
            for (size_t j_2__ = 0; j_2__ < X_sigma_bwt_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < X_sigma_bwt_j_1_max__; ++j_1__) {
                    X_sigma_bwt(j_1__, j_2__) = vals_r__[pos__++];
                }
            }

            current_statement_begin__ = 60;
            context__.validate_dims("data initialization", "N_tarsus", "int", context__.to_vec());
            N_tarsus = int(0);
            vals_i__ = context__.vals_i("N_tarsus");
            pos__ = 0;
            N_tarsus = vals_i__[pos__++];
            check_greater_or_equal(function__, "N_tarsus", N_tarsus, 1);

            current_statement_begin__ = 61;
            validate_non_negative_index("Y_tarsus", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Y_tarsus", "vector_d", context__.to_vec(N_tarsus));
            Y_tarsus = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Y_tarsus");
            pos__ = 0;
            size_t Y_tarsus_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Y_tarsus_j_1_max__; ++j_1__) {
                Y_tarsus(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 62;
            context__.validate_dims("data initialization", "K_tarsus", "int", context__.to_vec());
            K_tarsus = int(0);
            vals_i__ = context__.vals_i("K_tarsus");
            pos__ = 0;
            K_tarsus = vals_i__[pos__++];
            check_greater_or_equal(function__, "K_tarsus", K_tarsus, 1);

            current_statement_begin__ = 63;
            validate_non_negative_index("X_tarsus", "N_tarsus", N_tarsus);
            validate_non_negative_index("X_tarsus", "K_tarsus", K_tarsus);
            context__.validate_dims("data initialization", "X_tarsus", "matrix_d", context__.to_vec(N_tarsus,K_tarsus));
            X_tarsus = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_tarsus, K_tarsus);
            vals_r__ = context__.vals_r("X_tarsus");
            pos__ = 0;
            size_t X_tarsus_j_2_max__ = K_tarsus;
            size_t X_tarsus_j_1_max__ = N_tarsus;
            for (size_t j_2__ = 0; j_2__ < X_tarsus_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < X_tarsus_j_1_max__; ++j_1__) {
                    X_tarsus(j_1__, j_2__) = vals_r__[pos__++];
                }
            }

            current_statement_begin__ = 64;
            context__.validate_dims("data initialization", "K_sigma_tarsus", "int", context__.to_vec());
            K_sigma_tarsus = int(0);
            vals_i__ = context__.vals_i("K_sigma_tarsus");
            pos__ = 0;
            K_sigma_tarsus = vals_i__[pos__++];
            check_greater_or_equal(function__, "K_sigma_tarsus", K_sigma_tarsus, 1);

            current_statement_begin__ = 65;
            validate_non_negative_index("X_sigma_tarsus", "N_tarsus", N_tarsus);
            validate_non_negative_index("X_sigma_tarsus", "K_sigma_tarsus", K_sigma_tarsus);
            context__.validate_dims("data initialization", "X_sigma_tarsus", "matrix_d", context__.to_vec(N_tarsus,K_sigma_tarsus));
            X_sigma_tarsus = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_tarsus, K_sigma_tarsus);
            vals_r__ = context__.vals_r("X_sigma_tarsus");
            pos__ = 0;
            size_t X_sigma_tarsus_j_2_max__ = K_sigma_tarsus;
            size_t X_sigma_tarsus_j_1_max__ = N_tarsus;
            for (size_t j_2__ = 0; j_2__ < X_sigma_tarsus_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < X_sigma_tarsus_j_1_max__; ++j_1__) {
                    X_sigma_tarsus(j_1__, j_2__) = vals_r__[pos__++];
                }
            }

            current_statement_begin__ = 66;
            context__.validate_dims("data initialization", "nresp", "int", context__.to_vec());
            nresp = int(0);
            vals_i__ = context__.vals_i("nresp");
            pos__ = 0;
            nresp = vals_i__[pos__++];
            check_greater_or_equal(function__, "nresp", nresp, 1);

            current_statement_begin__ = 67;
            context__.validate_dims("data initialization", "nrescor", "int", context__.to_vec());
            nrescor = int(0);
            vals_i__ = context__.vals_i("nrescor");
            pos__ = 0;
            nrescor = vals_i__[pos__++];

            current_statement_begin__ = 69;
            context__.validate_dims("data initialization", "N_1", "int", context__.to_vec());
            N_1 = int(0);
            vals_i__ = context__.vals_i("N_1");
            pos__ = 0;
            N_1 = vals_i__[pos__++];
            check_greater_or_equal(function__, "N_1", N_1, 1);

            current_statement_begin__ = 70;
            context__.validate_dims("data initialization", "M_1", "int", context__.to_vec());
            M_1 = int(0);
            vals_i__ = context__.vals_i("M_1");
            pos__ = 0;
            M_1 = vals_i__[pos__++];
            check_greater_or_equal(function__, "M_1", M_1, 1);

            current_statement_begin__ = 71;
            validate_non_negative_index("J_1_bwt", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "J_1_bwt", "int", context__.to_vec(N_bwt));
            J_1_bwt = std::vector<int>(N_bwt, int(0));
            vals_i__ = context__.vals_i("J_1_bwt");
            pos__ = 0;
            size_t J_1_bwt_k_0_max__ = N_bwt;
            for (size_t k_0__ = 0; k_0__ < J_1_bwt_k_0_max__; ++k_0__) {
                J_1_bwt[k_0__] = vals_i__[pos__++];
            }
            size_t J_1_bwt_i_0_max__ = N_bwt;
            for (size_t i_0__ = 0; i_0__ < J_1_bwt_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_1_bwt[i_0__]", J_1_bwt[i_0__], 1);
            }

            current_statement_begin__ = 72;
            validate_non_negative_index("J_1_tarsus", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "J_1_tarsus", "int", context__.to_vec(N_tarsus));
            J_1_tarsus = std::vector<int>(N_tarsus, int(0));
            vals_i__ = context__.vals_i("J_1_tarsus");
            pos__ = 0;
            size_t J_1_tarsus_k_0_max__ = N_tarsus;
            for (size_t k_0__ = 0; k_0__ < J_1_tarsus_k_0_max__; ++k_0__) {
                J_1_tarsus[k_0__] = vals_i__[pos__++];
            }
            size_t J_1_tarsus_i_0_max__ = N_tarsus;
            for (size_t i_0__ = 0; i_0__ < J_1_tarsus_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_1_tarsus[i_0__]", J_1_tarsus[i_0__], 1);
            }

            current_statement_begin__ = 73;
            validate_non_negative_index("Lcov_1", "N_1", N_1);
            validate_non_negative_index("Lcov_1", "N_1", N_1);
            context__.validate_dims("data initialization", "Lcov_1", "matrix_d", context__.to_vec(N_1,N_1));
            Lcov_1 = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_1, N_1);
            vals_r__ = context__.vals_r("Lcov_1");
            pos__ = 0;
            size_t Lcov_1_j_2_max__ = N_1;
            size_t Lcov_1_j_1_max__ = N_1;
            for (size_t j_2__ = 0; j_2__ < Lcov_1_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < Lcov_1_j_1_max__; ++j_1__) {
                    Lcov_1(j_1__, j_2__) = vals_r__[pos__++];
                }
            }

            current_statement_begin__ = 75;
            validate_non_negative_index("Z_1_bwt_1", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_1_bwt_1", "vector_d", context__.to_vec(N_bwt));
            Z_1_bwt_1 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_1_bwt_1");
            pos__ = 0;
            size_t Z_1_bwt_1_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_1_bwt_1_j_1_max__; ++j_1__) {
                Z_1_bwt_1(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 76;
            validate_non_negative_index("Z_1_bwt_2", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_1_bwt_2", "vector_d", context__.to_vec(N_bwt));
            Z_1_bwt_2 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_1_bwt_2");
            pos__ = 0;
            size_t Z_1_bwt_2_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_1_bwt_2_j_1_max__; ++j_1__) {
                Z_1_bwt_2(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 77;
            validate_non_negative_index("Z_1_tarsus_3", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_1_tarsus_3", "vector_d", context__.to_vec(N_tarsus));
            Z_1_tarsus_3 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_1_tarsus_3");
            pos__ = 0;
            size_t Z_1_tarsus_3_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_1_tarsus_3_j_1_max__; ++j_1__) {
                Z_1_tarsus_3(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 78;
            validate_non_negative_index("Z_1_tarsus_4", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_1_tarsus_4", "vector_d", context__.to_vec(N_tarsus));
            Z_1_tarsus_4 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_1_tarsus_4");
            pos__ = 0;
            size_t Z_1_tarsus_4_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_1_tarsus_4_j_1_max__; ++j_1__) {
                Z_1_tarsus_4(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 79;
            context__.validate_dims("data initialization", "NC_1", "int", context__.to_vec());
            NC_1 = int(0);
            vals_i__ = context__.vals_i("NC_1");
            pos__ = 0;
            NC_1 = vals_i__[pos__++];
            check_greater_or_equal(function__, "NC_1", NC_1, 1);

            current_statement_begin__ = 81;
            context__.validate_dims("data initialization", "N_2", "int", context__.to_vec());
            N_2 = int(0);
            vals_i__ = context__.vals_i("N_2");
            pos__ = 0;
            N_2 = vals_i__[pos__++];
            check_greater_or_equal(function__, "N_2", N_2, 1);

            current_statement_begin__ = 82;
            context__.validate_dims("data initialization", "M_2", "int", context__.to_vec());
            M_2 = int(0);
            vals_i__ = context__.vals_i("M_2");
            pos__ = 0;
            M_2 = vals_i__[pos__++];
            check_greater_or_equal(function__, "M_2", M_2, 1);

            current_statement_begin__ = 83;
            validate_non_negative_index("J_2_bwt", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "J_2_bwt", "int", context__.to_vec(N_bwt));
            J_2_bwt = std::vector<int>(N_bwt, int(0));
            vals_i__ = context__.vals_i("J_2_bwt");
            pos__ = 0;
            size_t J_2_bwt_k_0_max__ = N_bwt;
            for (size_t k_0__ = 0; k_0__ < J_2_bwt_k_0_max__; ++k_0__) {
                J_2_bwt[k_0__] = vals_i__[pos__++];
            }
            size_t J_2_bwt_i_0_max__ = N_bwt;
            for (size_t i_0__ = 0; i_0__ < J_2_bwt_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_2_bwt[i_0__]", J_2_bwt[i_0__], 1);
            }

            current_statement_begin__ = 84;
            validate_non_negative_index("J_2_tarsus", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "J_2_tarsus", "int", context__.to_vec(N_tarsus));
            J_2_tarsus = std::vector<int>(N_tarsus, int(0));
            vals_i__ = context__.vals_i("J_2_tarsus");
            pos__ = 0;
            size_t J_2_tarsus_k_0_max__ = N_tarsus;
            for (size_t k_0__ = 0; k_0__ < J_2_tarsus_k_0_max__; ++k_0__) {
                J_2_tarsus[k_0__] = vals_i__[pos__++];
            }
            size_t J_2_tarsus_i_0_max__ = N_tarsus;
            for (size_t i_0__ = 0; i_0__ < J_2_tarsus_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_2_tarsus[i_0__]", J_2_tarsus[i_0__], 1);
            }

            current_statement_begin__ = 86;
            validate_non_negative_index("Z_2_bwt_1", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_2_bwt_1", "vector_d", context__.to_vec(N_bwt));
            Z_2_bwt_1 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_2_bwt_1");
            pos__ = 0;
            size_t Z_2_bwt_1_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_2_bwt_1_j_1_max__; ++j_1__) {
                Z_2_bwt_1(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 87;
            validate_non_negative_index("Z_2_bwt_2", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_2_bwt_2", "vector_d", context__.to_vec(N_bwt));
            Z_2_bwt_2 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_2_bwt_2");
            pos__ = 0;
            size_t Z_2_bwt_2_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_2_bwt_2_j_1_max__; ++j_1__) {
                Z_2_bwt_2(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 88;
            validate_non_negative_index("Z_2_tarsus_3", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_2_tarsus_3", "vector_d", context__.to_vec(N_tarsus));
            Z_2_tarsus_3 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_2_tarsus_3");
            pos__ = 0;
            size_t Z_2_tarsus_3_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_2_tarsus_3_j_1_max__; ++j_1__) {
                Z_2_tarsus_3(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 89;
            validate_non_negative_index("Z_2_tarsus_4", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_2_tarsus_4", "vector_d", context__.to_vec(N_tarsus));
            Z_2_tarsus_4 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_2_tarsus_4");
            pos__ = 0;
            size_t Z_2_tarsus_4_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_2_tarsus_4_j_1_max__; ++j_1__) {
                Z_2_tarsus_4(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 90;
            context__.validate_dims("data initialization", "NC_2", "int", context__.to_vec());
            NC_2 = int(0);
            vals_i__ = context__.vals_i("NC_2");
            pos__ = 0;
            NC_2 = vals_i__[pos__++];
            check_greater_or_equal(function__, "NC_2", NC_2, 1);

            current_statement_begin__ = 92;
            context__.validate_dims("data initialization", "N_3", "int", context__.to_vec());
            N_3 = int(0);
            vals_i__ = context__.vals_i("N_3");
            pos__ = 0;
            N_3 = vals_i__[pos__++];
            check_greater_or_equal(function__, "N_3", N_3, 1);

            current_statement_begin__ = 93;
            context__.validate_dims("data initialization", "M_3", "int", context__.to_vec());
            M_3 = int(0);
            vals_i__ = context__.vals_i("M_3");
            pos__ = 0;
            M_3 = vals_i__[pos__++];
            check_greater_or_equal(function__, "M_3", M_3, 1);

            current_statement_begin__ = 94;
            validate_non_negative_index("J_3_bwt", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "J_3_bwt", "int", context__.to_vec(N_bwt));
            J_3_bwt = std::vector<int>(N_bwt, int(0));
            vals_i__ = context__.vals_i("J_3_bwt");
            pos__ = 0;
            size_t J_3_bwt_k_0_max__ = N_bwt;
            for (size_t k_0__ = 0; k_0__ < J_3_bwt_k_0_max__; ++k_0__) {
                J_3_bwt[k_0__] = vals_i__[pos__++];
            }
            size_t J_3_bwt_i_0_max__ = N_bwt;
            for (size_t i_0__ = 0; i_0__ < J_3_bwt_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_3_bwt[i_0__]", J_3_bwt[i_0__], 1);
            }

            current_statement_begin__ = 95;
            validate_non_negative_index("J_3_tarsus", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "J_3_tarsus", "int", context__.to_vec(N_tarsus));
            J_3_tarsus = std::vector<int>(N_tarsus, int(0));
            vals_i__ = context__.vals_i("J_3_tarsus");
            pos__ = 0;
            size_t J_3_tarsus_k_0_max__ = N_tarsus;
            for (size_t k_0__ = 0; k_0__ < J_3_tarsus_k_0_max__; ++k_0__) {
                J_3_tarsus[k_0__] = vals_i__[pos__++];
            }
            size_t J_3_tarsus_i_0_max__ = N_tarsus;
            for (size_t i_0__ = 0; i_0__ < J_3_tarsus_i_0_max__; ++i_0__) {
                check_greater_or_equal(function__, "J_3_tarsus[i_0__]", J_3_tarsus[i_0__], 1);
            }

            current_statement_begin__ = 97;
            validate_non_negative_index("Z_3_bwt_1", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_3_bwt_1", "vector_d", context__.to_vec(N_bwt));
            Z_3_bwt_1 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_3_bwt_1");
            pos__ = 0;
            size_t Z_3_bwt_1_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_3_bwt_1_j_1_max__; ++j_1__) {
                Z_3_bwt_1(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 98;
            validate_non_negative_index("Z_3_bwt_2", "N_bwt", N_bwt);
            context__.validate_dims("data initialization", "Z_3_bwt_2", "vector_d", context__.to_vec(N_bwt));
            Z_3_bwt_2 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_bwt);
            vals_r__ = context__.vals_r("Z_3_bwt_2");
            pos__ = 0;
            size_t Z_3_bwt_2_j_1_max__ = N_bwt;
            for (size_t j_1__ = 0; j_1__ < Z_3_bwt_2_j_1_max__; ++j_1__) {
                Z_3_bwt_2(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 99;
            validate_non_negative_index("Z_3_tarsus_3", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_3_tarsus_3", "vector_d", context__.to_vec(N_tarsus));
            Z_3_tarsus_3 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_3_tarsus_3");
            pos__ = 0;
            size_t Z_3_tarsus_3_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_3_tarsus_3_j_1_max__; ++j_1__) {
                Z_3_tarsus_3(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 100;
            validate_non_negative_index("Z_3_tarsus_4", "N_tarsus", N_tarsus);
            context__.validate_dims("data initialization", "Z_3_tarsus_4", "vector_d", context__.to_vec(N_tarsus));
            Z_3_tarsus_4 = Eigen::Matrix<double, Eigen::Dynamic, 1>(N_tarsus);
            vals_r__ = context__.vals_r("Z_3_tarsus_4");
            pos__ = 0;
            size_t Z_3_tarsus_4_j_1_max__ = N_tarsus;
            for (size_t j_1__ = 0; j_1__ < Z_3_tarsus_4_j_1_max__; ++j_1__) {
                Z_3_tarsus_4(j_1__) = vals_r__[pos__++];
            }

            current_statement_begin__ = 101;
            context__.validate_dims("data initialization", "NC_3", "int", context__.to_vec());
            NC_3 = int(0);
            vals_i__ = context__.vals_i("NC_3");
            pos__ = 0;
            NC_3 = vals_i__[pos__++];
            check_greater_or_equal(function__, "NC_3", NC_3, 1);

            current_statement_begin__ = 102;
            context__.validate_dims("data initialization", "prior_only", "int", context__.to_vec());
            prior_only = int(0);
            vals_i__ = context__.vals_i("prior_only");
            pos__ = 0;
            prior_only = vals_i__[pos__++];


            // initialize transformed data variables
            current_statement_begin__ = 105;
            Kc_bwt = int(0);
            stan::math::fill(Kc_bwt, std::numeric_limits<int>::min());
            stan::math::assign(Kc_bwt,(K_bwt - 1));

            current_statement_begin__ = 106;
            validate_non_negative_index("Xc_bwt", "N_bwt", N_bwt);
            validate_non_negative_index("Xc_bwt", "Kc_bwt", Kc_bwt);
            Xc_bwt = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_bwt, Kc_bwt);
            stan::math::fill(Xc_bwt, DUMMY_VAR__);

            current_statement_begin__ = 107;
            validate_non_negative_index("means_X_bwt", "Kc_bwt", Kc_bwt);
            means_X_bwt = Eigen::Matrix<double, Eigen::Dynamic, 1>(Kc_bwt);
            stan::math::fill(means_X_bwt, DUMMY_VAR__);

            current_statement_begin__ = 108;
            Kc_tarsus = int(0);
            stan::math::fill(Kc_tarsus, std::numeric_limits<int>::min());
            stan::math::assign(Kc_tarsus,(K_tarsus - 1));

            current_statement_begin__ = 109;
            validate_non_negative_index("Xc_tarsus", "N_tarsus", N_tarsus);
            validate_non_negative_index("Xc_tarsus", "Kc_tarsus", Kc_tarsus);
            Xc_tarsus = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>(N_tarsus, Kc_tarsus);
            stan::math::fill(Xc_tarsus, DUMMY_VAR__);

            current_statement_begin__ = 110;
            validate_non_negative_index("means_X_tarsus", "Kc_tarsus", Kc_tarsus);
            means_X_tarsus = Eigen::Matrix<double, Eigen::Dynamic, 1>(Kc_tarsus);
            stan::math::fill(means_X_tarsus, DUMMY_VAR__);

            current_statement_begin__ = 111;
            validate_non_negative_index("Y", "nresp", nresp);
            validate_non_negative_index("Y", "N", N);
            Y = std::vector<Eigen::Matrix<double, Eigen::Dynamic, 1> >(N, Eigen::Matrix<double, Eigen::Dynamic, 1>(nresp));
            stan::math::fill(Y, DUMMY_VAR__);

            // execute transformed data statements
            current_statement_begin__ = 112;
            for (int i = 2; i <= K_bwt; ++i) {

                current_statement_begin__ = 113;
                stan::model::assign(means_X_bwt, 
                            stan::model::cons_list(stan::model::index_uni((i - 1)), stan::model::nil_index_list()), 
                            mean(stan::model::rvalue(X_bwt, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(i), stan::model::nil_index_list())), "X_bwt")), 
                            "assigning variable means_X_bwt");
                current_statement_begin__ = 114;
                stan::model::assign(Xc_bwt, 
                            stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni((i - 1)), stan::model::nil_index_list())), 
                            subtract(stan::model::rvalue(X_bwt, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(i), stan::model::nil_index_list())), "X_bwt"), get_base1(means_X_bwt, (i - 1), "means_X_bwt", 1)), 
                            "assigning variable Xc_bwt");
            }
            current_statement_begin__ = 116;
            for (int i = 2; i <= K_tarsus; ++i) {

                current_statement_begin__ = 117;
                stan::model::assign(means_X_tarsus, 
                            stan::model::cons_list(stan::model::index_uni((i - 1)), stan::model::nil_index_list()), 
                            mean(stan::model::rvalue(X_tarsus, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(i), stan::model::nil_index_list())), "X_tarsus")), 
                            "assigning variable means_X_tarsus");
                current_statement_begin__ = 118;
                stan::model::assign(Xc_tarsus, 
                            stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni((i - 1)), stan::model::nil_index_list())), 
                            subtract(stan::model::rvalue(X_tarsus, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(i), stan::model::nil_index_list())), "X_tarsus"), get_base1(means_X_tarsus, (i - 1), "means_X_tarsus", 1)), 
                            "assigning variable Xc_tarsus");
            }
            current_statement_begin__ = 120;
            for (int n = 1; n <= N; ++n) {

                current_statement_begin__ = 121;
                stan::model::assign(Y, 
                            stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                            transpose(stan::math::to_row_vector(stan::math::array_builder<double >().add(get_base1(Y_bwt, n, "Y_bwt", 1)).add(get_base1(Y_tarsus, n, "Y_tarsus", 1)).array())), 
                            "assigning variable Y");
            }

            // validate transformed data

            // validate, set parameter ranges
            num_params_r__ = 0U;
            param_ranges_i__.clear();
            current_statement_begin__ = 125;
            validate_non_negative_index("b_bwt", "Kc_bwt", Kc_bwt);
            num_params_r__ += Kc_bwt;
            current_statement_begin__ = 126;
            num_params_r__ += 1;
            current_statement_begin__ = 127;
            validate_non_negative_index("b_sigma_bwt", "K_sigma_bwt", K_sigma_bwt);
            num_params_r__ += K_sigma_bwt;
            current_statement_begin__ = 128;
            validate_non_negative_index("b_tarsus", "Kc_tarsus", Kc_tarsus);
            num_params_r__ += Kc_tarsus;
            current_statement_begin__ = 129;
            num_params_r__ += 1;
            current_statement_begin__ = 130;
            validate_non_negative_index("b_sigma_tarsus", "K_sigma_tarsus", K_sigma_tarsus);
            num_params_r__ += K_sigma_tarsus;
            current_statement_begin__ = 131;
            validate_non_negative_index("Lrescor", "nresp", nresp);
            validate_non_negative_index("Lrescor", "nresp", nresp);
            num_params_r__ += ((nresp * (nresp - 1)) / 2);
            current_statement_begin__ = 132;
            validate_non_negative_index("sd_1", "M_1", M_1);
            num_params_r__ += M_1;
            current_statement_begin__ = 133;
            validate_non_negative_index("z_1", "M_1", M_1);
            validate_non_negative_index("z_1", "N_1", N_1);
            num_params_r__ += (M_1 * N_1);
            current_statement_begin__ = 134;
            validate_non_negative_index("L_1", "M_1", M_1);
            validate_non_negative_index("L_1", "M_1", M_1);
            num_params_r__ += ((M_1 * (M_1 - 1)) / 2);
            current_statement_begin__ = 135;
            validate_non_negative_index("sd_2", "M_2", M_2);
            num_params_r__ += M_2;
            current_statement_begin__ = 136;
            validate_non_negative_index("z_2", "M_2", M_2);
            validate_non_negative_index("z_2", "N_2", N_2);
            num_params_r__ += (M_2 * N_2);
            current_statement_begin__ = 137;
            validate_non_negative_index("L_2", "M_2", M_2);
            validate_non_negative_index("L_2", "M_2", M_2);
            num_params_r__ += ((M_2 * (M_2 - 1)) / 2);
            current_statement_begin__ = 138;
            validate_non_negative_index("sd_3", "M_3", M_3);
            num_params_r__ += M_3;
            current_statement_begin__ = 139;
            validate_non_negative_index("z_3", "M_3", M_3);
            validate_non_negative_index("z_3", "N_3", N_3);
            num_params_r__ += (M_3 * N_3);
            current_statement_begin__ = 140;
            validate_non_negative_index("L_3", "M_3", M_3);
            validate_non_negative_index("L_3", "M_3", M_3);
            num_params_r__ += ((M_3 * (M_3 - 1)) / 2);
        } 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 ***");
        }
    }

    ~model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb() { }


    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__ = 125;
        if (!(context__.contains_r("b_bwt")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable b_bwt missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("b_bwt");
        pos__ = 0U;
        validate_non_negative_index("b_bwt", "Kc_bwt", Kc_bwt);
        context__.validate_dims("parameter initialization", "b_bwt", "vector_d", context__.to_vec(Kc_bwt));
        Eigen::Matrix<double, Eigen::Dynamic, 1> b_bwt(Kc_bwt);
        size_t b_bwt_j_1_max__ = Kc_bwt;
        for (size_t j_1__ = 0; j_1__ < b_bwt_j_1_max__; ++j_1__) {
            b_bwt(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_unconstrain(b_bwt);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable b_bwt: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 126;
        if (!(context__.contains_r("Intercept_bwt")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable Intercept_bwt missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("Intercept_bwt");
        pos__ = 0U;
        context__.validate_dims("parameter initialization", "Intercept_bwt", "double", context__.to_vec());
        double Intercept_bwt(0);
        Intercept_bwt = vals_r__[pos__++];
        try {
            writer__.scalar_unconstrain(Intercept_bwt);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable Intercept_bwt: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 127;
        if (!(context__.contains_r("b_sigma_bwt")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable b_sigma_bwt missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("b_sigma_bwt");
        pos__ = 0U;
        validate_non_negative_index("b_sigma_bwt", "K_sigma_bwt", K_sigma_bwt);
        context__.validate_dims("parameter initialization", "b_sigma_bwt", "vector_d", context__.to_vec(K_sigma_bwt));
        Eigen::Matrix<double, Eigen::Dynamic, 1> b_sigma_bwt(K_sigma_bwt);
        size_t b_sigma_bwt_j_1_max__ = K_sigma_bwt;
        for (size_t j_1__ = 0; j_1__ < b_sigma_bwt_j_1_max__; ++j_1__) {
            b_sigma_bwt(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_unconstrain(b_sigma_bwt);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable b_sigma_bwt: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 128;
        if (!(context__.contains_r("b_tarsus")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable b_tarsus missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("b_tarsus");
        pos__ = 0U;
        validate_non_negative_index("b_tarsus", "Kc_tarsus", Kc_tarsus);
        context__.validate_dims("parameter initialization", "b_tarsus", "vector_d", context__.to_vec(Kc_tarsus));
        Eigen::Matrix<double, Eigen::Dynamic, 1> b_tarsus(Kc_tarsus);
        size_t b_tarsus_j_1_max__ = Kc_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_tarsus_j_1_max__; ++j_1__) {
            b_tarsus(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_unconstrain(b_tarsus);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable b_tarsus: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 129;
        if (!(context__.contains_r("Intercept_tarsus")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable Intercept_tarsus missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("Intercept_tarsus");
        pos__ = 0U;
        context__.validate_dims("parameter initialization", "Intercept_tarsus", "double", context__.to_vec());
        double Intercept_tarsus(0);
        Intercept_tarsus = vals_r__[pos__++];
        try {
            writer__.scalar_unconstrain(Intercept_tarsus);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable Intercept_tarsus: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 130;
        if (!(context__.contains_r("b_sigma_tarsus")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable b_sigma_tarsus missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("b_sigma_tarsus");
        pos__ = 0U;
        validate_non_negative_index("b_sigma_tarsus", "K_sigma_tarsus", K_sigma_tarsus);
        context__.validate_dims("parameter initialization", "b_sigma_tarsus", "vector_d", context__.to_vec(K_sigma_tarsus));
        Eigen::Matrix<double, Eigen::Dynamic, 1> b_sigma_tarsus(K_sigma_tarsus);
        size_t b_sigma_tarsus_j_1_max__ = K_sigma_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_sigma_tarsus_j_1_max__; ++j_1__) {
            b_sigma_tarsus(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_unconstrain(b_sigma_tarsus);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable b_sigma_tarsus: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 131;
        if (!(context__.contains_r("Lrescor")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable Lrescor missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("Lrescor");
        pos__ = 0U;
        validate_non_negative_index("Lrescor", "nresp", nresp);
        validate_non_negative_index("Lrescor", "nresp", nresp);
        context__.validate_dims("parameter initialization", "Lrescor", "matrix_d", context__.to_vec(nresp,nresp));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Lrescor(nresp, nresp);
        size_t Lrescor_j_2_max__ = nresp;
        size_t Lrescor_j_1_max__ = nresp;
        for (size_t j_2__ = 0; j_2__ < Lrescor_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Lrescor_j_1_max__; ++j_1__) {
                Lrescor(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.cholesky_factor_corr_unconstrain(Lrescor);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable Lrescor: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 132;
        if (!(context__.contains_r("sd_1")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable sd_1 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("sd_1");
        pos__ = 0U;
        validate_non_negative_index("sd_1", "M_1", M_1);
        context__.validate_dims("parameter initialization", "sd_1", "vector_d", context__.to_vec(M_1));
        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_1(M_1);
        size_t sd_1_j_1_max__ = M_1;
        for (size_t j_1__ = 0; j_1__ < sd_1_j_1_max__; ++j_1__) {
            sd_1(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_lb_unconstrain(0, sd_1);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable sd_1: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 133;
        if (!(context__.contains_r("z_1")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable z_1 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("z_1");
        pos__ = 0U;
        validate_non_negative_index("z_1", "M_1", M_1);
        validate_non_negative_index("z_1", "N_1", N_1);
        context__.validate_dims("parameter initialization", "z_1", "matrix_d", context__.to_vec(M_1,N_1));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_1(M_1, N_1);
        size_t z_1_j_2_max__ = N_1;
        size_t z_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < z_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_1_j_1_max__; ++j_1__) {
                z_1(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.matrix_unconstrain(z_1);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable z_1: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 134;
        if (!(context__.contains_r("L_1")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable L_1 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("L_1");
        pos__ = 0U;
        validate_non_negative_index("L_1", "M_1", M_1);
        validate_non_negative_index("L_1", "M_1", M_1);
        context__.validate_dims("parameter initialization", "L_1", "matrix_d", context__.to_vec(M_1,M_1));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_1(M_1, M_1);
        size_t L_1_j_2_max__ = M_1;
        size_t L_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < L_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_1_j_1_max__; ++j_1__) {
                L_1(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.cholesky_factor_corr_unconstrain(L_1);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable L_1: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 135;
        if (!(context__.contains_r("sd_2")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable sd_2 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("sd_2");
        pos__ = 0U;
        validate_non_negative_index("sd_2", "M_2", M_2);
        context__.validate_dims("parameter initialization", "sd_2", "vector_d", context__.to_vec(M_2));
        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_2(M_2);
        size_t sd_2_j_1_max__ = M_2;
        for (size_t j_1__ = 0; j_1__ < sd_2_j_1_max__; ++j_1__) {
            sd_2(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_lb_unconstrain(0, sd_2);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable sd_2: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 136;
        if (!(context__.contains_r("z_2")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable z_2 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("z_2");
        pos__ = 0U;
        validate_non_negative_index("z_2", "M_2", M_2);
        validate_non_negative_index("z_2", "N_2", N_2);
        context__.validate_dims("parameter initialization", "z_2", "matrix_d", context__.to_vec(M_2,N_2));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_2(M_2, N_2);
        size_t z_2_j_2_max__ = N_2;
        size_t z_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < z_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_2_j_1_max__; ++j_1__) {
                z_2(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.matrix_unconstrain(z_2);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable z_2: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 137;
        if (!(context__.contains_r("L_2")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable L_2 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("L_2");
        pos__ = 0U;
        validate_non_negative_index("L_2", "M_2", M_2);
        validate_non_negative_index("L_2", "M_2", M_2);
        context__.validate_dims("parameter initialization", "L_2", "matrix_d", context__.to_vec(M_2,M_2));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_2(M_2, M_2);
        size_t L_2_j_2_max__ = M_2;
        size_t L_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < L_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_2_j_1_max__; ++j_1__) {
                L_2(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.cholesky_factor_corr_unconstrain(L_2);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable L_2: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 138;
        if (!(context__.contains_r("sd_3")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable sd_3 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("sd_3");
        pos__ = 0U;
        validate_non_negative_index("sd_3", "M_3", M_3);
        context__.validate_dims("parameter initialization", "sd_3", "vector_d", context__.to_vec(M_3));
        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_3(M_3);
        size_t sd_3_j_1_max__ = M_3;
        for (size_t j_1__ = 0; j_1__ < sd_3_j_1_max__; ++j_1__) {
            sd_3(j_1__) = vals_r__[pos__++];
        }
        try {
            writer__.vector_lb_unconstrain(0, sd_3);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable sd_3: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 139;
        if (!(context__.contains_r("z_3")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable z_3 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("z_3");
        pos__ = 0U;
        validate_non_negative_index("z_3", "M_3", M_3);
        validate_non_negative_index("z_3", "N_3", N_3);
        context__.validate_dims("parameter initialization", "z_3", "matrix_d", context__.to_vec(M_3,N_3));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_3(M_3, N_3);
        size_t z_3_j_2_max__ = N_3;
        size_t z_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < z_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_3_j_1_max__; ++j_1__) {
                z_3(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.matrix_unconstrain(z_3);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable z_3: ") + e.what()), current_statement_begin__, prog_reader__());
        }

        current_statement_begin__ = 140;
        if (!(context__.contains_r("L_3")))
            stan::lang::rethrow_located(std::runtime_error(std::string("Variable L_3 missing")), current_statement_begin__, prog_reader__());
        vals_r__ = context__.vals_r("L_3");
        pos__ = 0U;
        validate_non_negative_index("L_3", "M_3", M_3);
        validate_non_negative_index("L_3", "M_3", M_3);
        context__.validate_dims("parameter initialization", "L_3", "matrix_d", context__.to_vec(M_3,M_3));
        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_3(M_3, M_3);
        size_t L_3_j_2_max__ = M_3;
        size_t L_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < L_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_3_j_1_max__; ++j_1__) {
                L_3(j_1__, j_2__) = vals_r__[pos__++];
            }
        }
        try {
            writer__.cholesky_factor_corr_unconstrain(L_3);
        } catch (const std::exception& e) {
            stan::lang::rethrow_located(std::runtime_error(std::string("Error transforming variable L_3: ") + 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__ = 125;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> b_bwt;
            (void) b_bwt;  // dummy to suppress unused var warning
            if (jacobian__)
                b_bwt = in__.vector_constrain(Kc_bwt, lp__);
            else
                b_bwt = in__.vector_constrain(Kc_bwt);

            current_statement_begin__ = 126;
            local_scalar_t__ Intercept_bwt;
            (void) Intercept_bwt;  // dummy to suppress unused var warning
            if (jacobian__)
                Intercept_bwt = in__.scalar_constrain(lp__);
            else
                Intercept_bwt = in__.scalar_constrain();

            current_statement_begin__ = 127;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> b_sigma_bwt;
            (void) b_sigma_bwt;  // dummy to suppress unused var warning
            if (jacobian__)
                b_sigma_bwt = in__.vector_constrain(K_sigma_bwt, lp__);
            else
                b_sigma_bwt = in__.vector_constrain(K_sigma_bwt);

            current_statement_begin__ = 128;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> b_tarsus;
            (void) b_tarsus;  // dummy to suppress unused var warning
            if (jacobian__)
                b_tarsus = in__.vector_constrain(Kc_tarsus, lp__);
            else
                b_tarsus = in__.vector_constrain(Kc_tarsus);

            current_statement_begin__ = 129;
            local_scalar_t__ Intercept_tarsus;
            (void) Intercept_tarsus;  // dummy to suppress unused var warning
            if (jacobian__)
                Intercept_tarsus = in__.scalar_constrain(lp__);
            else
                Intercept_tarsus = in__.scalar_constrain();

            current_statement_begin__ = 130;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> b_sigma_tarsus;
            (void) b_sigma_tarsus;  // dummy to suppress unused var warning
            if (jacobian__)
                b_sigma_tarsus = in__.vector_constrain(K_sigma_tarsus, lp__);
            else
                b_sigma_tarsus = in__.vector_constrain(K_sigma_tarsus);

            current_statement_begin__ = 131;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> Lrescor;
            (void) Lrescor;  // dummy to suppress unused var warning
            if (jacobian__)
                Lrescor = in__.cholesky_factor_corr_constrain(nresp, lp__);
            else
                Lrescor = in__.cholesky_factor_corr_constrain(nresp);

            current_statement_begin__ = 132;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> sd_1;
            (void) sd_1;  // dummy to suppress unused var warning
            if (jacobian__)
                sd_1 = in__.vector_lb_constrain(0, M_1, lp__);
            else
                sd_1 = in__.vector_lb_constrain(0, M_1);

            current_statement_begin__ = 133;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> z_1;
            (void) z_1;  // dummy to suppress unused var warning
            if (jacobian__)
                z_1 = in__.matrix_constrain(M_1, N_1, lp__);
            else
                z_1 = in__.matrix_constrain(M_1, N_1);

            current_statement_begin__ = 134;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> L_1;
            (void) L_1;  // dummy to suppress unused var warning
            if (jacobian__)
                L_1 = in__.cholesky_factor_corr_constrain(M_1, lp__);
            else
                L_1 = in__.cholesky_factor_corr_constrain(M_1);

            current_statement_begin__ = 135;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> sd_2;
            (void) sd_2;  // dummy to suppress unused var warning
            if (jacobian__)
                sd_2 = in__.vector_lb_constrain(0, M_2, lp__);
            else
                sd_2 = in__.vector_lb_constrain(0, M_2);

            current_statement_begin__ = 136;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> z_2;
            (void) z_2;  // dummy to suppress unused var warning
            if (jacobian__)
                z_2 = in__.matrix_constrain(M_2, N_2, lp__);
            else
                z_2 = in__.matrix_constrain(M_2, N_2);

            current_statement_begin__ = 137;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> L_2;
            (void) L_2;  // dummy to suppress unused var warning
            if (jacobian__)
                L_2 = in__.cholesky_factor_corr_constrain(M_2, lp__);
            else
                L_2 = in__.cholesky_factor_corr_constrain(M_2);

            current_statement_begin__ = 138;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> sd_3;
            (void) sd_3;  // dummy to suppress unused var warning
            if (jacobian__)
                sd_3 = in__.vector_lb_constrain(0, M_3, lp__);
            else
                sd_3 = in__.vector_lb_constrain(0, M_3);

            current_statement_begin__ = 139;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> z_3;
            (void) z_3;  // dummy to suppress unused var warning
            if (jacobian__)
                z_3 = in__.matrix_constrain(M_3, N_3, lp__);
            else
                z_3 = in__.matrix_constrain(M_3, N_3);

            current_statement_begin__ = 140;
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> L_3;
            (void) L_3;  // dummy to suppress unused var warning
            if (jacobian__)
                L_3 = in__.cholesky_factor_corr_constrain(M_3, lp__);
            else
                L_3 = in__.cholesky_factor_corr_constrain(M_3);

            // transformed parameters
            current_statement_begin__ = 143;
            validate_non_negative_index("r_1", "N_1", N_1);
            validate_non_negative_index("r_1", "M_1", M_1);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> r_1(N_1, M_1);
            stan::math::initialize(r_1, DUMMY_VAR__);
            stan::math::fill(r_1, DUMMY_VAR__);

            current_statement_begin__ = 145;
            validate_non_negative_index("r_1_bwt_1", "N_1", N_1);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_1_bwt_1(N_1);
            stan::math::initialize(r_1_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_1_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 146;
            validate_non_negative_index("r_1_bwt_2", "N_1", N_1);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_1_bwt_2(N_1);
            stan::math::initialize(r_1_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_1_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 147;
            validate_non_negative_index("r_1_tarsus_3", "N_1", N_1);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_1_tarsus_3(N_1);
            stan::math::initialize(r_1_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_1_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 148;
            validate_non_negative_index("r_1_tarsus_4", "N_1", N_1);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_1_tarsus_4(N_1);
            stan::math::initialize(r_1_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_1_tarsus_4, DUMMY_VAR__);

            current_statement_begin__ = 149;
            validate_non_negative_index("r_2", "N_2", N_2);
            validate_non_negative_index("r_2", "M_2", M_2);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> r_2(N_2, M_2);
            stan::math::initialize(r_2, DUMMY_VAR__);
            stan::math::fill(r_2, DUMMY_VAR__);

            current_statement_begin__ = 151;
            validate_non_negative_index("r_2_bwt_1", "N_2", N_2);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_2_bwt_1(N_2);
            stan::math::initialize(r_2_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_2_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 152;
            validate_non_negative_index("r_2_bwt_2", "N_2", N_2);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_2_bwt_2(N_2);
            stan::math::initialize(r_2_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_2_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 153;
            validate_non_negative_index("r_2_tarsus_3", "N_2", N_2);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_2_tarsus_3(N_2);
            stan::math::initialize(r_2_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_2_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 154;
            validate_non_negative_index("r_2_tarsus_4", "N_2", N_2);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_2_tarsus_4(N_2);
            stan::math::initialize(r_2_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_2_tarsus_4, DUMMY_VAR__);

            current_statement_begin__ = 155;
            validate_non_negative_index("r_3", "N_3", N_3);
            validate_non_negative_index("r_3", "M_3", M_3);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic> r_3(N_3, M_3);
            stan::math::initialize(r_3, DUMMY_VAR__);
            stan::math::fill(r_3, DUMMY_VAR__);

            current_statement_begin__ = 157;
            validate_non_negative_index("r_3_bwt_1", "N_3", N_3);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_3_bwt_1(N_3);
            stan::math::initialize(r_3_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_3_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 158;
            validate_non_negative_index("r_3_bwt_2", "N_3", N_3);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_3_bwt_2(N_3);
            stan::math::initialize(r_3_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_3_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 159;
            validate_non_negative_index("r_3_tarsus_3", "N_3", N_3);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_3_tarsus_3(N_3);
            stan::math::initialize(r_3_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_3_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 160;
            validate_non_negative_index("r_3_tarsus_4", "N_3", N_3);
            Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> r_3_tarsus_4(N_3);
            stan::math::initialize(r_3_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_3_tarsus_4, DUMMY_VAR__);

            // transformed parameters block statements
            current_statement_begin__ = 162;
            stan::math::assign(r_1, scale_r_cor_cov(z_1, sd_1, L_1, Lcov_1, pstream__));
            current_statement_begin__ = 163;
            stan::math::assign(r_1_bwt_1, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 164;
            stan::math::assign(r_1_bwt_2, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 165;
            stan::math::assign(r_1_tarsus_3, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 166;
            stan::math::assign(r_1_tarsus_4, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 168;
            stan::math::assign(r_2, scale_r_cor(z_2, sd_2, L_2, pstream__));
            current_statement_begin__ = 169;
            stan::math::assign(r_2_bwt_1, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 170;
            stan::math::assign(r_2_bwt_2, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 171;
            stan::math::assign(r_2_tarsus_3, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 172;
            stan::math::assign(r_2_tarsus_4, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 174;
            stan::math::assign(r_3, scale_r_cor(z_3, sd_3, L_3, pstream__));
            current_statement_begin__ = 175;
            stan::math::assign(r_3_bwt_1, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 176;
            stan::math::assign(r_3_bwt_2, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 177;
            stan::math::assign(r_3_tarsus_3, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 178;
            stan::math::assign(r_3_tarsus_4, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_3"));

            // validate transformed parameters
            const char* function__ = "validate transformed params";
            (void) function__;  // dummy to suppress unused var warning

            current_statement_begin__ = 143;
            size_t r_1_j_1_max__ = N_1;
            size_t r_1_j_2_max__ = M_1;
            for (size_t j_1__ = 0; j_1__ < r_1_j_1_max__; ++j_1__) {
                for (size_t j_2__ = 0; j_2__ < r_1_j_2_max__; ++j_2__) {
                    if (stan::math::is_uninitialized(r_1(j_1__, j_2__))) {
                        std::stringstream msg__;
                        msg__ << "Undefined transformed parameter: r_1" << "(" << j_1__ << ", " << j_2__ << ")";
                        stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_1: ") + msg__.str()), current_statement_begin__, prog_reader__());
                    }
                }
            }
            current_statement_begin__ = 145;
            size_t r_1_bwt_1_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_1_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_1_bwt_1(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_1_bwt_1" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_1_bwt_1: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 146;
            size_t r_1_bwt_2_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_2_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_1_bwt_2(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_1_bwt_2" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_1_bwt_2: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 147;
            size_t r_1_tarsus_3_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_3_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_1_tarsus_3(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_1_tarsus_3" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_1_tarsus_3: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 148;
            size_t r_1_tarsus_4_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_4_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_1_tarsus_4(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_1_tarsus_4" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_1_tarsus_4: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 149;
            size_t r_2_j_1_max__ = N_2;
            size_t r_2_j_2_max__ = M_2;
            for (size_t j_1__ = 0; j_1__ < r_2_j_1_max__; ++j_1__) {
                for (size_t j_2__ = 0; j_2__ < r_2_j_2_max__; ++j_2__) {
                    if (stan::math::is_uninitialized(r_2(j_1__, j_2__))) {
                        std::stringstream msg__;
                        msg__ << "Undefined transformed parameter: r_2" << "(" << j_1__ << ", " << j_2__ << ")";
                        stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_2: ") + msg__.str()), current_statement_begin__, prog_reader__());
                    }
                }
            }
            current_statement_begin__ = 151;
            size_t r_2_bwt_1_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_1_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_2_bwt_1(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_2_bwt_1" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_2_bwt_1: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 152;
            size_t r_2_bwt_2_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_2_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_2_bwt_2(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_2_bwt_2" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_2_bwt_2: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 153;
            size_t r_2_tarsus_3_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_3_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_2_tarsus_3(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_2_tarsus_3" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_2_tarsus_3: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 154;
            size_t r_2_tarsus_4_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_4_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_2_tarsus_4(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_2_tarsus_4" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_2_tarsus_4: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 155;
            size_t r_3_j_1_max__ = N_3;
            size_t r_3_j_2_max__ = M_3;
            for (size_t j_1__ = 0; j_1__ < r_3_j_1_max__; ++j_1__) {
                for (size_t j_2__ = 0; j_2__ < r_3_j_2_max__; ++j_2__) {
                    if (stan::math::is_uninitialized(r_3(j_1__, j_2__))) {
                        std::stringstream msg__;
                        msg__ << "Undefined transformed parameter: r_3" << "(" << j_1__ << ", " << j_2__ << ")";
                        stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_3: ") + msg__.str()), current_statement_begin__, prog_reader__());
                    }
                }
            }
            current_statement_begin__ = 157;
            size_t r_3_bwt_1_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_1_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_3_bwt_1(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_3_bwt_1" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_3_bwt_1: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 158;
            size_t r_3_bwt_2_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_2_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_3_bwt_2(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_3_bwt_2" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_3_bwt_2: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 159;
            size_t r_3_tarsus_3_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_3_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_3_tarsus_3(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_3_tarsus_3" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_3_tarsus_3: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }
            current_statement_begin__ = 160;
            size_t r_3_tarsus_4_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_4_j_1_max__; ++j_1__) {
                if (stan::math::is_uninitialized(r_3_tarsus_4(j_1__))) {
                    std::stringstream msg__;
                    msg__ << "Undefined transformed parameter: r_3_tarsus_4" << "(" << j_1__ << ")";
                    stan::lang::rethrow_located(std::runtime_error(std::string("Error initializing variable r_3_tarsus_4: ") + msg__.str()), current_statement_begin__, prog_reader__());
                }
            }

            // model body

            current_statement_begin__ = 182;
            if (as_bool(logical_negation(prior_only))) {
                {
                current_statement_begin__ = 184;
                validate_non_negative_index("mu_bwt", "N_bwt", N_bwt);
                Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> mu_bwt(N_bwt);
                stan::math::initialize(mu_bwt, DUMMY_VAR__);
                stan::math::fill(mu_bwt, DUMMY_VAR__);
                stan::math::assign(mu_bwt,add(Intercept_bwt, multiply(Xc_bwt, b_bwt)));

                current_statement_begin__ = 186;
                validate_non_negative_index("sigma_bwt", "N_bwt", N_bwt);
                Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> sigma_bwt(N_bwt);
                stan::math::initialize(sigma_bwt, DUMMY_VAR__);
                stan::math::fill(sigma_bwt, DUMMY_VAR__);
                stan::math::assign(sigma_bwt,multiply(X_sigma_bwt, b_sigma_bwt));

                current_statement_begin__ = 188;
                validate_non_negative_index("mu_tarsus", "N_tarsus", N_tarsus);
                Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> mu_tarsus(N_tarsus);
                stan::math::initialize(mu_tarsus, DUMMY_VAR__);
                stan::math::fill(mu_tarsus, DUMMY_VAR__);
                stan::math::assign(mu_tarsus,add(Intercept_tarsus, multiply(Xc_tarsus, b_tarsus)));

                current_statement_begin__ = 190;
                validate_non_negative_index("sigma_tarsus", "N_tarsus", N_tarsus);
                Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1> sigma_tarsus(N_tarsus);
                stan::math::initialize(sigma_tarsus, DUMMY_VAR__);
                stan::math::fill(sigma_tarsus, DUMMY_VAR__);
                stan::math::assign(sigma_tarsus,multiply(X_sigma_tarsus, b_sigma_tarsus));

                current_statement_begin__ = 192;
                validate_non_negative_index("Mu", "nresp", nresp);
                validate_non_negative_index("Mu", "N", N);
                std::vector<Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1>  > Mu(N, Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1>(nresp));
                stan::math::initialize(Mu, DUMMY_VAR__);
                stan::math::fill(Mu, DUMMY_VAR__);

                current_statement_begin__ = 193;
                validate_non_negative_index("sigma", "nresp", nresp);
                validate_non_negative_index("sigma", "N", N);
                std::vector<Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1>  > sigma(N, Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, 1>(nresp));
                stan::math::initialize(sigma, DUMMY_VAR__);
                stan::math::fill(sigma, DUMMY_VAR__);

                current_statement_begin__ = 195;
                validate_non_negative_index("LSigma", "nresp", nresp);
                validate_non_negative_index("LSigma", "nresp", nresp);
                validate_non_negative_index("LSigma", "N", N);
                std::vector<Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic>  > LSigma(N, Eigen::Matrix<local_scalar_t__, Eigen::Dynamic, Eigen::Dynamic>(nresp, nresp));
                stan::math::initialize(LSigma, DUMMY_VAR__);
                stan::math::fill(LSigma, DUMMY_VAR__);


                current_statement_begin__ = 196;
                for (int n = 1; n <= N_bwt; ++n) {

                    current_statement_begin__ = 198;
                    stan::model::assign(mu_bwt, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                (stan::model::rvalue(mu_bwt, stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), "mu_bwt") + ((((((get_base1(r_1_bwt_1, get_base1(J_1_bwt, n, "J_1_bwt", 1), "r_1_bwt_1", 1) * get_base1(Z_1_bwt_1, n, "Z_1_bwt_1", 1)) + (get_base1(r_1_bwt_2, get_base1(J_1_bwt, n, "J_1_bwt", 1), "r_1_bwt_2", 1) * get_base1(Z_1_bwt_2, n, "Z_1_bwt_2", 1))) + (get_base1(r_2_bwt_1, get_base1(J_2_bwt, n, "J_2_bwt", 1), "r_2_bwt_1", 1) * get_base1(Z_2_bwt_1, n, "Z_2_bwt_1", 1))) + (get_base1(r_2_bwt_2, get_base1(J_2_bwt, n, "J_2_bwt", 1), "r_2_bwt_2", 1) * get_base1(Z_2_bwt_2, n, "Z_2_bwt_2", 1))) + (get_base1(r_3_bwt_1, get_base1(J_3_bwt, n, "J_3_bwt", 1), "r_3_bwt_1", 1) * get_base1(Z_3_bwt_1, n, "Z_3_bwt_1", 1))) + (get_base1(r_3_bwt_2, get_base1(J_3_bwt, n, "J_3_bwt", 1), "r_3_bwt_2", 1) * get_base1(Z_3_bwt_2, n, "Z_3_bwt_2", 1)))), 
                                "assigning variable mu_bwt");
                }
                current_statement_begin__ = 200;
                for (int n = 1; n <= N_tarsus; ++n) {

                    current_statement_begin__ = 202;
                    stan::model::assign(mu_tarsus, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                (stan::model::rvalue(mu_tarsus, stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), "mu_tarsus") + ((((((get_base1(r_1_tarsus_3, get_base1(J_1_tarsus, n, "J_1_tarsus", 1), "r_1_tarsus_3", 1) * get_base1(Z_1_tarsus_3, n, "Z_1_tarsus_3", 1)) + (get_base1(r_1_tarsus_4, get_base1(J_1_tarsus, n, "J_1_tarsus", 1), "r_1_tarsus_4", 1) * get_base1(Z_1_tarsus_4, n, "Z_1_tarsus_4", 1))) + (get_base1(r_2_tarsus_3, get_base1(J_2_tarsus, n, "J_2_tarsus", 1), "r_2_tarsus_3", 1) * get_base1(Z_2_tarsus_3, n, "Z_2_tarsus_3", 1))) + (get_base1(r_2_tarsus_4, get_base1(J_2_tarsus, n, "J_2_tarsus", 1), "r_2_tarsus_4", 1) * get_base1(Z_2_tarsus_4, n, "Z_2_tarsus_4", 1))) + (get_base1(r_3_tarsus_3, get_base1(J_3_tarsus, n, "J_3_tarsus", 1), "r_3_tarsus_3", 1) * get_base1(Z_3_tarsus_3, n, "Z_3_tarsus_3", 1))) + (get_base1(r_3_tarsus_4, get_base1(J_3_tarsus, n, "J_3_tarsus", 1), "r_3_tarsus_4", 1) * get_base1(Z_3_tarsus_4, n, "Z_3_tarsus_4", 1)))), 
                                "assigning variable mu_tarsus");
                }
                current_statement_begin__ = 204;
                for (int n = 1; n <= N_bwt; ++n) {

                    current_statement_begin__ = 206;
                    stan::model::assign(sigma_bwt, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                stan::math::exp(get_base1(sigma_bwt, n, "sigma_bwt", 1)), 
                                "assigning variable sigma_bwt");
                }
                current_statement_begin__ = 208;
                for (int n = 1; n <= N_tarsus; ++n) {

                    current_statement_begin__ = 210;
                    stan::model::assign(sigma_tarsus, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                stan::math::exp(get_base1(sigma_tarsus, n, "sigma_tarsus", 1)), 
                                "assigning variable sigma_tarsus");
                }
                current_statement_begin__ = 213;
                for (int n = 1; n <= N; ++n) {

                    current_statement_begin__ = 214;
                    stan::model::assign(Mu, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                transpose(stan::math::to_row_vector(stan::math::array_builder<local_scalar_t__ >().add(get_base1(mu_bwt, n, "mu_bwt", 1)).add(get_base1(mu_tarsus, n, "mu_tarsus", 1)).array())), 
                                "assigning variable Mu");
                    current_statement_begin__ = 215;
                    stan::model::assign(sigma, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                transpose(stan::math::to_row_vector(stan::math::array_builder<local_scalar_t__ >().add(get_base1(sigma_bwt, n, "sigma_bwt", 1)).add(get_base1(sigma_tarsus, n, "sigma_tarsus", 1)).array())), 
                                "assigning variable sigma");
                    current_statement_begin__ = 216;
                    stan::model::assign(LSigma, 
                                stan::model::cons_list(stan::model::index_uni(n), stan::model::nil_index_list()), 
                                diag_pre_multiply(get_base1(sigma, n, "sigma", 1), Lrescor), 
                                "assigning variable LSigma");
                }
                current_statement_begin__ = 218;
                for (int n = 1; n <= N; ++n) {

                    current_statement_begin__ = 219;
                    lp_accum__.add(multi_normal_cholesky_log(get_base1(Y, n, "Y", 1), get_base1(Mu, n, "Mu", 1), get_base1(LSigma, n, "LSigma", 1)));
                }
                }
            }
            current_statement_begin__ = 223;
            lp_accum__.add(student_t_log(Intercept_bwt, 3, 7.3, 2.8));
            current_statement_begin__ = 224;
            lp_accum__.add(student_t_log(Intercept_tarsus, 3, 20.4, 5.2));
            current_statement_begin__ = 225;
            lp_accum__.add(lkj_corr_cholesky_log(Lrescor, 1));
            current_statement_begin__ = 226;
            lp_accum__.add((student_t_log(get_base1(sd_1, 1, "sd_1", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 228;
            lp_accum__.add((student_t_log(get_base1(sd_1, 2, "sd_1", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 230;
            lp_accum__.add((student_t_log(get_base1(sd_1, 3, "sd_1", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 232;
            lp_accum__.add((student_t_log(get_base1(sd_1, 4, "sd_1", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 234;
            lp_accum__.add(std_normal_log(to_vector(z_1)));
            current_statement_begin__ = 235;
            lp_accum__.add(lkj_corr_cholesky_log(L_1, 1));
            current_statement_begin__ = 236;
            lp_accum__.add((student_t_log(get_base1(sd_2, 1, "sd_2", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 238;
            lp_accum__.add((student_t_log(get_base1(sd_2, 2, "sd_2", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 240;
            lp_accum__.add((student_t_log(get_base1(sd_2, 3, "sd_2", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 242;
            lp_accum__.add((student_t_log(get_base1(sd_2, 4, "sd_2", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 244;
            lp_accum__.add(std_normal_log(to_vector(z_2)));
            current_statement_begin__ = 245;
            lp_accum__.add(lkj_corr_cholesky_log(L_2, 1));
            current_statement_begin__ = 246;
            lp_accum__.add((student_t_log(get_base1(sd_3, 1, "sd_3", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 248;
            lp_accum__.add((student_t_log(get_base1(sd_3, 2, "sd_3", 1), 3, 0, 2.8) - (1 * student_t_ccdf_log(0, 3, 0, 2.8))));
            current_statement_begin__ = 250;
            lp_accum__.add((student_t_log(get_base1(sd_3, 3, "sd_3", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 252;
            lp_accum__.add((student_t_log(get_base1(sd_3, 4, "sd_3", 1), 3, 0, 5.2) - (1 * student_t_ccdf_log(0, 3, 0, 5.2))));
            current_statement_begin__ = 254;
            lp_accum__.add(std_normal_log(to_vector(z_3)));
            current_statement_begin__ = 255;
            lp_accum__.add(lkj_corr_cholesky_log(L_3, 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("b_bwt");
        names__.push_back("Intercept_bwt");
        names__.push_back("b_sigma_bwt");
        names__.push_back("b_tarsus");
        names__.push_back("Intercept_tarsus");
        names__.push_back("b_sigma_tarsus");
        names__.push_back("Lrescor");
        names__.push_back("sd_1");
        names__.push_back("z_1");
        names__.push_back("L_1");
        names__.push_back("sd_2");
        names__.push_back("z_2");
        names__.push_back("L_2");
        names__.push_back("sd_3");
        names__.push_back("z_3");
        names__.push_back("L_3");
        names__.push_back("r_1");
        names__.push_back("r_1_bwt_1");
        names__.push_back("r_1_bwt_2");
        names__.push_back("r_1_tarsus_3");
        names__.push_back("r_1_tarsus_4");
        names__.push_back("r_2");
        names__.push_back("r_2_bwt_1");
        names__.push_back("r_2_bwt_2");
        names__.push_back("r_2_tarsus_3");
        names__.push_back("r_2_tarsus_4");
        names__.push_back("r_3");
        names__.push_back("r_3_bwt_1");
        names__.push_back("r_3_bwt_2");
        names__.push_back("r_3_tarsus_3");
        names__.push_back("r_3_tarsus_4");
        names__.push_back("b_bwt_Intercept");
        names__.push_back("b_tarsus_Intercept");
        names__.push_back("Rescor");
        names__.push_back("rescor");
        names__.push_back("Cor_1");
        names__.push_back("cor_1");
        names__.push_back("Cor_2");
        names__.push_back("cor_2");
        names__.push_back("Cor_3");
        names__.push_back("cor_3");
    }


    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(Kc_bwt);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(K_sigma_bwt);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(Kc_tarsus);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(K_sigma_tarsus);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(nresp);
        dims__.push_back(nresp);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_1);
        dims__.push_back(N_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_1);
        dims__.push_back(M_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_2);
        dims__.push_back(N_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_2);
        dims__.push_back(M_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_3);
        dims__.push_back(N_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_3);
        dims__.push_back(M_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_1);
        dims__.push_back(M_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_2);
        dims__.push_back(M_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_3);
        dims__.push_back(M_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(N_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(nresp);
        dims__.push_back(nresp);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(nrescor);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_1);
        dims__.push_back(M_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(NC_1);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_2);
        dims__.push_back(M_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(NC_2);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(M_3);
        dims__.push_back(M_3);
        dimss__.push_back(dims__);
        dims__.resize(0);
        dims__.push_back(NC_3);
        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__ = "model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb_namespace::write_array";
        (void) function__;  // dummy to suppress unused var warning

        // read-transform, write parameters
        Eigen::Matrix<double, Eigen::Dynamic, 1> b_bwt = in__.vector_constrain(Kc_bwt);
        size_t b_bwt_j_1_max__ = Kc_bwt;
        for (size_t j_1__ = 0; j_1__ < b_bwt_j_1_max__; ++j_1__) {
            vars__.push_back(b_bwt(j_1__));
        }

        double Intercept_bwt = in__.scalar_constrain();
        vars__.push_back(Intercept_bwt);

        Eigen::Matrix<double, Eigen::Dynamic, 1> b_sigma_bwt = in__.vector_constrain(K_sigma_bwt);
        size_t b_sigma_bwt_j_1_max__ = K_sigma_bwt;
        for (size_t j_1__ = 0; j_1__ < b_sigma_bwt_j_1_max__; ++j_1__) {
            vars__.push_back(b_sigma_bwt(j_1__));
        }

        Eigen::Matrix<double, Eigen::Dynamic, 1> b_tarsus = in__.vector_constrain(Kc_tarsus);
        size_t b_tarsus_j_1_max__ = Kc_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_tarsus_j_1_max__; ++j_1__) {
            vars__.push_back(b_tarsus(j_1__));
        }

        double Intercept_tarsus = in__.scalar_constrain();
        vars__.push_back(Intercept_tarsus);

        Eigen::Matrix<double, Eigen::Dynamic, 1> b_sigma_tarsus = in__.vector_constrain(K_sigma_tarsus);
        size_t b_sigma_tarsus_j_1_max__ = K_sigma_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_sigma_tarsus_j_1_max__; ++j_1__) {
            vars__.push_back(b_sigma_tarsus(j_1__));
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Lrescor = in__.cholesky_factor_corr_constrain(nresp);
        size_t Lrescor_j_2_max__ = nresp;
        size_t Lrescor_j_1_max__ = nresp;
        for (size_t j_2__ = 0; j_2__ < Lrescor_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Lrescor_j_1_max__; ++j_1__) {
                vars__.push_back(Lrescor(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_1 = in__.vector_lb_constrain(0, M_1);
        size_t sd_1_j_1_max__ = M_1;
        for (size_t j_1__ = 0; j_1__ < sd_1_j_1_max__; ++j_1__) {
            vars__.push_back(sd_1(j_1__));
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_1 = in__.matrix_constrain(M_1, N_1);
        size_t z_1_j_2_max__ = N_1;
        size_t z_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < z_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_1_j_1_max__; ++j_1__) {
                vars__.push_back(z_1(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_1 = in__.cholesky_factor_corr_constrain(M_1);
        size_t L_1_j_2_max__ = M_1;
        size_t L_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < L_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_1_j_1_max__; ++j_1__) {
                vars__.push_back(L_1(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_2 = in__.vector_lb_constrain(0, M_2);
        size_t sd_2_j_1_max__ = M_2;
        for (size_t j_1__ = 0; j_1__ < sd_2_j_1_max__; ++j_1__) {
            vars__.push_back(sd_2(j_1__));
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_2 = in__.matrix_constrain(M_2, N_2);
        size_t z_2_j_2_max__ = N_2;
        size_t z_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < z_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_2_j_1_max__; ++j_1__) {
                vars__.push_back(z_2(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_2 = in__.cholesky_factor_corr_constrain(M_2);
        size_t L_2_j_2_max__ = M_2;
        size_t L_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < L_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_2_j_1_max__; ++j_1__) {
                vars__.push_back(L_2(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, 1> sd_3 = in__.vector_lb_constrain(0, M_3);
        size_t sd_3_j_1_max__ = M_3;
        for (size_t j_1__ = 0; j_1__ < sd_3_j_1_max__; ++j_1__) {
            vars__.push_back(sd_3(j_1__));
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> z_3 = in__.matrix_constrain(M_3, N_3);
        size_t z_3_j_2_max__ = N_3;
        size_t z_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < z_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_3_j_1_max__; ++j_1__) {
                vars__.push_back(z_3(j_1__, j_2__));
            }
        }

        Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> L_3 = in__.cholesky_factor_corr_constrain(M_3);
        size_t L_3_j_2_max__ = M_3;
        size_t L_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < L_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_3_j_1_max__; ++j_1__) {
                vars__.push_back(L_3(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__ = 143;
            validate_non_negative_index("r_1", "N_1", N_1);
            validate_non_negative_index("r_1", "M_1", M_1);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> r_1(N_1, M_1);
            stan::math::initialize(r_1, DUMMY_VAR__);
            stan::math::fill(r_1, DUMMY_VAR__);

            current_statement_begin__ = 145;
            validate_non_negative_index("r_1_bwt_1", "N_1", N_1);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_1_bwt_1(N_1);
            stan::math::initialize(r_1_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_1_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 146;
            validate_non_negative_index("r_1_bwt_2", "N_1", N_1);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_1_bwt_2(N_1);
            stan::math::initialize(r_1_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_1_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 147;
            validate_non_negative_index("r_1_tarsus_3", "N_1", N_1);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_1_tarsus_3(N_1);
            stan::math::initialize(r_1_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_1_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 148;
            validate_non_negative_index("r_1_tarsus_4", "N_1", N_1);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_1_tarsus_4(N_1);
            stan::math::initialize(r_1_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_1_tarsus_4, DUMMY_VAR__);

            current_statement_begin__ = 149;
            validate_non_negative_index("r_2", "N_2", N_2);
            validate_non_negative_index("r_2", "M_2", M_2);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> r_2(N_2, M_2);
            stan::math::initialize(r_2, DUMMY_VAR__);
            stan::math::fill(r_2, DUMMY_VAR__);

            current_statement_begin__ = 151;
            validate_non_negative_index("r_2_bwt_1", "N_2", N_2);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_2_bwt_1(N_2);
            stan::math::initialize(r_2_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_2_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 152;
            validate_non_negative_index("r_2_bwt_2", "N_2", N_2);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_2_bwt_2(N_2);
            stan::math::initialize(r_2_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_2_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 153;
            validate_non_negative_index("r_2_tarsus_3", "N_2", N_2);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_2_tarsus_3(N_2);
            stan::math::initialize(r_2_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_2_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 154;
            validate_non_negative_index("r_2_tarsus_4", "N_2", N_2);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_2_tarsus_4(N_2);
            stan::math::initialize(r_2_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_2_tarsus_4, DUMMY_VAR__);

            current_statement_begin__ = 155;
            validate_non_negative_index("r_3", "N_3", N_3);
            validate_non_negative_index("r_3", "M_3", M_3);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> r_3(N_3, M_3);
            stan::math::initialize(r_3, DUMMY_VAR__);
            stan::math::fill(r_3, DUMMY_VAR__);

            current_statement_begin__ = 157;
            validate_non_negative_index("r_3_bwt_1", "N_3", N_3);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_3_bwt_1(N_3);
            stan::math::initialize(r_3_bwt_1, DUMMY_VAR__);
            stan::math::fill(r_3_bwt_1, DUMMY_VAR__);

            current_statement_begin__ = 158;
            validate_non_negative_index("r_3_bwt_2", "N_3", N_3);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_3_bwt_2(N_3);
            stan::math::initialize(r_3_bwt_2, DUMMY_VAR__);
            stan::math::fill(r_3_bwt_2, DUMMY_VAR__);

            current_statement_begin__ = 159;
            validate_non_negative_index("r_3_tarsus_3", "N_3", N_3);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_3_tarsus_3(N_3);
            stan::math::initialize(r_3_tarsus_3, DUMMY_VAR__);
            stan::math::fill(r_3_tarsus_3, DUMMY_VAR__);

            current_statement_begin__ = 160;
            validate_non_negative_index("r_3_tarsus_4", "N_3", N_3);
            Eigen::Matrix<double, Eigen::Dynamic, 1> r_3_tarsus_4(N_3);
            stan::math::initialize(r_3_tarsus_4, DUMMY_VAR__);
            stan::math::fill(r_3_tarsus_4, DUMMY_VAR__);

            // do transformed parameters statements
            current_statement_begin__ = 162;
            stan::math::assign(r_1, scale_r_cor_cov(z_1, sd_1, L_1, Lcov_1, pstream__));
            current_statement_begin__ = 163;
            stan::math::assign(r_1_bwt_1, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 164;
            stan::math::assign(r_1_bwt_2, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 165;
            stan::math::assign(r_1_tarsus_3, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 166;
            stan::math::assign(r_1_tarsus_4, stan::model::rvalue(r_1, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_1"));
            current_statement_begin__ = 168;
            stan::math::assign(r_2, scale_r_cor(z_2, sd_2, L_2, pstream__));
            current_statement_begin__ = 169;
            stan::math::assign(r_2_bwt_1, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 170;
            stan::math::assign(r_2_bwt_2, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 171;
            stan::math::assign(r_2_tarsus_3, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 172;
            stan::math::assign(r_2_tarsus_4, stan::model::rvalue(r_2, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_2"));
            current_statement_begin__ = 174;
            stan::math::assign(r_3, scale_r_cor(z_3, sd_3, L_3, pstream__));
            current_statement_begin__ = 175;
            stan::math::assign(r_3_bwt_1, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(1), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 176;
            stan::math::assign(r_3_bwt_2, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(2), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 177;
            stan::math::assign(r_3_tarsus_3, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(3), stan::model::nil_index_list())), "r_3"));
            current_statement_begin__ = 178;
            stan::math::assign(r_3_tarsus_4, stan::model::rvalue(r_3, stan::model::cons_list(stan::model::index_omni(), stan::model::cons_list(stan::model::index_uni(4), stan::model::nil_index_list())), "r_3"));

            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 r_1_j_2_max__ = M_1;
                size_t r_1_j_1_max__ = N_1;
                for (size_t j_2__ = 0; j_2__ < r_1_j_2_max__; ++j_2__) {
                    for (size_t j_1__ = 0; j_1__ < r_1_j_1_max__; ++j_1__) {
                        vars__.push_back(r_1(j_1__, j_2__));
                    }
                }
                size_t r_1_bwt_1_j_1_max__ = N_1;
                for (size_t j_1__ = 0; j_1__ < r_1_bwt_1_j_1_max__; ++j_1__) {
                    vars__.push_back(r_1_bwt_1(j_1__));
                }
                size_t r_1_bwt_2_j_1_max__ = N_1;
                for (size_t j_1__ = 0; j_1__ < r_1_bwt_2_j_1_max__; ++j_1__) {
                    vars__.push_back(r_1_bwt_2(j_1__));
                }
                size_t r_1_tarsus_3_j_1_max__ = N_1;
                for (size_t j_1__ = 0; j_1__ < r_1_tarsus_3_j_1_max__; ++j_1__) {
                    vars__.push_back(r_1_tarsus_3(j_1__));
                }
                size_t r_1_tarsus_4_j_1_max__ = N_1;
                for (size_t j_1__ = 0; j_1__ < r_1_tarsus_4_j_1_max__; ++j_1__) {
                    vars__.push_back(r_1_tarsus_4(j_1__));
                }
                size_t r_2_j_2_max__ = M_2;
                size_t r_2_j_1_max__ = N_2;
                for (size_t j_2__ = 0; j_2__ < r_2_j_2_max__; ++j_2__) {
                    for (size_t j_1__ = 0; j_1__ < r_2_j_1_max__; ++j_1__) {
                        vars__.push_back(r_2(j_1__, j_2__));
                    }
                }
                size_t r_2_bwt_1_j_1_max__ = N_2;
                for (size_t j_1__ = 0; j_1__ < r_2_bwt_1_j_1_max__; ++j_1__) {
                    vars__.push_back(r_2_bwt_1(j_1__));
                }
                size_t r_2_bwt_2_j_1_max__ = N_2;
                for (size_t j_1__ = 0; j_1__ < r_2_bwt_2_j_1_max__; ++j_1__) {
                    vars__.push_back(r_2_bwt_2(j_1__));
                }
                size_t r_2_tarsus_3_j_1_max__ = N_2;
                for (size_t j_1__ = 0; j_1__ < r_2_tarsus_3_j_1_max__; ++j_1__) {
                    vars__.push_back(r_2_tarsus_3(j_1__));
                }
                size_t r_2_tarsus_4_j_1_max__ = N_2;
                for (size_t j_1__ = 0; j_1__ < r_2_tarsus_4_j_1_max__; ++j_1__) {
                    vars__.push_back(r_2_tarsus_4(j_1__));
                }
                size_t r_3_j_2_max__ = M_3;
                size_t r_3_j_1_max__ = N_3;
                for (size_t j_2__ = 0; j_2__ < r_3_j_2_max__; ++j_2__) {
                    for (size_t j_1__ = 0; j_1__ < r_3_j_1_max__; ++j_1__) {
                        vars__.push_back(r_3(j_1__, j_2__));
                    }
                }
                size_t r_3_bwt_1_j_1_max__ = N_3;
                for (size_t j_1__ = 0; j_1__ < r_3_bwt_1_j_1_max__; ++j_1__) {
                    vars__.push_back(r_3_bwt_1(j_1__));
                }
                size_t r_3_bwt_2_j_1_max__ = N_3;
                for (size_t j_1__ = 0; j_1__ < r_3_bwt_2_j_1_max__; ++j_1__) {
                    vars__.push_back(r_3_bwt_2(j_1__));
                }
                size_t r_3_tarsus_3_j_1_max__ = N_3;
                for (size_t j_1__ = 0; j_1__ < r_3_tarsus_3_j_1_max__; ++j_1__) {
                    vars__.push_back(r_3_tarsus_3(j_1__));
                }
                size_t r_3_tarsus_4_j_1_max__ = N_3;
                for (size_t j_1__ = 0; j_1__ < r_3_tarsus_4_j_1_max__; ++j_1__) {
                    vars__.push_back(r_3_tarsus_4(j_1__));
                }
            }
            if (!include_gqs__) return;
            // declare and define generated quantities
            current_statement_begin__ = 259;
            double b_bwt_Intercept;
            (void) b_bwt_Intercept;  // dummy to suppress unused var warning
            stan::math::initialize(b_bwt_Intercept, DUMMY_VAR__);
            stan::math::fill(b_bwt_Intercept, DUMMY_VAR__);
            stan::math::assign(b_bwt_Intercept,(Intercept_bwt - dot_product(means_X_bwt, b_bwt)));

            current_statement_begin__ = 261;
            double b_tarsus_Intercept;
            (void) b_tarsus_Intercept;  // dummy to suppress unused var warning
            stan::math::initialize(b_tarsus_Intercept, DUMMY_VAR__);
            stan::math::fill(b_tarsus_Intercept, DUMMY_VAR__);
            stan::math::assign(b_tarsus_Intercept,(Intercept_tarsus - dot_product(means_X_tarsus, b_tarsus)));

            current_statement_begin__ = 263;
            validate_non_negative_index("Rescor", "nresp", nresp);
            validate_non_negative_index("Rescor", "nresp", nresp);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Rescor(nresp, nresp);
            stan::math::initialize(Rescor, DUMMY_VAR__);
            stan::math::fill(Rescor, DUMMY_VAR__);
            stan::math::assign(Rescor,multiply_lower_tri_self_transpose(Lrescor));

            current_statement_begin__ = 264;
            validate_non_negative_index("rescor", "nrescor", nrescor);
            Eigen::Matrix<double, Eigen::Dynamic, 1> rescor(nrescor);
            stan::math::initialize(rescor, DUMMY_VAR__);
            stan::math::fill(rescor, DUMMY_VAR__);

            current_statement_begin__ = 266;
            validate_non_negative_index("Cor_1", "M_1", M_1);
            validate_non_negative_index("Cor_1", "M_1", M_1);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Cor_1(M_1, M_1);
            stan::math::initialize(Cor_1, DUMMY_VAR__);
            stan::math::fill(Cor_1, DUMMY_VAR__);
            stan::math::assign(Cor_1,multiply_lower_tri_self_transpose(L_1));

            current_statement_begin__ = 267;
            validate_non_negative_index("cor_1", "NC_1", NC_1);
            Eigen::Matrix<double, Eigen::Dynamic, 1> cor_1(NC_1);
            stan::math::initialize(cor_1, DUMMY_VAR__);
            stan::math::fill(cor_1, DUMMY_VAR__);

            current_statement_begin__ = 269;
            validate_non_negative_index("Cor_2", "M_2", M_2);
            validate_non_negative_index("Cor_2", "M_2", M_2);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Cor_2(M_2, M_2);
            stan::math::initialize(Cor_2, DUMMY_VAR__);
            stan::math::fill(Cor_2, DUMMY_VAR__);
            stan::math::assign(Cor_2,multiply_lower_tri_self_transpose(L_2));

            current_statement_begin__ = 270;
            validate_non_negative_index("cor_2", "NC_2", NC_2);
            Eigen::Matrix<double, Eigen::Dynamic, 1> cor_2(NC_2);
            stan::math::initialize(cor_2, DUMMY_VAR__);
            stan::math::fill(cor_2, DUMMY_VAR__);

            current_statement_begin__ = 272;
            validate_non_negative_index("Cor_3", "M_3", M_3);
            validate_non_negative_index("Cor_3", "M_3", M_3);
            Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> Cor_3(M_3, M_3);
            stan::math::initialize(Cor_3, DUMMY_VAR__);
            stan::math::fill(Cor_3, DUMMY_VAR__);
            stan::math::assign(Cor_3,multiply_lower_tri_self_transpose(L_3));

            current_statement_begin__ = 273;
            validate_non_negative_index("cor_3", "NC_3", NC_3);
            Eigen::Matrix<double, Eigen::Dynamic, 1> cor_3(NC_3);
            stan::math::initialize(cor_3, DUMMY_VAR__);
            stan::math::fill(cor_3, DUMMY_VAR__);

            // generated quantities statements
            current_statement_begin__ = 275;
            for (int k = 1; k <= nresp; ++k) {

                current_statement_begin__ = 276;
                for (int j = 1; j <= (k - 1); ++j) {

                    current_statement_begin__ = 277;
                    stan::model::assign(rescor, 
                                stan::model::cons_list(stan::model::index_uni((choose((k - 1), 2) + j)), stan::model::nil_index_list()), 
                                get_base1(Rescor, j, k, "Rescor", 1), 
                                "assigning variable rescor");
                }
            }
            current_statement_begin__ = 281;
            for (int k = 1; k <= M_1; ++k) {

                current_statement_begin__ = 282;
                for (int j = 1; j <= (k - 1); ++j) {

                    current_statement_begin__ = 283;
                    stan::model::assign(cor_1, 
                                stan::model::cons_list(stan::model::index_uni((choose((k - 1), 2) + j)), stan::model::nil_index_list()), 
                                get_base1(Cor_1, j, k, "Cor_1", 1), 
                                "assigning variable cor_1");
                }
            }
            current_statement_begin__ = 287;
            for (int k = 1; k <= M_2; ++k) {

                current_statement_begin__ = 288;
                for (int j = 1; j <= (k - 1); ++j) {

                    current_statement_begin__ = 289;
                    stan::model::assign(cor_2, 
                                stan::model::cons_list(stan::model::index_uni((choose((k - 1), 2) + j)), stan::model::nil_index_list()), 
                                get_base1(Cor_2, j, k, "Cor_2", 1), 
                                "assigning variable cor_2");
                }
            }
            current_statement_begin__ = 293;
            for (int k = 1; k <= M_3; ++k) {

                current_statement_begin__ = 294;
                for (int j = 1; j <= (k - 1); ++j) {

                    current_statement_begin__ = 295;
                    stan::model::assign(cor_3, 
                                stan::model::cons_list(stan::model::index_uni((choose((k - 1), 2) + j)), stan::model::nil_index_list()), 
                                get_base1(Cor_3, j, k, "Cor_3", 1), 
                                "assigning variable cor_3");
                }
            }

            // validate, write generated quantities
            current_statement_begin__ = 259;
            vars__.push_back(b_bwt_Intercept);

            current_statement_begin__ = 261;
            vars__.push_back(b_tarsus_Intercept);

            current_statement_begin__ = 263;
            stan::math::check_corr_matrix(function__, "Rescor", Rescor);

            size_t Rescor_j_2_max__ = nresp;
            size_t Rescor_j_1_max__ = nresp;
            for (size_t j_2__ = 0; j_2__ < Rescor_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < Rescor_j_1_max__; ++j_1__) {
                    vars__.push_back(Rescor(j_1__, j_2__));
                }
            }

            current_statement_begin__ = 264;
            check_greater_or_equal(function__, "rescor", rescor, -(1));
            check_less_or_equal(function__, "rescor", rescor, 1);

            size_t rescor_j_1_max__ = nrescor;
            for (size_t j_1__ = 0; j_1__ < rescor_j_1_max__; ++j_1__) {
                vars__.push_back(rescor(j_1__));
            }

            current_statement_begin__ = 266;
            stan::math::check_corr_matrix(function__, "Cor_1", Cor_1);

            size_t Cor_1_j_2_max__ = M_1;
            size_t Cor_1_j_1_max__ = M_1;
            for (size_t j_2__ = 0; j_2__ < Cor_1_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < Cor_1_j_1_max__; ++j_1__) {
                    vars__.push_back(Cor_1(j_1__, j_2__));
                }
            }

            current_statement_begin__ = 267;
            check_greater_or_equal(function__, "cor_1", cor_1, -(1));
            check_less_or_equal(function__, "cor_1", cor_1, 1);

            size_t cor_1_j_1_max__ = NC_1;
            for (size_t j_1__ = 0; j_1__ < cor_1_j_1_max__; ++j_1__) {
                vars__.push_back(cor_1(j_1__));
            }

            current_statement_begin__ = 269;
            stan::math::check_corr_matrix(function__, "Cor_2", Cor_2);

            size_t Cor_2_j_2_max__ = M_2;
            size_t Cor_2_j_1_max__ = M_2;
            for (size_t j_2__ = 0; j_2__ < Cor_2_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < Cor_2_j_1_max__; ++j_1__) {
                    vars__.push_back(Cor_2(j_1__, j_2__));
                }
            }

            current_statement_begin__ = 270;
            check_greater_or_equal(function__, "cor_2", cor_2, -(1));
            check_less_or_equal(function__, "cor_2", cor_2, 1);

            size_t cor_2_j_1_max__ = NC_2;
            for (size_t j_1__ = 0; j_1__ < cor_2_j_1_max__; ++j_1__) {
                vars__.push_back(cor_2(j_1__));
            }

            current_statement_begin__ = 272;
            stan::math::check_corr_matrix(function__, "Cor_3", Cor_3);

            size_t Cor_3_j_2_max__ = M_3;
            size_t Cor_3_j_1_max__ = M_3;
            for (size_t j_2__ = 0; j_2__ < Cor_3_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < Cor_3_j_1_max__; ++j_1__) {
                    vars__.push_back(Cor_3(j_1__, j_2__));
                }
            }

            current_statement_begin__ = 273;
            check_greater_or_equal(function__, "cor_3", cor_3, -(1));
            check_less_or_equal(function__, "cor_3", cor_3, 1);

            size_t cor_3_j_1_max__ = NC_3;
            for (size_t j_1__ = 0; j_1__ < cor_3_j_1_max__; ++j_1__) {
                vars__.push_back(cor_3(j_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 ***");
        }
    }

    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 "model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb";
    }


    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 b_bwt_j_1_max__ = Kc_bwt;
        for (size_t j_1__ = 0; j_1__ < b_bwt_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_bwt" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        param_name_stream__.str(std::string());
        param_name_stream__ << "Intercept_bwt";
        param_names__.push_back(param_name_stream__.str());
        size_t b_sigma_bwt_j_1_max__ = K_sigma_bwt;
        for (size_t j_1__ = 0; j_1__ < b_sigma_bwt_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_sigma_bwt" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t b_tarsus_j_1_max__ = Kc_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_tarsus_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_tarsus" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        param_name_stream__.str(std::string());
        param_name_stream__ << "Intercept_tarsus";
        param_names__.push_back(param_name_stream__.str());
        size_t b_sigma_tarsus_j_1_max__ = K_sigma_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_sigma_tarsus_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_sigma_tarsus" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Lrescor_j_2_max__ = nresp;
        size_t Lrescor_j_1_max__ = nresp;
        for (size_t j_2__ = 0; j_2__ < Lrescor_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Lrescor_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "Lrescor" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t sd_1_j_1_max__ = M_1;
        for (size_t j_1__ = 0; j_1__ < sd_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_1_j_2_max__ = N_1;
        size_t z_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < z_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_1_j_2_max__ = M_1;
        size_t L_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < L_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "L_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t sd_2_j_1_max__ = M_2;
        for (size_t j_1__ = 0; j_1__ < sd_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_2_j_2_max__ = N_2;
        size_t z_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < z_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_2_j_2_max__ = M_2;
        size_t L_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < L_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "L_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t sd_3_j_1_max__ = M_3;
        for (size_t j_1__ = 0; j_1__ < sd_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_3_j_2_max__ = N_3;
        size_t z_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < z_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_3" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_3_j_2_max__ = M_3;
        size_t L_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < L_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < L_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "L_3" << '.' << 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 r_1_j_2_max__ = M_1;
            size_t r_1_j_1_max__ = N_1;
            for (size_t j_2__ = 0; j_2__ < r_1_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_1_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_1_bwt_1_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_bwt_2_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_tarsus_3_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_tarsus_4_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_j_2_max__ = M_2;
            size_t r_2_j_1_max__ = N_2;
            for (size_t j_2__ = 0; j_2__ < r_2_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_2_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_2_bwt_1_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_bwt_2_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_tarsus_3_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_tarsus_4_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_j_2_max__ = M_3;
            size_t r_3_j_1_max__ = N_3;
            for (size_t j_2__ = 0; j_2__ < r_3_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_3_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_3" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_3_bwt_1_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_bwt_2_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_tarsus_3_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_tarsus_4_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }

        if (!include_gqs__) return;
        param_name_stream__.str(std::string());
        param_name_stream__ << "b_bwt_Intercept";
        param_names__.push_back(param_name_stream__.str());
        param_name_stream__.str(std::string());
        param_name_stream__ << "b_tarsus_Intercept";
        param_names__.push_back(param_name_stream__.str());
        size_t Rescor_j_2_max__ = nresp;
        size_t Rescor_j_1_max__ = nresp;
        for (size_t j_2__ = 0; j_2__ < Rescor_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Rescor_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "Rescor" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t rescor_j_1_max__ = nrescor;
        for (size_t j_1__ = 0; j_1__ < rescor_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "rescor" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_1_j_2_max__ = M_1;
        size_t Cor_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < Cor_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Cor_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "Cor_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t cor_1_j_1_max__ = NC_1;
        for (size_t j_1__ = 0; j_1__ < cor_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_2_j_2_max__ = M_2;
        size_t Cor_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < Cor_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Cor_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "Cor_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t cor_2_j_1_max__ = NC_2;
        for (size_t j_1__ = 0; j_1__ < cor_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_3_j_2_max__ = M_3;
        size_t Cor_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < Cor_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < Cor_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "Cor_3" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t cor_3_j_1_max__ = NC_3;
        for (size_t j_1__ = 0; j_1__ < cor_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
    }


    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 b_bwt_j_1_max__ = Kc_bwt;
        for (size_t j_1__ = 0; j_1__ < b_bwt_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_bwt" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        param_name_stream__.str(std::string());
        param_name_stream__ << "Intercept_bwt";
        param_names__.push_back(param_name_stream__.str());
        size_t b_sigma_bwt_j_1_max__ = K_sigma_bwt;
        for (size_t j_1__ = 0; j_1__ < b_sigma_bwt_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_sigma_bwt" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t b_tarsus_j_1_max__ = Kc_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_tarsus_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_tarsus" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        param_name_stream__.str(std::string());
        param_name_stream__ << "Intercept_tarsus";
        param_names__.push_back(param_name_stream__.str());
        size_t b_sigma_tarsus_j_1_max__ = K_sigma_tarsus;
        for (size_t j_1__ = 0; j_1__ < b_sigma_tarsus_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "b_sigma_tarsus" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Lrescor_j_1_max__ = ((nresp * (nresp - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < Lrescor_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "Lrescor" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t sd_1_j_1_max__ = M_1;
        for (size_t j_1__ = 0; j_1__ < sd_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_1_j_2_max__ = N_1;
        size_t z_1_j_1_max__ = M_1;
        for (size_t j_2__ = 0; j_2__ < z_1_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_1_j_1_max__ = ((M_1 * (M_1 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < L_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "L_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t sd_2_j_1_max__ = M_2;
        for (size_t j_1__ = 0; j_1__ < sd_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_2_j_2_max__ = N_2;
        size_t z_2_j_1_max__ = M_2;
        for (size_t j_2__ = 0; j_2__ < z_2_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_2_j_1_max__ = ((M_2 * (M_2 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < L_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "L_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t sd_3_j_1_max__ = M_3;
        for (size_t j_1__ = 0; j_1__ < sd_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "sd_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t z_3_j_2_max__ = N_3;
        size_t z_3_j_1_max__ = M_3;
        for (size_t j_2__ = 0; j_2__ < z_3_j_2_max__; ++j_2__) {
            for (size_t j_1__ = 0; j_1__ < z_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "z_3" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }
        size_t L_3_j_1_max__ = ((M_3 * (M_3 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < L_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "L_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }

        if (!include_gqs__ && !include_tparams__) return;

        if (include_tparams__) {
            size_t r_1_j_2_max__ = M_1;
            size_t r_1_j_1_max__ = N_1;
            for (size_t j_2__ = 0; j_2__ < r_1_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_1_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_1" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_1_bwt_1_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_bwt_2_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_tarsus_3_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_1_tarsus_4_j_1_max__ = N_1;
            for (size_t j_1__ = 0; j_1__ < r_1_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_1_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_j_2_max__ = M_2;
            size_t r_2_j_1_max__ = N_2;
            for (size_t j_2__ = 0; j_2__ < r_2_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_2_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_2" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_2_bwt_1_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_bwt_2_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_tarsus_3_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_2_tarsus_4_j_1_max__ = N_2;
            for (size_t j_1__ = 0; j_1__ < r_2_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_2_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_j_2_max__ = M_3;
            size_t r_3_j_1_max__ = N_3;
            for (size_t j_2__ = 0; j_2__ < r_3_j_2_max__; ++j_2__) {
                for (size_t j_1__ = 0; j_1__ < r_3_j_1_max__; ++j_1__) {
                    param_name_stream__.str(std::string());
                    param_name_stream__ << "r_3" << '.' << j_1__ + 1 << '.' << j_2__ + 1;
                    param_names__.push_back(param_name_stream__.str());
                }
            }
            size_t r_3_bwt_1_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_1_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_bwt_1" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_bwt_2_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_bwt_2_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_bwt_2" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_tarsus_3_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_3_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_tarsus_3" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
            size_t r_3_tarsus_4_j_1_max__ = N_3;
            for (size_t j_1__ = 0; j_1__ < r_3_tarsus_4_j_1_max__; ++j_1__) {
                param_name_stream__.str(std::string());
                param_name_stream__ << "r_3_tarsus_4" << '.' << j_1__ + 1;
                param_names__.push_back(param_name_stream__.str());
            }
        }

        if (!include_gqs__) return;
        param_name_stream__.str(std::string());
        param_name_stream__ << "b_bwt_Intercept";
        param_names__.push_back(param_name_stream__.str());
        param_name_stream__.str(std::string());
        param_name_stream__ << "b_tarsus_Intercept";
        param_names__.push_back(param_name_stream__.str());
        size_t Rescor_j_1_max__ = ((nresp * (nresp - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < Rescor_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "Rescor" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t rescor_j_1_max__ = nrescor;
        for (size_t j_1__ = 0; j_1__ < rescor_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "rescor" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_1_j_1_max__ = ((M_1 * (M_1 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < Cor_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "Cor_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t cor_1_j_1_max__ = NC_1;
        for (size_t j_1__ = 0; j_1__ < cor_1_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_1" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_2_j_1_max__ = ((M_2 * (M_2 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < Cor_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "Cor_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t cor_2_j_1_max__ = NC_2;
        for (size_t j_1__ = 0; j_1__ < cor_2_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_2" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t Cor_3_j_1_max__ = ((M_3 * (M_3 - 1)) / 2);
        for (size_t j_1__ = 0; j_1__ < Cor_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "Cor_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
        size_t cor_3_j_1_max__ = NC_3;
        for (size_t j_1__ = 0; j_1__ < cor_3_j_1_max__; ++j_1__) {
            param_name_stream__.str(std::string());
            param_name_stream__ << "cor_3" << '.' << j_1__ + 1;
            param_names__.push_back(param_name_stream__.str());
        }
    }

}; // model

}  // namespace

typedef model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb_namespace::model1a4c67bd7c3f_9e1f937d1252c795015f7237580a6dfb 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