From a3b2a5a59d3319454ac69bfd6c38026ff4a9b8b0 Mon Sep 17 00:00:00 2001 From: Fonti Kar Date: Mon, 9 Oct 2023 15:07:35 +1100 Subject: [PATCH] Implemented two models, and workflow to how to add data to template to running model --- DESCRIPTION | 5 +- NAMESPACE | 4 +- R/rmot.R | 16 - R/rmot_assign_data.R | 31 + R/rmot_models.R | 80 +- R/rmot_run.R | 20 + R/stanmodels.R | 5 +- R/utils.R | 11 + inst/stan/constant_single.stan | 83 ++ inst/stan/{lm.stan => linear.stan} | 6 +- man/rmot_assign_data.Rd | 22 + man/rmot_model.Rd | 20 + man/{rmot_lm.Rd => rmot_run.Rd} | 14 +- src/RcppExports.cpp | 6 +- src/stanExports_constant_single.cc | 32 + src/stanExports_constant_single.h | 812 ++++++++++++++++++ ...tanExports_lm.cc => stanExports_linear.cc} | 6 +- ...{stanExports_lm.h => stanExports_linear.h} | 83 +- 18 files changed, 1134 insertions(+), 122 deletions(-) delete mode 100644 R/rmot.R create mode 100644 R/rmot_assign_data.R create mode 100644 R/rmot_run.R create mode 100644 R/utils.R create mode 100644 inst/stan/constant_single.stan rename inst/stan/{lm.stan => linear.stan} (66%) create mode 100644 man/rmot_assign_data.Rd create mode 100644 man/rmot_model.Rd rename man/{rmot_lm.Rd => rmot_run.Rd} (54%) create mode 100644 src/stanExports_constant_single.cc create mode 100644 src/stanExports_constant_single.h rename src/{stanExports_lm.cc => stanExports_linear.cc} (95%) rename src/{stanExports_lm.h => stanExports_linear.h} (86%) diff --git a/DESCRIPTION b/DESCRIPTION index 33ea190..3d20acc 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -15,10 +15,13 @@ Depends: R (>= 3.4.0) Imports: methods, + purrr, Rcpp (>= 0.12.0), RcppParallel (>= 5.0.1), + rlang, rstan (>= 2.18.1), - rstantools (>= 2.3.1.1) + rstantools (>= 2.3.1.1), + stats LinkingTo: BH (>= 1.66.0), Rcpp (>= 0.12.0), diff --git a/NAMESPACE b/NAMESPACE index 61840e3..c3c42cf 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,6 +1,8 @@ # Generated by roxygen2: do not edit by hand -export(rmot_lm) +export(rmot_assign_data) +export(rmot_model) +export(rmot_run) import(Rcpp) import(methods) importFrom(RcppParallel,RcppParallelLibs) diff --git a/R/rmot.R b/R/rmot.R deleted file mode 100644 index 16892a1..0000000 --- a/R/rmot.R +++ /dev/null @@ -1,16 +0,0 @@ -#' Run a linear model in Stan -#' -#' @param x vector of numeric or integer predictor variable -#' @param y vector of numeric or integer response variable -#' @param ... additional arguments passed to rstan::sampling -#' -#' @return stanfit model output -#' @export -#' -#' @examples -#' mtcars -#' rmot_lm(mtcars$mpg, mtcars$disp) -rmot <- function(x, y, ...) { - out <- rstan::sampling(stanmodels$model, data = standata, ...) - return(out) -} diff --git a/R/rmot_assign_data.R b/R/rmot_assign_data.R new file mode 100644 index 0000000..80b7028 --- /dev/null +++ b/R/rmot_assign_data.R @@ -0,0 +1,31 @@ +#' Assign data to template +#' +#' @param model_template output from rmot_model +#' @param ... data-masking name-value pairs +#' +#' @return updated named list with your data assigned to Stan model parameters +#' @export +#' +#' @examples +#' rmot_assign_data(X = Loblolly$age, Y = Loblolly$height) +rmot_assign_data <- function(model_template, ...){ + # Grab user expressions + user_code <- rlang::enexprs(..., .check_assign = TRUE) + + # Grab the names + fields <- names(user_code) + + #TODO: Check names are in model_template + + # Evaluate the RHS of expressions (the values) + data <- purrr::map(user_code, + eval) + + for(i in fields){ + model_template <- purrr::list_modify(model_template, !!!data[i]) + } + + #TODO: Check if N is supplied, if not, assign by default to length(X) and give warning + + return(model_template) +} diff --git a/R/rmot_models.R b/R/rmot_models.R index aa94057..512fe5b 100644 --- a/R/rmot_models.R +++ b/R/rmot_models.R @@ -1,61 +1,51 @@ -# Set list structures for different models -# An example for lm - -rmot_lm <- function(){ -list(X = NULL, - Y = NULL, - N = NULL, - model = "linear") -} +#' Select data configuration template for rmot supported model +#' +#' @param model model name character string +#' +#' @return named list that matches Stan model parameters +#' @export +#' +#' @examples +#' rmot_model("linear") +rmot_model <- function(model=NULL){ -# Need a mechanism to select models -# rmot_config(model = "linear") + #TODO: Need a mechanism to check model requested in one that is supported by rmot -rmot_config <- function(model=NULL){ output <- switch(model, - linear = rmot_lm()) + linear = rmot_lm(), + constant_single = rmot_cgs()) class(output) <- "rmot_object" return(output) } -# Need a mechanism to take user data and assign to slots in list -rmot_assign_data <- function(model_template, field, data){ - purrr::assign_in(model_template, field, data) -} - - -rmot_assign_data <- function(model_template, ...){ - # Grab user expressions - user_code <- rlang::enexprs(..., .check_assign = TRUE) +#' Data configuration template for linear model +#' @keywords internal +#' @noRd - # Evaluate the RHS of expressions (the values) - data <- purrr::map(user_code, - eval) - - # Grab the names - fields <- names(user_code) - - for(i in fields){ - model_template <- purrr::list_modify(model_template, !!!data[i]) - } - - return(model_template) +rmot_lm <- function(){ + list(X = NULL, + Y = NULL, + N = NULL, + model = "linear") } +#' Data configuration template for constant growth single species model +#' @keywords internal +#' @noRd + +rmot_cgs <- function(){ + list(N_obs = NULL, + N_ind = NULL, + S_obs = NULL, + census = NULL, + census_interval = NULL, + id_factor = NULL, + S_0_obs = NULL, + model = "constant_single") +} -list_rename = function(data, ...) { - mapping = sapply( - rlang::enquos(...), - rlang::as_name - ) - new_names = stats::setNames(nm=names(data)) - # `new_name = old_name` for consistency with `dplyr::rename` - new_names[mapping] = names(mapping) - # for `old_name = new_name` use: `new_names[names(mapping)] = mapping` - stats::setNames(data, new_names) -} diff --git a/R/rmot_run.R b/R/rmot_run.R new file mode 100644 index 0000000..16ed4d2 --- /dev/null +++ b/R/rmot_run.R @@ -0,0 +1,20 @@ +#' Run a linear model in Stan +#' +#' @param model_template model template generated by rmot_model and updated by rmot_assign_data +#' @param ... additional arguments passed to rstan::sampling +#' +#' @return Stanfit model output +#' @export +#' +#' @examples +#' mtcars +#' rmot_lm(mtcars$mpg, mtcars$disp) +rmot_run <- function(model_template, ...) { + + # Detect model + out <- switch(model_template$model, + linear = rstan::sampling(stanmodels$linear, data = model_template, ...), + constant_single = rstan::sampling(stanmodels$constant_single, data = model_template, ...)) + + return(out) +} diff --git a/R/stanmodels.R b/R/stanmodels.R index 78f7f7b..5a58cd8 100644 --- a/R/stanmodels.R +++ b/R/stanmodels.R @@ -1,10 +1,11 @@ # Generated by rstantools. Do not edit by hand. # names of stan models -stanmodels <- c("lm") +stanmodels <- c("constant_single", "linear") # load each stan module -Rcpp::loadModule("stan_fit4lm_mod", what = TRUE) +Rcpp::loadModule("stan_fit4constant_single_mod", what = TRUE) +Rcpp::loadModule("stan_fit4linear_mod", what = TRUE) # instantiate each stanmodel object stanmodels <- sapply(stanmodels, function(model_name) { diff --git a/R/utils.R b/R/utils.R new file mode 100644 index 0000000..8dee17a --- /dev/null +++ b/R/utils.R @@ -0,0 +1,11 @@ +list_rename = function(data, ...) { + mapping = sapply( + rlang::enquos(...), + rlang::as_name + ) + new_names = stats::setNames(nm=names(data)) + # `new_name = old_name` for consistency with `dplyr::rename` + new_names[mapping] = names(mapping) + # for `old_name = new_name` use: `new_names[names(mapping)] = mapping` + stats::setNames(data, new_names) +} diff --git a/inst/stan/constant_single.stan b/inst/stan/constant_single.stan new file mode 100644 index 0000000..24a984c --- /dev/null +++ b/inst/stan/constant_single.stan @@ -0,0 +1,83 @@ +//Constant Growth - Single species + +// Data structure +data { + int N_obs; + int N_ind; + real S_obs[N_obs]; + int census[N_obs]; + real census_interval[N_obs]; + int id_factor[N_obs]; + real S_0_obs[N_ind]; +} + +// The parameters accepted by the model. +parameters { + //Individual level + real ind_S_0[N_ind]; + real ind_beta[N_ind]; + + real species_beta_mu; + real species_beta_sigma; + + //Global level + real global_error_sigma; +} + +// The model to be estimated. +model { + real S_hat[N_obs]; + real G_hat[N_obs]; + + for(i in 1:N_obs){ + if(id_factor[i+1]==id_factor[i]){ + if(census[i]==1){//Fits the first size + S_hat[i] = ind_S_0[id_factor[i]]; + } + + if(i < N_obs){ //Analytic solution + G_hat[i] = ind_beta[id_factor[i]]; + S_hat[i+1] = S_hat[i] + G_hat[i]*census_interval[i+1]; + } + } else { + G_hat[i] = 0; //Gives 0 as the growth estimate for the last data point. + } + } + + //Likelihood + S_obs ~ normal(S_hat, global_error_sigma); + + //Priors + //Individual level + ind_S_0 ~ normal(S_0_obs, global_error_sigma); + ind_beta ~ lognormal(species_beta_mu, + species_beta_sigma); + + //Species level + species_beta_mu ~ normal(0.1, 1); + species_beta_sigma ~cauchy(0.1, 1); + + //Global level + global_error_sigma ~cauchy(0.1, 1); +} + +// The output +generated quantities { + real S_hat[N_obs]; + real G_hat[N_obs]; + + for(i in 1:N_obs){ + if(id_factor[i+1]==id_factor[i]){ + if(census[i]==1){//Fits the first size + S_hat[i] = ind_S_0[id_factor[i]]; + } + + if(i < N_obs){ //Analytic solution + G_hat[i] = ind_beta[id_factor[i]]; + S_hat[i+1] = S_hat[i] + G_hat[i]*census_interval[i+1]; + } + } else { + G_hat[i] = 0; //Gives 0 as the growth estimate for the last data point. + } + } +} diff --git a/inst/stan/lm.stan b/inst/stan/linear.stan similarity index 66% rename from inst/stan/lm.stan rename to inst/stan/linear.stan index 7ff45f3..5291f5e 100644 --- a/inst/stan/lm.stan +++ b/inst/stan/linear.stan @@ -1,8 +1,8 @@ // A linear model data { int N; - vector[N] x; - vector[N] y; + vector[N] X; + vector[N] Y; } parameters { real intercept; @@ -12,5 +12,5 @@ parameters { model { // ... priors, etc. - y ~ normal(intercept + beta * x, sigma); + Y ~ normal(intercept + beta * X, sigma); } diff --git a/man/rmot_assign_data.Rd b/man/rmot_assign_data.Rd new file mode 100644 index 0000000..23857f6 --- /dev/null +++ b/man/rmot_assign_data.Rd @@ -0,0 +1,22 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/rmot_assign_data.R +\name{rmot_assign_data} +\alias{rmot_assign_data} +\title{Assign data to template} +\usage{ +rmot_assign_data(model_template, ...) +} +\arguments{ +\item{model_template}{output from rmot_model} + +\item{...}{data-masking name-value pairs} +} +\value{ +updated named list with your data assigned to Stan model parameters +} +\description{ +Assign data to template +} +\examples{ +rmot_assign_data(X = Loblolly$age, Y = Loblolly$height) +} diff --git a/man/rmot_model.Rd b/man/rmot_model.Rd new file mode 100644 index 0000000..a20dd53 --- /dev/null +++ b/man/rmot_model.Rd @@ -0,0 +1,20 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/rmot_models.R +\name{rmot_model} +\alias{rmot_model} +\title{Select data configuration template for rmot supported model} +\usage{ +rmot_model(model = NULL) +} +\arguments{ +\item{model}{model name character string} +} +\value{ +named list that matches Stan model parameters +} +\description{ +Select data configuration template for rmot supported model +} +\examples{ +rmot_model("linear") +} diff --git a/man/rmot_lm.Rd b/man/rmot_run.Rd similarity index 54% rename from man/rmot_lm.Rd rename to man/rmot_run.Rd index aac125b..c399654 100644 --- a/man/rmot_lm.Rd +++ b/man/rmot_run.Rd @@ -1,20 +1,18 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/rmot_lm.R -\name{rmot_lm} -\alias{rmot_lm} +% Please edit documentation in R/rmot_run.R +\name{rmot_run} +\alias{rmot_run} \title{Run a linear model in Stan} \usage{ -rmot_lm(x, y, ...) +rmot_run(model_template, ...) } \arguments{ -\item{x}{vector of numeric or integer predictor variable} - -\item{y}{vector of numeric or integer response variable} +\item{model_template}{model template generated by rmot_model and updated by rmot_assign_data} \item{...}{additional arguments passed to rstan::sampling} } \value{ -Stan model output +Stanfit model output } \description{ Run a linear model in Stan diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 18a3e48..c924a64 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -12,10 +12,12 @@ Rcpp::Rostream& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get(); #endif -RcppExport SEXP _rcpp_module_boot_stan_fit4lm_mod(); +RcppExport SEXP _rcpp_module_boot_stan_fit4constant_single_mod(); +RcppExport SEXP _rcpp_module_boot_stan_fit4linear_mod(); static const R_CallMethodDef CallEntries[] = { - {"_rcpp_module_boot_stan_fit4lm_mod", (DL_FUNC) &_rcpp_module_boot_stan_fit4lm_mod, 0}, + {"_rcpp_module_boot_stan_fit4constant_single_mod", (DL_FUNC) &_rcpp_module_boot_stan_fit4constant_single_mod, 0}, + {"_rcpp_module_boot_stan_fit4linear_mod", (DL_FUNC) &_rcpp_module_boot_stan_fit4linear_mod, 0}, {NULL, NULL, 0} }; diff --git a/src/stanExports_constant_single.cc b/src/stanExports_constant_single.cc new file mode 100644 index 0000000..961bd38 --- /dev/null +++ b/src/stanExports_constant_single.cc @@ -0,0 +1,32 @@ +// Generated by rstantools. Do not edit by hand. + +#include +using namespace Rcpp ; +#include "stanExports_constant_single.h" + +RCPP_MODULE(stan_fit4constant_single_mod) { + + + class_ >("rstantools_model_constant_single") + + .constructor() + + + .method("call_sampler", &rstan::stan_fit ::call_sampler) + .method("param_names", &rstan::stan_fit ::param_names) + .method("param_names_oi", &rstan::stan_fit ::param_names_oi) + .method("param_fnames_oi", &rstan::stan_fit ::param_fnames_oi) + .method("param_dims", &rstan::stan_fit ::param_dims) + .method("param_dims_oi", &rstan::stan_fit ::param_dims_oi) + .method("update_param_oi", &rstan::stan_fit ::update_param_oi) + .method("param_oi_tidx", &rstan::stan_fit ::param_oi_tidx) + .method("grad_log_prob", &rstan::stan_fit ::grad_log_prob) + .method("log_prob", &rstan::stan_fit ::log_prob) + .method("unconstrain_pars", &rstan::stan_fit ::unconstrain_pars) + .method("constrain_pars", &rstan::stan_fit ::constrain_pars) + .method("num_pars_unconstrained", &rstan::stan_fit ::num_pars_unconstrained) + .method("unconstrained_param_names", &rstan::stan_fit ::unconstrained_param_names) + .method("constrained_param_names", &rstan::stan_fit ::constrained_param_names) + .method("standalone_gqs", &rstan::stan_fit ::standalone_gqs) + ; +} diff --git a/src/stanExports_constant_single.h b/src/stanExports_constant_single.h new file mode 100644 index 0000000..a617591 --- /dev/null +++ b/src/stanExports_constant_single.h @@ -0,0 +1,812 @@ +// Generated by rstantools. Do not edit by hand. + +/* + rmot 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. + + rmot 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 rmot. If not, see . +*/ +#ifndef MODELS_HPP +#define MODELS_HPP +#define STAN__SERVICES__COMMAND_HPP +#ifndef USE_STANC3 +#define USE_STANC3 +#endif +#include +// Code generated by stanc v2.26.1-4-gd72b68b7-dirty +#include +namespace model_constant_single_namespace { +inline void validate_positive_index(const char* var_name, const char* expr, + int val) { + if (val < 1) { + std::stringstream msg; + msg << "Found dimension size less than one in simplex declaration" + << "; variable=" << var_name << "; dimension size expression=" << expr + << "; expression value=" << val; + std::string msg_str(msg.str()); + throw std::invalid_argument(msg_str.c_str()); + } +} +inline void validate_unit_vector_index(const char* var_name, const char* expr, + int val) { + if (val <= 1) { + std::stringstream msg; + if (val == 1) { + msg << "Found dimension size one in unit vector declaration." + << " One-dimensional unit vector is discrete" + << " but the target distribution must be continuous." + << " variable=" << var_name << "; dimension size expression=" << expr; + } else { + msg << "Found dimension size less than one in unit vector declaration" + << "; variable=" << var_name << "; dimension size expression=" << expr + << "; expression value=" << val; + } + std::string msg_str(msg.str()); + throw std::invalid_argument(msg_str.c_str()); + } +} +using std::istream; +using std::string; +using std::stringstream; +using std::vector; +using std::pow; +using stan::io::dump; +using stan::math::lgamma; +using stan::model::model_base_crtp; +using stan::model::rvalue; +using stan::model::cons_list; +using stan::model::index_uni; +using stan::model::index_max; +using stan::model::index_min; +using stan::model::index_min_max; +using stan::model::index_multi; +using stan::model::index_omni; +using stan::model::nil_index_list; +using namespace stan::math; +using stan::math::pow; +stan::math::profile_map profiles__; +static int current_statement__= 0; +static const std::vector locations_array__ = {" (found before start of program)", + " (in 'constant_single', line 15, column 2 to column 31)", + " (in 'constant_single', line 16, column 2 to column 32)", + " (in 'constant_single', line 17, column 2 to column 23)", + " (in 'constant_single', line 18, column 2 to column 35)", + " (in 'constant_single', line 20, column 2 to column 35)", + " (in 'constant_single', line 54, column 2 to column 20)", + " (in 'constant_single', line 55, column 2 to column 20)", + " (in 'constant_single', line 66, column 6 to column 19)", + " (in 'constant_single', line 65, column 11 to line 67, column 5)", + " (in 'constant_single', line 59, column 8 to column 41)", + " (in 'constant_single', line 58, column 22 to line 60, column 7)", + " (in 'constant_single', line 58, column 6 to line 60, column 7)", + " (in 'constant_single', line 62, column 8 to column 42)", + " (in 'constant_single', line 63, column 8 to column 62)", + " (in 'constant_single', line 61, column 19 to line 64, column 7)", + " (in 'constant_single', line 61, column 6 to line 64, column 7)", + " (in 'constant_single', line 57, column 36 to line 65, column 5)", + " (in 'constant_single', line 57, column 4 to line 67, column 5)", + " (in 'constant_single', line 56, column 19 to line 68, column 3)", + " (in 'constant_single', line 56, column 2 to line 68, column 3)", + " (in 'constant_single', line 24, column 13 to column 18)", + " (in 'constant_single', line 24, column 2 to column 20)", + " (in 'constant_single', line 25, column 13 to column 18)", + " (in 'constant_single', line 25, column 2 to column 20)", + " (in 'constant_single', line 36, column 6 to column 19)", + " (in 'constant_single', line 35, column 11 to line 37, column 5)", + " (in 'constant_single', line 29, column 8 to column 41)", + " (in 'constant_single', line 28, column 22 to line 30, column 7)", + " (in 'constant_single', line 28, column 6 to line 30, column 7)", + " (in 'constant_single', line 32, column 8 to column 42)", + " (in 'constant_single', line 33, column 8 to column 62)", + " (in 'constant_single', line 31, column 19 to line 34, column 7)", + " (in 'constant_single', line 31, column 6 to line 34, column 7)", + " (in 'constant_single', line 27, column 36 to line 35, column 5)", + " (in 'constant_single', line 27, column 4 to line 37, column 5)", + " (in 'constant_single', line 26, column 19 to line 38, column 3)", + " (in 'constant_single', line 26, column 2 to line 38, column 3)", + " (in 'constant_single', line 40, column 2 to column 44)", + " (in 'constant_single', line 43, column 2 to column 48)", + " (in 'constant_single', line 44, column 2 to line 45, column 40)", + " (in 'constant_single', line 47, column 2 to column 35)", + " (in 'constant_single', line 48, column 2 to column 37)", + " (in 'constant_single', line 50, column 2 to column 37)", + " (in 'constant_single', line 4, column 2 to column 12)", + " (in 'constant_single', line 5, column 2 to column 12)", + " (in 'constant_single', line 6, column 13 to column 18)", + " (in 'constant_single', line 6, column 2 to column 20)", + " (in 'constant_single', line 7, column 13 to column 18)", + " (in 'constant_single', line 7, column 2 to column 20)", + " (in 'constant_single', line 8, column 23 to column 28)", + " (in 'constant_single', line 8, column 2 to column 30)", + " (in 'constant_single', line 9, column 16 to column 21)", + " (in 'constant_single', line 9, column 2 to column 23)", + " (in 'constant_single', line 10, column 15 to column 20)", + " (in 'constant_single', line 10, column 2 to column 22)", + " (in 'constant_single', line 15, column 24 to column 29)", + " (in 'constant_single', line 16, column 25 to column 30)", + " (in 'constant_single', line 54, column 13 to column 18)", + " (in 'constant_single', line 55, column 13 to column 18)"}; +#include +class model_constant_single final : public model_base_crtp { +private: + int N_obs; + int N_ind; + std::vector S_obs; + std::vector census; + std::vector census_interval; + std::vector id_factor; + std::vector S_0_obs; + +public: + ~model_constant_single() { } + + inline std::string model_name() const final { return "model_constant_single"; } + inline std::vector model_compile_info() const noexcept { + return std::vector{"stanc_version = stanc3 v2.26.1-4-gd72b68b7-dirty", "stancflags = "}; + } + + + model_constant_single(stan::io::var_context& context__, + unsigned int random_seed__ = 0, + std::ostream* pstream__ = nullptr) : model_base_crtp(0) { + using local_scalar_t__ = double ; + boost::ecuyer1988 base_rng__ = + stan::services::util::create_rng(random_seed__, 0); + (void) base_rng__; // suppress unused var warning + static const char* function__ = "model_constant_single_namespace::model_constant_single"; + (void) function__; // suppress unused var warning + local_scalar_t__ DUMMY_VAR__(std::numeric_limits::quiet_NaN()); + (void) DUMMY_VAR__; // suppress unused var warning + + try { + int pos__; + pos__ = std::numeric_limits::min(); + + pos__ = 1; + current_statement__ = 44; + context__.validate_dims("data initialization","N_obs","int", + context__.to_vec()); + N_obs = std::numeric_limits::min(); + + current_statement__ = 44; + N_obs = context__.vals_i("N_obs")[(1 - 1)]; + current_statement__ = 45; + context__.validate_dims("data initialization","N_ind","int", + context__.to_vec()); + N_ind = std::numeric_limits::min(); + + current_statement__ = 45; + N_ind = context__.vals_i("N_ind")[(1 - 1)]; + current_statement__ = 46; + validate_non_negative_index("S_obs", "N_obs", N_obs); + current_statement__ = 47; + context__.validate_dims("data initialization","S_obs","double", + context__.to_vec(N_obs)); + S_obs = std::vector(N_obs, std::numeric_limits::quiet_NaN()); + + current_statement__ = 47; + assign(S_obs, nil_index_list(), context__.vals_r("S_obs"), + "assigning variable S_obs"); + current_statement__ = 48; + validate_non_negative_index("census", "N_obs", N_obs); + current_statement__ = 49; + context__.validate_dims("data initialization","census","int", + context__.to_vec(N_obs)); + census = std::vector(N_obs, std::numeric_limits::min()); + + current_statement__ = 49; + assign(census, nil_index_list(), context__.vals_i("census"), + "assigning variable census"); + current_statement__ = 50; + validate_non_negative_index("census_interval", "N_obs", N_obs); + current_statement__ = 51; + context__.validate_dims("data initialization","census_interval", + "double",context__.to_vec(N_obs)); + census_interval = std::vector(N_obs, std::numeric_limits::quiet_NaN()); + + current_statement__ = 51; + assign(census_interval, nil_index_list(), + context__.vals_r("census_interval"), + "assigning variable census_interval"); + current_statement__ = 52; + validate_non_negative_index("id_factor", "N_obs", N_obs); + current_statement__ = 53; + context__.validate_dims("data initialization","id_factor","int", + context__.to_vec(N_obs)); + id_factor = std::vector(N_obs, std::numeric_limits::min()); + + current_statement__ = 53; + assign(id_factor, nil_index_list(), context__.vals_i("id_factor"), + "assigning variable id_factor"); + current_statement__ = 54; + validate_non_negative_index("S_0_obs", "N_ind", N_ind); + current_statement__ = 55; + context__.validate_dims("data initialization","S_0_obs","double", + context__.to_vec(N_ind)); + S_0_obs = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 55; + assign(S_0_obs, nil_index_list(), context__.vals_r("S_0_obs"), + "assigning variable S_0_obs"); + current_statement__ = 56; + validate_non_negative_index("ind_S_0", "N_ind", N_ind); + current_statement__ = 57; + validate_non_negative_index("ind_beta", "N_ind", N_ind); + current_statement__ = 58; + validate_non_negative_index("S_hat", "N_obs", N_obs); + current_statement__ = 59; + validate_non_negative_index("G_hat", "N_obs", N_obs); + } catch (const std::exception& e) { + stan::lang::rethrow_located(e, locations_array__[current_statement__]); + // Next line prevents compiler griping about no return + throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***"); + } + num_params_r__ = 0U; + + try { + num_params_r__ += N_ind; + num_params_r__ += N_ind; + num_params_r__ += 1; + num_params_r__ += 1; + num_params_r__ += 1; + } catch (const std::exception& e) { + stan::lang::rethrow_located(e, locations_array__[current_statement__]); + // Next line prevents compiler griping about no return + throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***"); + } + } + template * = nullptr, stan::require_vector_like_vt* = nullptr> + inline stan::scalar_type_t log_prob_impl(VecR& params_r__, + VecI& params_i__, + std::ostream* pstream__ = nullptr) const { + using T__ = stan::scalar_type_t; + using local_scalar_t__ = T__; + T__ lp__(0.0); + stan::math::accumulator lp_accum__; + static const char* function__ = "model_constant_single_namespace::log_prob"; +(void) function__; // suppress unused var warning + stan::io::reader in__(params_r__, params_i__); + local_scalar_t__ DUMMY_VAR__(std::numeric_limits::quiet_NaN()); + (void) DUMMY_VAR__; // suppress unused var warning + + try { + std::vector ind_S_0; + ind_S_0 = std::vector(N_ind, DUMMY_VAR__); + + current_statement__ = 1; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 1; + assign(ind_S_0, cons_list(index_uni(sym1__), nil_index_list()), + in__.scalar(), "assigning variable ind_S_0");} + current_statement__ = 1; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 1; + if (jacobian__) { + current_statement__ = 1; + assign(ind_S_0, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_S_0[(sym1__ - 1)], 0, lp__), + "assigning variable ind_S_0"); + } else { + current_statement__ = 1; + assign(ind_S_0, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_S_0[(sym1__ - 1)], 0), + "assigning variable ind_S_0"); + }} + std::vector ind_beta; + ind_beta = std::vector(N_ind, DUMMY_VAR__); + + current_statement__ = 2; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 2; + assign(ind_beta, cons_list(index_uni(sym1__), nil_index_list()), + in__.scalar(), "assigning variable ind_beta");} + current_statement__ = 2; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 2; + if (jacobian__) { + current_statement__ = 2; + assign(ind_beta, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_beta[(sym1__ - 1)], 0, lp__), + "assigning variable ind_beta"); + } else { + current_statement__ = 2; + assign(ind_beta, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_beta[(sym1__ - 1)], 0), + "assigning variable ind_beta"); + }} + local_scalar_t__ species_beta_mu; + species_beta_mu = DUMMY_VAR__; + + current_statement__ = 3; + species_beta_mu = in__.scalar(); + local_scalar_t__ species_beta_sigma; + species_beta_sigma = DUMMY_VAR__; + + current_statement__ = 4; + species_beta_sigma = in__.scalar(); + current_statement__ = 4; + if (jacobian__) { + current_statement__ = 4; + species_beta_sigma = stan::math::lb_constrain(species_beta_sigma, 0, + lp__); + } else { + current_statement__ = 4; + species_beta_sigma = stan::math::lb_constrain(species_beta_sigma, 0); + } + local_scalar_t__ global_error_sigma; + global_error_sigma = DUMMY_VAR__; + + current_statement__ = 5; + global_error_sigma = in__.scalar(); + current_statement__ = 5; + if (jacobian__) { + current_statement__ = 5; + global_error_sigma = stan::math::lb_constrain(global_error_sigma, 0, + lp__); + } else { + current_statement__ = 5; + global_error_sigma = stan::math::lb_constrain(global_error_sigma, 0); + } + { + current_statement__ = 21; + validate_non_negative_index("S_hat", "N_obs", N_obs); + std::vector S_hat; + S_hat = std::vector(N_obs, DUMMY_VAR__); + + current_statement__ = 23; + validate_non_negative_index("G_hat", "N_obs", N_obs); + std::vector G_hat; + G_hat = std::vector(N_obs, DUMMY_VAR__); + + current_statement__ = 37; + for (int i = 1; i <= N_obs; ++i) { + current_statement__ = 35; + if (logical_eq(id_factor[((i + 1) - 1)], id_factor[(i - 1)])) { + current_statement__ = 29; + if (logical_eq(census[(i - 1)], 1)) { + current_statement__ = 27; + assign(S_hat, cons_list(index_uni(i), nil_index_list()), + ind_S_0[(id_factor[(i - 1)] - 1)], "assigning variable S_hat"); + } + current_statement__ = 33; + if (logical_lt(i, N_obs)) { + current_statement__ = 30; + assign(G_hat, cons_list(index_uni(i), nil_index_list()), + ind_beta[(id_factor[(i - 1)] - 1)], + "assigning variable G_hat"); + current_statement__ = 31; + assign(S_hat, cons_list(index_uni((i + 1)), nil_index_list()), + (S_hat[(i - 1)] + + (G_hat[(i - 1)] * census_interval[((i + 1) - 1)])), + "assigning variable S_hat"); + } + } else { + current_statement__ = 25; + assign(G_hat, cons_list(index_uni(i), nil_index_list()), 0, + "assigning variable G_hat"); + }} + current_statement__ = 38; + lp_accum__.add( + normal_lpdf(S_obs, S_hat, global_error_sigma)); + current_statement__ = 39; + lp_accum__.add( + normal_lpdf(ind_S_0, S_0_obs, global_error_sigma)); + current_statement__ = 40; + lp_accum__.add( + lognormal_lpdf(ind_beta, species_beta_mu, + species_beta_sigma)); + current_statement__ = 41; + lp_accum__.add(normal_lpdf(species_beta_mu, 0.1, 1)); + current_statement__ = 42; + lp_accum__.add(cauchy_lpdf(species_beta_sigma, 0.1, 1)); + current_statement__ = 43; + lp_accum__.add(cauchy_lpdf(global_error_sigma, 0.1, 1)); + } + } catch (const std::exception& e) { + stan::lang::rethrow_located(e, locations_array__[current_statement__]); + // Next line prevents compiler griping about no return + throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***"); + } + lp_accum__.add(lp__); + return lp_accum__.sum(); + } // log_prob_impl() + + template * = nullptr, stan::require_vector_like_vt* = nullptr, stan::require_std_vector_vt* = nullptr> + inline void write_array_impl(RNG& base_rng__, VecR& params_r__, + VecI& params_i__, VecVar& vars__, + const bool emit_transformed_parameters__ = true, + const bool emit_generated_quantities__ = true, + std::ostream* pstream__ = nullptr) const { + using local_scalar_t__ = double; + vars__.resize(0); + stan::io::reader in__(params_r__, params_i__); + static const char* function__ = "model_constant_single_namespace::write_array"; +(void) function__; // suppress unused var warning + (void) function__; // suppress unused var warning + double lp__ = 0.0; + (void) lp__; // dummy to suppress unused var warning + stan::math::accumulator lp_accum__; + local_scalar_t__ DUMMY_VAR__(std::numeric_limits::quiet_NaN()); + (void) DUMMY_VAR__; // suppress unused var warning + + try { + std::vector ind_S_0; + ind_S_0 = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 1; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 1; + assign(ind_S_0, cons_list(index_uni(sym1__), nil_index_list()), + in__.scalar(), "assigning variable ind_S_0");} + current_statement__ = 1; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 1; + assign(ind_S_0, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_S_0[(sym1__ - 1)], 0), + "assigning variable ind_S_0");} + std::vector ind_beta; + ind_beta = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 2; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 2; + assign(ind_beta, cons_list(index_uni(sym1__), nil_index_list()), + in__.scalar(), "assigning variable ind_beta");} + current_statement__ = 2; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 2; + assign(ind_beta, cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_constrain(ind_beta[(sym1__ - 1)], 0), + "assigning variable ind_beta");} + double species_beta_mu; + species_beta_mu = std::numeric_limits::quiet_NaN(); + + current_statement__ = 3; + species_beta_mu = in__.scalar(); + double species_beta_sigma; + species_beta_sigma = std::numeric_limits::quiet_NaN(); + + current_statement__ = 4; + species_beta_sigma = in__.scalar(); + current_statement__ = 4; + species_beta_sigma = stan::math::lb_constrain(species_beta_sigma, 0); + double global_error_sigma; + global_error_sigma = std::numeric_limits::quiet_NaN(); + + current_statement__ = 5; + global_error_sigma = in__.scalar(); + current_statement__ = 5; + global_error_sigma = stan::math::lb_constrain(global_error_sigma, 0); + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + vars__.emplace_back(ind_S_0[(sym1__ - 1)]);} + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + vars__.emplace_back(ind_beta[(sym1__ - 1)]);} + vars__.emplace_back(species_beta_mu); + vars__.emplace_back(species_beta_sigma); + vars__.emplace_back(global_error_sigma); + if (logical_negation((primitive_value(emit_transformed_parameters__) || + primitive_value(emit_generated_quantities__)))) { + return ; + } + if (logical_negation(emit_generated_quantities__)) { + return ; + } + std::vector S_hat; + S_hat = std::vector(N_obs, std::numeric_limits::quiet_NaN()); + + std::vector G_hat; + G_hat = std::vector(N_obs, std::numeric_limits::quiet_NaN()); + + current_statement__ = 20; + for (int i = 1; i <= N_obs; ++i) { + current_statement__ = 18; + if (logical_eq(id_factor[((i + 1) - 1)], id_factor[(i - 1)])) { + current_statement__ = 12; + if (logical_eq(census[(i - 1)], 1)) { + current_statement__ = 10; + assign(S_hat, cons_list(index_uni(i), nil_index_list()), + ind_S_0[(id_factor[(i - 1)] - 1)], "assigning variable S_hat"); + } + current_statement__ = 16; + if (logical_lt(i, N_obs)) { + current_statement__ = 13; + assign(G_hat, cons_list(index_uni(i), nil_index_list()), + ind_beta[(id_factor[(i - 1)] - 1)], "assigning variable G_hat"); + current_statement__ = 14; + assign(S_hat, cons_list(index_uni((i + 1)), nil_index_list()), + (S_hat[(i - 1)] + + (G_hat[(i - 1)] * census_interval[((i + 1) - 1)])), + "assigning variable S_hat"); + } + } else { + current_statement__ = 8; + assign(G_hat, cons_list(index_uni(i), nil_index_list()), 0, + "assigning variable G_hat"); + }} + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + vars__.emplace_back(S_hat[(sym1__ - 1)]);} + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + vars__.emplace_back(G_hat[(sym1__ - 1)]);} + } catch (const std::exception& e) { + stan::lang::rethrow_located(e, locations_array__[current_statement__]); + // Next line prevents compiler griping about no return + throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***"); + } + } // write_array_impl() + + template * = nullptr, stan::require_vector_like_vt* = nullptr> + inline void transform_inits_impl(const stan::io::var_context& context__, + VecI& params_i__, VecVar& vars__, + std::ostream* pstream__ = nullptr) const { + using local_scalar_t__ = double; + vars__.clear(); + vars__.reserve(num_params_r__); + + try { + int pos__; + pos__ = std::numeric_limits::min(); + + pos__ = 1; + std::vector ind_S_0; + ind_S_0 = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 1; + assign(ind_S_0, nil_index_list(), context__.vals_r("ind_S_0"), + "assigning variable ind_S_0"); + std::vector ind_S_0_free__; + ind_S_0_free__ = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 1; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 1; + assign(ind_S_0_free__, + cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_free(ind_S_0[(sym1__ - 1)], 0), + "assigning variable ind_S_0_free__");} + std::vector ind_beta; + ind_beta = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 2; + assign(ind_beta, nil_index_list(), context__.vals_r("ind_beta"), + "assigning variable ind_beta"); + std::vector ind_beta_free__; + ind_beta_free__ = std::vector(N_ind, std::numeric_limits::quiet_NaN()); + + current_statement__ = 2; + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + current_statement__ = 2; + assign(ind_beta_free__, + cons_list(index_uni(sym1__), nil_index_list()), + stan::math::lb_free(ind_beta[(sym1__ - 1)], 0), + "assigning variable ind_beta_free__");} + double species_beta_mu; + species_beta_mu = std::numeric_limits::quiet_NaN(); + + current_statement__ = 3; + species_beta_mu = context__.vals_r("species_beta_mu")[(1 - 1)]; + double species_beta_sigma; + species_beta_sigma = std::numeric_limits::quiet_NaN(); + + current_statement__ = 4; + species_beta_sigma = context__.vals_r("species_beta_sigma")[(1 - 1)]; + double species_beta_sigma_free__; + species_beta_sigma_free__ = std::numeric_limits::quiet_NaN(); + + current_statement__ = 4; + species_beta_sigma_free__ = stan::math::lb_free(species_beta_sigma, 0); + double global_error_sigma; + global_error_sigma = std::numeric_limits::quiet_NaN(); + + current_statement__ = 5; + global_error_sigma = context__.vals_r("global_error_sigma")[(1 - 1)]; + double global_error_sigma_free__; + global_error_sigma_free__ = std::numeric_limits::quiet_NaN(); + + current_statement__ = 5; + global_error_sigma_free__ = stan::math::lb_free(global_error_sigma, 0); + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + vars__.emplace_back(ind_S_0_free__[(sym1__ - 1)]);} + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + vars__.emplace_back(ind_beta_free__[(sym1__ - 1)]);} + vars__.emplace_back(species_beta_mu); + vars__.emplace_back(species_beta_sigma_free__); + vars__.emplace_back(global_error_sigma_free__); + } catch (const std::exception& e) { + stan::lang::rethrow_located(e, locations_array__[current_statement__]); + // Next line prevents compiler griping about no return + throw std::runtime_error("*** IF YOU SEE THIS, PLEASE REPORT A BUG ***"); + } + } // transform_inits_impl() + + inline void get_param_names(std::vector& names__) const { + + names__.clear(); + names__.emplace_back("ind_S_0"); + names__.emplace_back("ind_beta"); + names__.emplace_back("species_beta_mu"); + names__.emplace_back("species_beta_sigma"); + names__.emplace_back("global_error_sigma"); + names__.emplace_back("S_hat"); + names__.emplace_back("G_hat"); + } // get_param_names() + + inline void get_dims(std::vector>& dimss__) const { + dimss__.clear(); + dimss__.emplace_back(std::vector{static_cast(N_ind)}); + + dimss__.emplace_back(std::vector{static_cast(N_ind)}); + + dimss__.emplace_back(std::vector{}); + + dimss__.emplace_back(std::vector{}); + + dimss__.emplace_back(std::vector{}); + + dimss__.emplace_back(std::vector{static_cast(N_obs)}); + + dimss__.emplace_back(std::vector{static_cast(N_obs)}); + + } // get_dims() + + inline void constrained_param_names( + std::vector& param_names__, + bool emit_transformed_parameters__ = true, + bool emit_generated_quantities__ = true) const + final { + + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + { + param_names__.emplace_back(std::string() + "ind_S_0" + '.' + std::to_string(sym1__)); + }} + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + { + param_names__.emplace_back(std::string() + "ind_beta" + '.' + std::to_string(sym1__)); + }} + param_names__.emplace_back(std::string() + "species_beta_mu"); + param_names__.emplace_back(std::string() + "species_beta_sigma"); + param_names__.emplace_back(std::string() + "global_error_sigma"); + if (emit_transformed_parameters__) { + + } + + if (emit_generated_quantities__) { + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + { + param_names__.emplace_back(std::string() + "S_hat" + '.' + std::to_string(sym1__)); + }} + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + { + param_names__.emplace_back(std::string() + "G_hat" + '.' + std::to_string(sym1__)); + }} + } + + } // constrained_param_names() + + inline void unconstrained_param_names( + std::vector& param_names__, + bool emit_transformed_parameters__ = true, + bool emit_generated_quantities__ = true) const + final { + + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + { + param_names__.emplace_back(std::string() + "ind_S_0" + '.' + std::to_string(sym1__)); + }} + for (int sym1__ = 1; sym1__ <= N_ind; ++sym1__) { + { + param_names__.emplace_back(std::string() + "ind_beta" + '.' + std::to_string(sym1__)); + }} + param_names__.emplace_back(std::string() + "species_beta_mu"); + param_names__.emplace_back(std::string() + "species_beta_sigma"); + param_names__.emplace_back(std::string() + "global_error_sigma"); + if (emit_transformed_parameters__) { + + } + + if (emit_generated_quantities__) { + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + { + param_names__.emplace_back(std::string() + "S_hat" + '.' + std::to_string(sym1__)); + }} + for (int sym1__ = 1; sym1__ <= N_obs; ++sym1__) { + { + param_names__.emplace_back(std::string() + "G_hat" + '.' + std::to_string(sym1__)); + }} + } + + } // unconstrained_param_names() + + inline std::string get_constrained_sizedtypes() const { + stringstream s__; + s__ << "[{\"name\":\"ind_S_0\",\"type\":{\"name\":\"array\",\"length\":" << N_ind << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"parameters\"},{\"name\":\"ind_beta\",\"type\":{\"name\":\"array\",\"length\":" << N_ind << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"parameters\"},{\"name\":\"species_beta_mu\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"species_beta_sigma\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"global_error_sigma\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"S_hat\",\"type\":{\"name\":\"array\",\"length\":" << N_obs << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"generated_quantities\"},{\"name\":\"G_hat\",\"type\":{\"name\":\"array\",\"length\":" << N_obs << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"generated_quantities\"}]"; + return s__.str(); + } // get_constrained_sizedtypes() + + inline std::string get_unconstrained_sizedtypes() const { + stringstream s__; + s__ << "[{\"name\":\"ind_S_0\",\"type\":{\"name\":\"array\",\"length\":" << N_ind << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"parameters\"},{\"name\":\"ind_beta\",\"type\":{\"name\":\"array\",\"length\":" << N_ind << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"parameters\"},{\"name\":\"species_beta_mu\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"species_beta_sigma\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"global_error_sigma\",\"type\":{\"name\":\"real\"},\"block\":\"parameters\"},{\"name\":\"S_hat\",\"type\":{\"name\":\"array\",\"length\":" << N_obs << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"generated_quantities\"},{\"name\":\"G_hat\",\"type\":{\"name\":\"array\",\"length\":" << N_obs << ",\"element_type\":{\"name\":\"real\"}},\"block\":\"generated_quantities\"}]"; + return s__.str(); + } // get_unconstrained_sizedtypes() + + + // Begin method overload boilerplate + template + inline void write_array(RNG& base_rng, + Eigen::Matrix& params_r, + Eigen::Matrix& vars, + const bool emit_transformed_parameters = true, + const bool emit_generated_quantities = true, + std::ostream* pstream = nullptr) const { + std::vector vars_vec(vars.size()); + std::vector params_i; + write_array_impl(base_rng, params_r, params_i, vars_vec, + emit_transformed_parameters, emit_generated_quantities, pstream); + vars.resize(vars_vec.size()); + for (int i = 0; i < vars.size(); ++i) { + vars.coeffRef(i) = vars_vec[i]; + } + } + template + inline void write_array(RNG& base_rng, std::vector& params_r, + std::vector& params_i, + std::vector& vars, + bool emit_transformed_parameters = true, + bool emit_generated_quantities = true, + std::ostream* pstream = nullptr) const { + write_array_impl(base_rng, params_r, params_i, vars, emit_transformed_parameters, emit_generated_quantities, pstream); + } + template + inline T_ log_prob(Eigen::Matrix& params_r, + std::ostream* pstream = nullptr) const { + Eigen::Matrix params_i; + return log_prob_impl(params_r, params_i, pstream); + } + template + inline T__ log_prob(std::vector& params_r, + std::vector& params_i, + std::ostream* pstream = nullptr) const { + return log_prob_impl(params_r, params_i, pstream); + } + + inline void transform_inits(const stan::io::var_context& context, + Eigen::Matrix& params_r, + std::ostream* pstream = nullptr) const final { + std::vector params_r_vec(params_r.size()); + std::vector params_i; + transform_inits_impl(context, params_i, params_r_vec, pstream); + params_r.resize(params_r_vec.size()); + for (int i = 0; i < params_r.size(); ++i) { + params_r.coeffRef(i) = params_r_vec[i]; + } + } + inline void transform_inits(const stan::io::var_context& context, + std::vector& params_i, + std::vector& vars, + std::ostream* pstream = nullptr) const final { + transform_inits_impl(context, params_i, vars, pstream); + } +}; +} +using stan_model = model_constant_single_namespace::model_constant_single; +#ifndef USING_R +// Boilerplate +stan::model::model_base& new_model( + stan::io::var_context& data_context, + unsigned int seed, + std::ostream* msg_stream) { + stan_model* m = new stan_model(data_context, seed, msg_stream); + return *m; +} +stan::math::profile_map& get_stan_profile_data() { + return model_constant_single_namespace::profiles__; +} +#endif +#endif diff --git a/src/stanExports_lm.cc b/src/stanExports_linear.cc similarity index 95% rename from src/stanExports_lm.cc rename to src/stanExports_linear.cc index 873a18d..a6c0a65 100644 --- a/src/stanExports_lm.cc +++ b/src/stanExports_linear.cc @@ -2,12 +2,12 @@ #include using namespace Rcpp ; -#include "stanExports_lm.h" +#include "stanExports_linear.h" -RCPP_MODULE(stan_fit4lm_mod) { +RCPP_MODULE(stan_fit4linear_mod) { - class_ >("rstantools_model_lm") + class_ >("rstantools_model_linear") .constructor() diff --git a/src/stanExports_lm.h b/src/stanExports_linear.h similarity index 86% rename from src/stanExports_lm.h rename to src/stanExports_linear.h index 8ea21b3..b072723 100644 --- a/src/stanExports_lm.h +++ b/src/stanExports_linear.h @@ -23,7 +23,7 @@ #include // Code generated by stanc v2.26.1-4-gd72b68b7-dirty #include -namespace model_lm_namespace { +namespace model_linear_namespace { inline void validate_positive_index(const char* var_name, const char* expr, int val) { if (val < 1) { @@ -75,38 +75,39 @@ using stan::math::pow; stan::math::profile_map profiles__; static int current_statement__= 0; static const std::vector locations_array__ = {" (found before start of program)", - " (in 'lm', line 8, column 2 to column 17)", - " (in 'lm', line 9, column 2 to column 12)", - " (in 'lm', line 10, column 2 to column 22)", - " (in 'lm', line 14, column 2 to column 42)", - " (in 'lm', line 3, column 2 to column 17)", - " (in 'lm', line 4, column 9 to column 10)", - " (in 'lm', line 4, column 2 to column 14)", - " (in 'lm', line 5, column 9 to column 10)", - " (in 'lm', line 5, column 2 to column 14)"}; + " (in 'linear', line 8, column 2 to column 17)", + " (in 'linear', line 9, column 2 to column 12)", + " (in 'linear', line 10, column 2 to column 22)", + " (in 'linear', line 14, column 2 to column 42)", + " (in 'linear', line 3, column 2 to column 17)", + " (in 'linear', line 4, column 9 to column 10)", + " (in 'linear', line 4, column 2 to column 14)", + " (in 'linear', line 5, column 9 to column 10)", + " (in 'linear', line 5, column 2 to column 14)"}; #include -class model_lm final : public model_base_crtp { +class model_linear final : public model_base_crtp { private: int N; - Eigen::Matrix x; - Eigen::Matrix y; + Eigen::Matrix X; + Eigen::Matrix Y; public: - ~model_lm() { } + ~model_linear() { } - inline std::string model_name() const final { return "model_lm"; } + inline std::string model_name() const final { return "model_linear"; } inline std::vector model_compile_info() const noexcept { return std::vector{"stanc_version = stanc3 v2.26.1-4-gd72b68b7-dirty", "stancflags = "}; } - model_lm(stan::io::var_context& context__, unsigned int random_seed__ = 0, - std::ostream* pstream__ = nullptr) : model_base_crtp(0) { + model_linear(stan::io::var_context& context__, + unsigned int random_seed__ = 0, + std::ostream* pstream__ = nullptr) : model_base_crtp(0) { using local_scalar_t__ = double ; boost::ecuyer1988 base_rng__ = stan::services::util::create_rng(random_seed__, 0); (void) base_rng__; // suppress unused var warning - static const char* function__ = "model_lm_namespace::model_lm"; + static const char* function__ = "model_linear_namespace::model_linear"; (void) function__; // suppress unused var warning local_scalar_t__ DUMMY_VAR__(std::numeric_limits::quiet_NaN()); (void) DUMMY_VAR__; // suppress unused var warning @@ -127,48 +128,48 @@ class model_lm final : public model_base_crtp { current_statement__ = 5; check_greater_or_equal(function__, "N", N, 1); current_statement__ = 6; - validate_non_negative_index("x", "N", N); + validate_non_negative_index("X", "N", N); current_statement__ = 7; - context__.validate_dims("data initialization","x","double", + context__.validate_dims("data initialization","X","double", context__.to_vec(N)); - x = Eigen::Matrix(N); - stan::math::fill(x, std::numeric_limits::quiet_NaN()); + X = Eigen::Matrix(N); + stan::math::fill(X, std::numeric_limits::quiet_NaN()); { - std::vector x_flat__; + std::vector X_flat__; current_statement__ = 7; - assign(x_flat__, nil_index_list(), context__.vals_r("x"), - "assigning variable x_flat__"); + assign(X_flat__, nil_index_list(), context__.vals_r("X"), + "assigning variable X_flat__"); current_statement__ = 7; pos__ = 1; current_statement__ = 7; for (int sym1__ = 1; sym1__ <= N; ++sym1__) { current_statement__ = 7; - assign(x, cons_list(index_uni(sym1__), nil_index_list()), - x_flat__[(pos__ - 1)], "assigning variable x"); + assign(X, cons_list(index_uni(sym1__), nil_index_list()), + X_flat__[(pos__ - 1)], "assigning variable X"); current_statement__ = 7; pos__ = (pos__ + 1);} } current_statement__ = 8; - validate_non_negative_index("y", "N", N); + validate_non_negative_index("Y", "N", N); current_statement__ = 9; - context__.validate_dims("data initialization","y","double", + context__.validate_dims("data initialization","Y","double", context__.to_vec(N)); - y = Eigen::Matrix(N); - stan::math::fill(y, std::numeric_limits::quiet_NaN()); + Y = Eigen::Matrix(N); + stan::math::fill(Y, std::numeric_limits::quiet_NaN()); { - std::vector y_flat__; + std::vector Y_flat__; current_statement__ = 9; - assign(y_flat__, nil_index_list(), context__.vals_r("y"), - "assigning variable y_flat__"); + assign(Y_flat__, nil_index_list(), context__.vals_r("Y"), + "assigning variable Y_flat__"); current_statement__ = 9; pos__ = 1; current_statement__ = 9; for (int sym1__ = 1; sym1__ <= N; ++sym1__) { current_statement__ = 9; - assign(y, cons_list(index_uni(sym1__), nil_index_list()), - y_flat__[(pos__ - 1)], "assigning variable y"); + assign(Y, cons_list(index_uni(sym1__), nil_index_list()), + Y_flat__[(pos__ - 1)], "assigning variable Y"); current_statement__ = 9; pos__ = (pos__ + 1);} } @@ -197,7 +198,7 @@ class model_lm final : public model_base_crtp { using local_scalar_t__ = T__; T__ lp__(0.0); stan::math::accumulator lp_accum__; - static const char* function__ = "model_lm_namespace::log_prob"; + static const char* function__ = "model_linear_namespace::log_prob"; (void) function__; // suppress unused var warning stan::io::reader in__(params_r__, params_i__); local_scalar_t__ DUMMY_VAR__(std::numeric_limits::quiet_NaN()); @@ -230,7 +231,7 @@ class model_lm final : public model_base_crtp { { current_statement__ = 4; lp_accum__.add( - normal_lpdf(y, add(intercept, multiply(beta, x)), sigma)); + normal_lpdf(Y, add(intercept, multiply(beta, X)), sigma)); } } catch (const std::exception& e) { stan::lang::rethrow_located(e, locations_array__[current_statement__]); @@ -250,7 +251,7 @@ class model_lm final : public model_base_crtp { using local_scalar_t__ = double; vars__.resize(0); stan::io::reader in__(params_r__, params_i__); - static const char* function__ = "model_lm_namespace::write_array"; + static const char* function__ = "model_linear_namespace::write_array"; (void) function__; // suppress unused var warning (void) function__; // suppress unused var warning double lp__ = 0.0; @@ -464,7 +465,7 @@ class model_lm final : public model_base_crtp { } }; } -using stan_model = model_lm_namespace::model_lm; +using stan_model = model_linear_namespace::model_linear; #ifndef USING_R // Boilerplate stan::model::model_base& new_model( @@ -475,7 +476,7 @@ stan::model::model_base& new_model( return *m; } stan::math::profile_map& get_stan_profile_data() { - return model_lm_namespace::profiles__; + return model_linear_namespace::profiles__; } #endif #endif