text: /*
      bayesreg is free software: you can redistribute it and/or modify
      it under the terms of the GNU General Public License as published by
      the Free Software Foundation, either version 3 of the License, or
      (at your option) any later version.
  
      bayesreg is distributed in the hope that it will be useful,
      but WITHOUT ANY WARRANTY; without even the implied warranty of
      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
      GNU General Public License for more details.
  
      You should have received a copy of the GNU General Public License
      along with bayesreg.  If not, see <http://www.gnu.org/licenses/>.
  */
  #ifndef MODELS_HPP
  #define MODELS_HPP
  #define STAN__SERVICES__COMMAND_HPP
  #include <rstan/rstaninc.hpp>
  // Code generated by Stan version 2.18.0
  
  #include <stan/model/model_header.hpp>
  
  namespace model_lm_noNA_regularizedHorseshoe_FB_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", "model_lm_noNA_regularizedHorseshoe_FB");
      reader.add_event(66, 64, "end", "model_lm_noNA_regularizedHorseshoe_FB");
      return reader;
  }
  
  #include <meta_header.hpp>
   class model_lm_noNA_regularizedHorseshoe_FB : public prob_grad {
  private:
      int N_train;
      int p;
      vector<double> y_train;
      matrix_d X_train;
      int N_test;
      matrix_d X_test;
      double scale_global;
      double nu_global;
      double nu_local;
      double slab_scale;
      double slab_df;
  public:
      model_lm_noNA_regularizedHorseshoe_FB(stan::io::var_context& context__,
          std::ostream* pstream__ = 0)
          : prob_grad(0) {
          ctor_body(context__, 0, pstream__);
      }
  
      model_lm_noNA_regularizedHorseshoe_FB(stan::io::var_context& context__,
          unsigned int random_seed__,
          std::ostream* pstream__ = 0)
          : prob_grad(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__ = "model_lm_noNA_regularizedHorseshoe_FB_namespace::model_lm_noNA_regularizedHorseshoe_FB";
          (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
  
          // initialize member variables
          try {
              current_statement_begin__ = 3;
              context__.validate_dims("data initialization", "N_train", "int", context__.to_vec());
              N_train = int(0);
              vals_i__ = context__.vals_i("N_train");
              pos__ = 0;
              N_train = vals_i__[pos__++];
              current_statement_begin__ = 4;
              context__.validate_dims("data initialization", "p", "int", context__.to_vec());
              p = int(0);
              vals_i__ = context__.vals_i("p");
              pos__ = 0;
              p = vals_i__[pos__++];
              current_statement_begin__ = 5;
              validate_non_negative_index("y_train", "N_train", N_train);
              context__.validate_dims("data initialization", "y_train", "double", context__.to_vec(N_train));
              validate_non_negative_index("y_train", "N_train", N_train);
              y_train = std::vector<double>(N_train,double(0));
              vals_r__ = context__.vals_r("y_train");
              pos__ = 0;
              size_t y_train_limit_0__ = N_train;
              for (size_t i_0__ = 0; i_0__ < y_train_limit_0__; ++i_0__) {
                  y_train[i_0__] = vals_r__[pos__++];
              }
              current_statement_begin__ = 6;
              validate_non_negative_index("X_train", "N_train", N_train);
              validate_non_negative_index("X_train", "p", p);
              context__.validate_dims("data initialization", "X_train", "matrix_d", context__.to_vec(N_train,p));
              validate_non_negative_index("X_train", "N_train", N_train);
              validate_non_negative_index("X_train", "p", p);
              X_train = matrix_d(static_cast<Eigen::VectorXd::Index>(N_train),static_cast<Eigen::VectorXd::Index>(p));
              vals_r__ = context__.vals_r("X_train");
              pos__ = 0;
              size_t X_train_m_mat_lim__ = N_train;
              size_t X_train_n_mat_lim__ = p;
              for (size_t n_mat__ = 0; n_mat__ < X_train_n_mat_lim__; ++n_mat__) {
                  for (size_t m_mat__ = 0; m_mat__ < X_train_m_mat_lim__; ++m_mat__) {
                      X_train(m_mat__,n_mat__) = vals_r__[pos__++];
                  }
              }
              current_statement_begin__ = 8;
              context__.validate_dims("data initialization", "N_test", "int", context__.to_vec());
              N_test = int(0);
              vals_i__ = context__.vals_i("N_test");
              pos__ = 0;
              N_test = vals_i__[pos__++];
              current_statement_begin__ = 9;
              validate_non_negative_index("X_test", "N_test", N_test);
              validate_non_negative_index("X_test", "p", p);
              context__.validate_dims("data initialization", "X_test", "matrix_d", context__.to_vec(N_test,p));
              validate_non_negative_index("X_test", "N_test", N_test);
              validate_non_negative_index("X_test", "p", p);
              X_test = matrix_d(static_cast<Eigen::VectorXd::Index>(N_test),static_cast<Eigen::VectorXd::Index>(p));
              vals_r__ = context__.vals_r("X_test");
              pos__ = 0;
              size_t X_test_m_mat_lim__ = N_test;
              size_t X_test_n_mat_lim__ = p;
              for (size_t n_mat__ = 0; n_mat__ < X_test_n_mat_lim__; ++n_mat__) {
                  for (size_t m_mat__ = 0; m_mat__ < X_test_m_mat_lim__; ++m_mat__) {
                      X_test(m_mat__,n_mat__) = vals_r__[pos__++];
                  }
              }
              current_statement_begin__ = 10;
              context__.validate_dims("data initialization", "scale_global", "double", context__.to_vec());
              scale_global = double(0);
              vals_r__ = context__.vals_r("scale_global");
              pos__ = 0;
              scale_global = vals_r__[pos__++];
              current_statement_begin__ = 11;
              context__.validate_dims("data initialization", "nu_global", "double", context__.to_vec());
              nu_global = double(0);
              vals_r__ = context__.vals_r("nu_global");
              pos__ = 0;
              nu_global = vals_r__[pos__++];
              current_statement_begin__ = 12;
              context__.validate_dims("data initialization", "nu_local", "double", context__.to_vec());
              nu_local = double(0);
              vals_r__ = context__.vals_r("nu_local");
              pos__ = 0;
              nu_local = vals_r__[pos__++];
              current_statement_begin__ = 13;
              context__.validate_dims("data initialization", "slab_scale", "double", context__.to_vec());
              slab_scale = double(0);
              vals_r__ = context__.vals_r("slab_scale");
              pos__ = 0;
              slab_scale = vals_r__[pos__++];
              current_statement_begin__ = 14;
              context__.validate_dims("data initialization", "slab_df", "double", context__.to_vec());
              slab_df = double(0);
              vals_r__ = context__.vals_r("slab_df");
              pos__ = 0;
              slab_df = vals_r__[pos__++];
  
              // validate, data variables
              current_statement_begin__ = 3;
              current_statement_begin__ = 4;
              current_statement_begin__ = 5;
              current_statement_begin__ = 6;
              current_statement_begin__ = 8;
              current_statement_begin__ = 9;
              current_statement_begin__ = 10;
              check_greater_or_equal(function__,"scale_global",scale_global,0);
              current_statement_begin__ = 11;
              check_greater_or_equal(function__,"nu_global",nu_global,1);
              current_statement_begin__ = 12;
              check_greater_or_equal(function__,"nu_local",nu_local,1);
              current_statement_begin__ = 13;
              check_greater_or_equal(function__,"slab_scale",slab_scale,0);
              current_statement_begin__ = 14;
              check_greater_or_equal(function__,"slab_df",slab_df,0);
              // initialize data variables
  
  
              // validate transformed data
  
              // validate, set parameter ranges
              num_params_r__ = 0U;
              param_ranges_i__.clear();
              current_statement_begin__ = 17;
              ++num_params_r__;
              current_statement_begin__ = 18;
              ++num_params_r__;
              current_statement_begin__ = 20;
              validate_non_negative_index("z", "p", p);
              num_params_r__ += p;
              current_statement_begin__ = 21;
              ++num_params_r__;
              current_statement_begin__ = 22;
              ++num_params_r__;
              current_statement_begin__ = 23;
              validate_non_negative_index("aux1_local", "p", p);
              num_params_r__ += p;
              current_statement_begin__ = 24;
              validate_non_negative_index("aux2_local", "p", p);
              num_params_r__ += p;
              current_statement_begin__ = 25;
              ++num_params_r__;
          } 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 ***");
          }
      }
  
      ~model_lm_noNA_regularizedHorseshoe_FB() { }
  
  
      void transform_inits(const stan::io::var_context& context__,
                           std::vector<int>& params_i__,
                           std::vector<double>& params_r__,
                           std::ostream* pstream__) const {
          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__;
  
          if (!(context__.contains_r("mu")))
              throw std::runtime_error("variable mu missing");
          vals_r__ = context__.vals_r("mu");
          pos__ = 0U;
          context__.validate_dims("initialization", "mu", "double", context__.to_vec());
          double mu(0);
          mu = vals_r__[pos__++];
          try {
              writer__.scalar_unconstrain(mu);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable mu: ") + e.what());
          }
  
          if (!(context__.contains_r("sigma2")))
              throw std::runtime_error("variable sigma2 missing");
          vals_r__ = context__.vals_r("sigma2");
          pos__ = 0U;
          context__.validate_dims("initialization", "sigma2", "double", context__.to_vec());
          double sigma2(0);
          sigma2 = vals_r__[pos__++];
          try {
              writer__.scalar_lb_unconstrain(0,sigma2);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable sigma2: ") + e.what());
          }
  
          if (!(context__.contains_r("z")))
              throw std::runtime_error("variable z missing");
          vals_r__ = context__.vals_r("z");
          pos__ = 0U;
          validate_non_negative_index("z", "p", p);
          context__.validate_dims("initialization", "z", "vector_d", context__.to_vec(p));
          vector_d z(static_cast<Eigen::VectorXd::Index>(p));
          for (int j1__ = 0U; j1__ < p; ++j1__)
              z(j1__) = vals_r__[pos__++];
          try {
              writer__.vector_unconstrain(z);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable z: ") + e.what());
          }
  
          if (!(context__.contains_r("aux1_global")))
              throw std::runtime_error("variable aux1_global missing");
          vals_r__ = context__.vals_r("aux1_global");
          pos__ = 0U;
          context__.validate_dims("initialization", "aux1_global", "double", context__.to_vec());
          double aux1_global(0);
          aux1_global = vals_r__[pos__++];
          try {
              writer__.scalar_lb_unconstrain(0,aux1_global);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable aux1_global: ") + e.what());
          }
  
          if (!(context__.contains_r("aux2_global")))
              throw std::runtime_error("variable aux2_global missing");
          vals_r__ = context__.vals_r("aux2_global");
          pos__ = 0U;
          context__.validate_dims("initialization", "aux2_global", "double", context__.to_vec());
          double aux2_global(0);
          aux2_global = vals_r__[pos__++];
          try {
              writer__.scalar_lb_unconstrain(0,aux2_global);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable aux2_global: ") + e.what());
          }
  
          if (!(context__.contains_r("aux1_local")))
              throw std::runtime_error("variable aux1_local missing");
          vals_r__ = context__.vals_r("aux1_local");
          pos__ = 0U;
          validate_non_negative_index("aux1_local", "p", p);
          context__.validate_dims("initialization", "aux1_local", "vector_d", context__.to_vec(p));
          vector_d aux1_local(static_cast<Eigen::VectorXd::Index>(p));
          for (int j1__ = 0U; j1__ < p; ++j1__)
              aux1_local(j1__) = vals_r__[pos__++];
          try {
              writer__.vector_lb_unconstrain(0,aux1_local);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable aux1_local: ") + e.what());
          }
  
          if (!(context__.contains_r("aux2_local")))
              throw std::runtime_error("variable aux2_local missing");
          vals_r__ = context__.vals_r("aux2_local");
          pos__ = 0U;
          validate_non_negative_index("aux2_local", "p", p);
          context__.validate_dims("initialization", "aux2_local", "vector_d", context__.to_vec(p));
          vector_d aux2_local(static_cast<Eigen::VectorXd::Index>(p));
          for (int j1__ = 0U; j1__ < p; ++j1__)
              aux2_local(j1__) = vals_r__[pos__++];
          try {
              writer__.vector_lb_unconstrain(0,aux2_local);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable aux2_local: ") + e.what());
          }
  
          if (!(context__.contains_r("caux")))
              throw std::runtime_error("variable caux missing");
          vals_r__ = context__.vals_r("caux");
          pos__ = 0U;
          context__.validate_dims("initialization", "caux", "double", context__.to_vec());
          double caux(0);
          caux = vals_r__[pos__++];
          try {
              writer__.scalar_lb_unconstrain(0,caux);
          } catch (const std::exception& e) { 
              throw std::runtime_error(std::string("Error transforming variable caux: ") + e.what());
          }
  
          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(vector<T__>& params_r__,
                   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__;  // suppress unused var warning
  
          T__ lp__(0.0);
          stan::math::accumulator<T__> lp_accum__;
  
          try {
              // model parameters
              stan::io::reader<local_scalar_t__> in__(params_r__,params_i__);
  
              local_scalar_t__ mu;
              (void) mu;  // dummy to suppress unused var warning
              if (jacobian__)
                  mu = in__.scalar_constrain(lp__);
              else
                  mu = in__.scalar_constrain();
  
              local_scalar_t__ sigma2;
              (void) sigma2;  // dummy to suppress unused var warning
              if (jacobian__)
                  sigma2 = in__.scalar_lb_constrain(0,lp__);
              else
                  sigma2 = in__.scalar_lb_constrain(0);
  
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  z;
              (void) z;  // dummy to suppress unused var warning
              if (jacobian__)
                  z = in__.vector_constrain(p,lp__);
              else
                  z = in__.vector_constrain(p);
  
              local_scalar_t__ aux1_global;
              (void) aux1_global;  // dummy to suppress unused var warning
              if (jacobian__)
                  aux1_global = in__.scalar_lb_constrain(0,lp__);
              else
                  aux1_global = in__.scalar_lb_constrain(0);
  
              local_scalar_t__ aux2_global;
              (void) aux2_global;  // dummy to suppress unused var warning
              if (jacobian__)
                  aux2_global = in__.scalar_lb_constrain(0,lp__);
              else
                  aux2_global = in__.scalar_lb_constrain(0);
  
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  aux1_local;
              (void) aux1_local;  // dummy to suppress unused var warning
              if (jacobian__)
                  aux1_local = in__.vector_lb_constrain(0,p,lp__);
              else
                  aux1_local = in__.vector_lb_constrain(0,p);
  
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  aux2_local;
              (void) aux2_local;  // dummy to suppress unused var warning
              if (jacobian__)
                  aux2_local = in__.vector_lb_constrain(0,p,lp__);
              else
                  aux2_local = in__.vector_lb_constrain(0,p);
  
              local_scalar_t__ caux;
              (void) caux;  // dummy to suppress unused var warning
              if (jacobian__)
                  caux = in__.scalar_lb_constrain(0,lp__);
              else
                  caux = in__.scalar_lb_constrain(0);
  
  
              // transformed parameters
              current_statement_begin__ = 28;
              local_scalar_t__ sigma;
              (void) sigma;  // dummy to suppress unused var warning
  
              stan::math::initialize(sigma, DUMMY_VAR__);
              stan::math::fill(sigma,DUMMY_VAR__);
              current_statement_begin__ = 29;
              validate_non_negative_index("linpred", "N_train", N_train);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  linpred(static_cast<Eigen::VectorXd::Index>(N_train));
              (void) linpred;  // dummy to suppress unused var warning
  
              stan::math::initialize(linpred, DUMMY_VAR__);
              stan::math::fill(linpred,DUMMY_VAR__);
              current_statement_begin__ = 30;
              local_scalar_t__ tau;
              (void) tau;  // dummy to suppress unused var warning
  
              stan::math::initialize(tau, DUMMY_VAR__);
              stan::math::fill(tau,DUMMY_VAR__);
              current_statement_begin__ = 31;
              validate_non_negative_index("lambda", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  lambda(static_cast<Eigen::VectorXd::Index>(p));
              (void) lambda;  // dummy to suppress unused var warning
  
              stan::math::initialize(lambda, DUMMY_VAR__);
              stan::math::fill(lambda,DUMMY_VAR__);
              current_statement_begin__ = 32;
              validate_non_negative_index("lambda_tilde", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  lambda_tilde(static_cast<Eigen::VectorXd::Index>(p));
              (void) lambda_tilde;  // dummy to suppress unused var warning
  
              stan::math::initialize(lambda_tilde, DUMMY_VAR__);
              stan::math::fill(lambda_tilde,DUMMY_VAR__);
              current_statement_begin__ = 33;
              local_scalar_t__ c;
              (void) c;  // dummy to suppress unused var warning
  
              stan::math::initialize(c, DUMMY_VAR__);
              stan::math::fill(c,DUMMY_VAR__);
              current_statement_begin__ = 34;
              validate_non_negative_index("beta", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  beta(static_cast<Eigen::VectorXd::Index>(p));
              (void) beta;  // dummy to suppress unused var warning
  
              stan::math::initialize(beta, DUMMY_VAR__);
              stan::math::fill(beta,DUMMY_VAR__);
  
  
              current_statement_begin__ = 35;
              stan::math::assign(sigma, stan::math::sqrt(sigma2));
              current_statement_begin__ = 36;
              stan::math::assign(lambda, elt_multiply(aux1_local,stan::math::sqrt(aux2_local)));
              current_statement_begin__ = 37;
              stan::math::assign(tau, (((aux1_global * stan::math::sqrt(aux2_global)) * scale_global) * sigma));
              current_statement_begin__ = 38;
              stan::math::assign(c, (slab_scale * stan::math::sqrt(caux)));
              current_statement_begin__ = 39;
              stan::math::assign(lambda_tilde, stan::math::sqrt(elt_divide(multiply(pow(c,2),square(lambda)),add(pow(c,2),multiply(pow(tau,2),square(lambda))))));
              current_statement_begin__ = 40;
              stan::math::assign(beta, multiply(elt_multiply(z,lambda_tilde),tau));
              current_statement_begin__ = 41;
              stan::math::assign(linpred, add(mu,multiply(X_train,beta)));
  
              // validate transformed parameters
              if (stan::math::is_uninitialized(sigma)) {
                  std::stringstream msg__;
                  msg__ << "Undefined transformed parameter: sigma";
                  throw std::runtime_error(msg__.str());
              }
              for (int i0__ = 0; i0__ < N_train; ++i0__) {
                  if (stan::math::is_uninitialized(linpred(i0__))) {
                      std::stringstream msg__;
                      msg__ << "Undefined transformed parameter: linpred" << '[' << i0__ << ']';
                      throw std::runtime_error(msg__.str());
                  }
              }
              if (stan::math::is_uninitialized(tau)) {
                  std::stringstream msg__;
                  msg__ << "Undefined transformed parameter: tau";
                  throw std::runtime_error(msg__.str());
              }
              for (int i0__ = 0; i0__ < p; ++i0__) {
                  if (stan::math::is_uninitialized(lambda(i0__))) {
                      std::stringstream msg__;
                      msg__ << "Undefined transformed parameter: lambda" << '[' << i0__ << ']';
                      throw std::runtime_error(msg__.str());
                  }
              }
              for (int i0__ = 0; i0__ < p; ++i0__) {
                  if (stan::math::is_uninitialized(lambda_tilde(i0__))) {
                      std::stringstream msg__;
                      msg__ << "Undefined transformed parameter: lambda_tilde" << '[' << i0__ << ']';
                      throw std::runtime_error(msg__.str());
                  }
              }
              if (stan::math::is_uninitialized(c)) {
                  std::stringstream msg__;
                  msg__ << "Undefined transformed parameter: c";
                  throw std::runtime_error(msg__.str());
              }
              for (int i0__ = 0; i0__ < p; ++i0__) {
                  if (stan::math::is_uninitialized(beta(i0__))) {
                      std::stringstream msg__;
                      msg__ << "Undefined transformed parameter: beta" << '[' << i0__ << ']';
                      throw std::runtime_error(msg__.str());
                  }
              }
  
              const char* function__ = "validate transformed params";
              (void) function__;  // dummy to suppress unused var warning
              current_statement_begin__ = 28;
              check_greater_or_equal(function__,"sigma",sigma,0);
              current_statement_begin__ = 29;
              current_statement_begin__ = 30;
              check_greater_or_equal(function__,"tau",tau,0);
              current_statement_begin__ = 31;
              check_greater_or_equal(function__,"lambda",lambda,0);
              current_statement_begin__ = 32;
              check_greater_or_equal(function__,"lambda_tilde",lambda_tilde,0);
              current_statement_begin__ = 33;
              check_greater_or_equal(function__,"c",c,0);
              current_statement_begin__ = 34;
  
              // model body
  
              current_statement_begin__ = 46;
              lp_accum__.add(normal_log<propto__>(z, 0, 1));
              current_statement_begin__ = 47;
              lp_accum__.add(normal_log<propto__>(aux1_local, 0, 1));
              current_statement_begin__ = 48;
              lp_accum__.add(inv_gamma_log<propto__>(aux2_local, (0.5 * nu_local), (0.5 * nu_local)));
              current_statement_begin__ = 49;
              lp_accum__.add(normal_log<propto__>(aux1_global, 0, 1));
              current_statement_begin__ = 50;
              lp_accum__.add(inv_gamma_log<propto__>(aux2_global, (0.5 * nu_global), (0.5 * nu_global)));
              current_statement_begin__ = 51;
              lp_accum__.add(inv_gamma_log<propto__>(caux, (0.5 * slab_df), (0.5 * slab_df)));
              current_statement_begin__ = 54;
              lp_accum__.add((-(2) * stan::math::log(sigma)));
              current_statement_begin__ = 57;
              lp_accum__.add(normal_log<propto__>(y_train, linpred, sigma));
  
          } 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("mu");
          names__.push_back("sigma2");
          names__.push_back("z");
          names__.push_back("aux1_global");
          names__.push_back("aux2_global");
          names__.push_back("aux1_local");
          names__.push_back("aux2_local");
          names__.push_back("caux");
          names__.push_back("sigma");
          names__.push_back("linpred");
          names__.push_back("tau");
          names__.push_back("lambda");
          names__.push_back("lambda_tilde");
          names__.push_back("c");
          names__.push_back("beta");
          names__.push_back("y_test");
      }
  
  
      void get_dims(std::vector<std::vector<size_t> >& dimss__) const {
          dimss__.resize(0);
          std::vector<size_t> dims__;
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(N_train);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(p);
          dimss__.push_back(dims__);
          dims__.resize(0);
          dims__.push_back(N_test);
          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__ = "model_lm_noNA_regularizedHorseshoe_FB_namespace::write_array";
          (void) function__;  // dummy to suppress unused var warning
          // read-transform, write parameters
          double mu = in__.scalar_constrain();
          double sigma2 = in__.scalar_lb_constrain(0);
          vector_d z = in__.vector_constrain(p);
          double aux1_global = in__.scalar_lb_constrain(0);
          double aux2_global = in__.scalar_lb_constrain(0);
          vector_d aux1_local = in__.vector_lb_constrain(0,p);
          vector_d aux2_local = in__.vector_lb_constrain(0,p);
          double caux = in__.scalar_lb_constrain(0);
          vars__.push_back(mu);
          vars__.push_back(sigma2);
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(z[k_0__]);
              }
          vars__.push_back(aux1_global);
          vars__.push_back(aux2_global);
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(aux1_local[k_0__]);
              }
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(aux2_local[k_0__]);
              }
          vars__.push_back(caux);
  
          // declare and define transformed parameters
          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
  
          try {
              current_statement_begin__ = 28;
              local_scalar_t__ sigma;
              (void) sigma;  // dummy to suppress unused var warning
  
              stan::math::initialize(sigma, DUMMY_VAR__);
              stan::math::fill(sigma,DUMMY_VAR__);
              current_statement_begin__ = 29;
              validate_non_negative_index("linpred", "N_train", N_train);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  linpred(static_cast<Eigen::VectorXd::Index>(N_train));
              (void) linpred;  // dummy to suppress unused var warning
  
              stan::math::initialize(linpred, DUMMY_VAR__);
              stan::math::fill(linpred,DUMMY_VAR__);
              current_statement_begin__ = 30;
              local_scalar_t__ tau;
              (void) tau;  // dummy to suppress unused var warning
  
              stan::math::initialize(tau, DUMMY_VAR__);
              stan::math::fill(tau,DUMMY_VAR__);
              current_statement_begin__ = 31;
              validate_non_negative_index("lambda", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  lambda(static_cast<Eigen::VectorXd::Index>(p));
              (void) lambda;  // dummy to suppress unused var warning
  
              stan::math::initialize(lambda, DUMMY_VAR__);
              stan::math::fill(lambda,DUMMY_VAR__);
              current_statement_begin__ = 32;
              validate_non_negative_index("lambda_tilde", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  lambda_tilde(static_cast<Eigen::VectorXd::Index>(p));
              (void) lambda_tilde;  // dummy to suppress unused var warning
  
              stan::math::initialize(lambda_tilde, DUMMY_VAR__);
              stan::math::fill(lambda_tilde,DUMMY_VAR__);
              current_statement_begin__ = 33;
              local_scalar_t__ c;
              (void) c;  // dummy to suppress unused var warning
  
              stan::math::initialize(c, DUMMY_VAR__);
              stan::math::fill(c,DUMMY_VAR__);
              current_statement_begin__ = 34;
              validate_non_negative_index("beta", "p", p);
              Eigen::Matrix<local_scalar_t__,Eigen::Dynamic,1>  beta(static_cast<Eigen::VectorXd::Index>(p));
              (void) beta;  // dummy to suppress unused var warning
  
              stan::math::initialize(beta, DUMMY_VAR__);
              stan::math::fill(beta,DUMMY_VAR__);
  
  
              current_statement_begin__ = 35;
              stan::math::assign(sigma, stan::math::sqrt(sigma2));
              current_statement_begin__ = 36;
              stan::math::assign(lambda, elt_multiply(aux1_local,stan::math::sqrt(aux2_local)));
              current_statement_begin__ = 37;
              stan::math::assign(tau, (((aux1_global * stan::math::sqrt(aux2_global)) * scale_global) * sigma));
              current_statement_begin__ = 38;
              stan::math::assign(c, (slab_scale * stan::math::sqrt(caux)));
              current_statement_begin__ = 39;
              stan::math::assign(lambda_tilde, stan::math::sqrt(elt_divide(multiply(pow(c,2),square(lambda)),add(pow(c,2),multiply(pow(tau,2),square(lambda))))));
              current_statement_begin__ = 40;
              stan::math::assign(beta, multiply(elt_multiply(z,lambda_tilde),tau));
              current_statement_begin__ = 41;
              stan::math::assign(linpred, add(mu,multiply(X_train,beta)));
  
              // validate transformed parameters
              current_statement_begin__ = 28;
              check_greater_or_equal(function__,"sigma",sigma,0);
              current_statement_begin__ = 29;
              current_statement_begin__ = 30;
              check_greater_or_equal(function__,"tau",tau,0);
              current_statement_begin__ = 31;
              check_greater_or_equal(function__,"lambda",lambda,0);
              current_statement_begin__ = 32;
              check_greater_or_equal(function__,"lambda_tilde",lambda_tilde,0);
              current_statement_begin__ = 33;
              check_greater_or_equal(function__,"c",c,0);
              current_statement_begin__ = 34;
  
              // write transformed parameters
              if (include_tparams__) {
          vars__.push_back(sigma);
              for (int k_0__ = 0; k_0__ < N_train; ++k_0__) {
              vars__.push_back(linpred[k_0__]);
              }
          vars__.push_back(tau);
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(lambda[k_0__]);
              }
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(lambda_tilde[k_0__]);
              }
          vars__.push_back(c);
              for (int k_0__ = 0; k_0__ < p; ++k_0__) {
              vars__.push_back(beta[k_0__]);
              }
              }
              if (!include_gqs__) return;
              // declare and define generated quantities
              current_statement_begin__ = 60;
              validate_non_negative_index("y_test", "N_test", N_test);
              vector<local_scalar_t__> y_test(N_test);
              stan::math::initialize(y_test, DUMMY_VAR__);
              stan::math::fill(y_test,DUMMY_VAR__);
  
  
              current_statement_begin__ = 61;
              for (int i = 1; i <= N_test; ++i) {
  
                  current_statement_begin__ = 62;
                  stan::model::assign(y_test, 
                              stan::model::cons_list(stan::model::index_uni(i), stan::model::nil_index_list()), 
                              normal_rng((mu + multiply(stan::model::rvalue(X_test, stan::model::cons_list(stan::model::index_uni(i), stan::model::cons_list(stan::model::index_omni(), stan::model::nil_index_list())), "X_test"),beta)),sigma, base_rng__), 
                              "assigning variable y_test");
              }
  
              // validate generated quantities
              current_statement_begin__ = 60;
  
              // write generated quantities
              for (int k_0__ = 0; k_0__ < N_test; ++k_0__) {
              vars__.push_back(y_test[k_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 ***");
          }
      }
  
      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];
      }
  
      static std::string model_name() {
          return "model_lm_noNA_regularizedHorseshoe_FB";
      }
  
  
      void constrained_param_names(std::vector<std::string>& param_names__,
                                   bool include_tparams__ = true,
                                   bool include_gqs__ = true) const {
          std::stringstream param_name_stream__;
          param_name_stream__.str(std::string());
          param_name_stream__ << "mu";
          param_names__.push_back(param_name_stream__.str());
          param_name_stream__.str(std::string());
          param_name_stream__ << "sigma2";
          param_names__.push_back(param_name_stream__.str());
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "z" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          param_name_stream__.str(std::string());
          param_name_stream__ << "aux1_global";
          param_names__.push_back(param_name_stream__.str());
          param_name_stream__.str(std::string());
          param_name_stream__ << "aux2_global";
          param_names__.push_back(param_name_stream__.str());
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "aux1_local" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "aux2_local" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          param_name_stream__.str(std::string());
          param_name_stream__ << "caux";
          param_names__.push_back(param_name_stream__.str());
  
          if (!include_gqs__ && !include_tparams__) return;
  
          if (include_tparams__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "sigma";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= N_train; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "linpred" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              param_name_stream__.str(std::string());
              param_name_stream__ << "tau";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "lambda" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "lambda_tilde" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              param_name_stream__.str(std::string());
              param_name_stream__ << "c";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "beta" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
          }
  
  
          if (!include_gqs__) return;
          for (int k_0__ = 1; k_0__ <= N_test; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "y_test" << '.' << k_0__;
              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__;
          param_name_stream__.str(std::string());
          param_name_stream__ << "mu";
          param_names__.push_back(param_name_stream__.str());
          param_name_stream__.str(std::string());
          param_name_stream__ << "sigma2";
          param_names__.push_back(param_name_stream__.str());
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "z" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          param_name_stream__.str(std::string());
          param_name_stream__ << "aux1_global";
          param_names__.push_back(param_name_stream__.str());
          param_name_stream__.str(std::string());
          param_name_stream__ << "aux2_global";
          param_names__.push_back(param_name_stream__.str());
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "aux1_local" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "aux2_local" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
          param_name_stream__.str(std::string());
          param_name_stream__ << "caux";
          param_names__.push_back(param_name_stream__.str());
  
          if (!include_gqs__ && !include_tparams__) return;
  
          if (include_tparams__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "sigma";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= N_train; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "linpred" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              param_name_stream__.str(std::string());
              param_name_stream__ << "tau";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "lambda" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "lambda_tilde" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
              param_name_stream__.str(std::string());
              param_name_stream__ << "c";
              param_names__.push_back(param_name_stream__.str());
              for (int k_0__ = 1; k_0__ <= p; ++k_0__) {
                  param_name_stream__.str(std::string());
                  param_name_stream__ << "beta" << '.' << k_0__;
                  param_names__.push_back(param_name_stream__.str());
              }
          }
  
  
          if (!include_gqs__) return;
          for (int k_0__ = 1; k_0__ <= N_test; ++k_0__) {
              param_name_stream__.str(std::string());
              param_name_stream__ << "y_test" << '.' << k_0__;
              param_names__.push_back(param_name_stream__.str());
          }
      }
  
  }; // model
  
  }
  
  typedef model_lm_noNA_regularizedHorseshoe_FB_namespace::model_lm_noNA_regularizedHorseshoe_FB stan_model;
  
  
  #endif
  