Skip to content

Instantly share code, notes, and snippets.

@seabbs

seabbs/pomp-tb-model

Created August 8, 2019 11:05
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 1 You must be signed in to fork a gist
  • Save seabbs/f08a8a46b1342b8649df963ac015ea31 to your computer and use it in GitHub Desktop.
Save seabbs/f08a8a46b1342b8649df963ac015ea31 to your computer and use it in GitHub Desktop.
Download ZIP
An example POMP TB transmission model
Raw
pomp-tb-model
#' Disease Model, with Vaccination, Data Based Population Demographics and Ageing as a Pomp Object
#'
#' @param df A data frame of data to be passed to the pomp model.
#' @param t0 An integer value denoting the starting time for the model, this may be negative.
#' @param nstageA The number of age compartments to model.
#' @param timestep Value indicating the timestep over which to solve the model.
#' @param time The variable to use the time parameter from the data.
#' @param covars Additional data to use within the model.
#' @param tcovar A character string indicating the time variable in covar.
#' @param obsnames A character string of the observable variables.
#' @param ... Pass additional arguements to pomp function.
#' @import magrittr
#' @import pomp
#' @importFrom dplyr data_frame
#' @importFrom idmodelr add_pointer_struct
#' @return A pomp model object.
#' @export
#'
#' @examples
#'
tb_model <- function(df, t0 = 0, nstageA = 18, timestep = 1/365, time = "year",
covars, tcovar = "time", obsnames = paste0("incidence", 1:18), ...) {
rskel <- Csnippet("
// Set up populations
const double *S = &S1; const double *E_1 = &E_11;
const double *E_2 = &E_21; const double *I = &I1;
const double *E_R = &E_R1;
double *DS = &DS1; double *DE_1 = &DE_11;
double *DE_2 = &DE_21; double *DI = &DI1;
double *DE_R = &DE_R1;
const double *S_v = &S_v1; const double *E_1v = &E_1v1;
const double *E_2v = &E_2v1; const double *I_v = &I_v1;
const double *E_Rv = &E_Rv1;
double *DS_v = &DS_v1; double *DE_1v = &DE_1v1;
double *DE_2v = &DE_2v1; double *DI_v = &DI_v1;
double *DE_Rv = &DE_Rv1;
// Set up tracking populations
double *DN_E_1I = &DN_E_1I1; double *DN_E_2I = &DN_E_2I1;
double *DN_E_RI = &DN_E_RI1;
double *DN_E_1I_v = &DN_E_1I_v1; double *DN_E_2I_v = &DN_E_2I_v1;
double *DN_E_RI_v = &DN_E_RI_v1;
double DDeaths = &DDeaths;
// Model parameters
double *gamma = &gamma1;
double *alpha = &alpha1;
const double *theta = &theta1;
const double *C = &C1;
const double *iota = &iota1;
const double *mu = &mu1;
const double *delta = &delta1;
const double *chi = &chi1;
//Function variables
int i; int j;
double trans[18];
double foi;
double mu_t;
double epsilon_1; double epsilon_2; double kappa;
// Override scenario definition if using school age override, force to school aged
if (school_override == 1) {
scenario_no = 1;
}
//Set up scenarios
// School-aged vaccination
if (scenario_no == 1) {
// Add vaccine coverage
gamma[3] = gamma_school;
// Add vaccine effectiveness at reducing activation
alpha[3] = alpha_school[0];
alpha[4] = alpha_school[1];
alpha[5] = alpha_school[2];
alpha[6] = alpha_shool[3];
// Add vaccine effectiveness at reducing initial infection
chi[3] = chi_school;
chi[4] = chi_school;
chi[5] = chi_school;
chi[6] = chi_school;
// Neonatal vaccination
}else if (scenario_no == 2) {
// Add vaccine coverage
gamma[0] = gamma_neo;
// Add vaccine effectiveness at reducing activation
alpha[0] = alpha_neo[0];
alpha[1] = alpha_neo[1];
alpha[2] = alpha_neo[2];
alpha[3] = alpha_neo[3];
// Add vaccine effectiveness at reducing initial infection
chi[0] = chi_neo;
chi[1] = chi_neo;
chi[2] = chi_neo;
chi[3] = chi_neo;
}else if (scenario_no == 3) {
// No vaccination make no alterations
}
// Assign proportion of cases that are smear positive
double tau[nstageA];
for (i=0; i < 3; i++) {
tau[i] = tau_child;
}
for (i=3; i < 12; i++) {
tau[i] = tau_adult;
}
for (i=12; i < nstageA; i++) {
tau[i] = tau_older_adult;
}
// Population
double N = 0;
double N_age[nstageA];
double replacement_people;
for (i=0; i < nstageA; i++) {
N_age[i] = S[i] + E_1[i] + E_2[i] + E_R[i] + I[i] + S_v[i] + E_1v[i] + E_2v[i] + E_Rv[i] + I_v[i];
N += N_age[i];
}
//Dynamic Equations birth and ageing for first age group
DS[0] = (1 - gamma[0]) * omega - (theta[0] + mu[0]) * S[0];
DE_1[0] = - (theta[0] + mu[0]) * E_1[0];
DE_2[0] = - (theta[0] + mu[0]) * E_2[0];
DE_R[0] = - (theta[0] + mu[0]) * E_R[0];
DI[0] = - (theta[0] + mu[0])* I[0];
DS_v[0] = gamma[0] * omega - (theta[0] + mu[0]) * S_v[0];
DE_1v[0] = - (theta[0] + mu[0]) * E_1v[0];
DE_2v[0] = - (theta[0] + mu[0]) * E_2v[0];
DE_Rv[0] = - (theta[0] + mu[0]) * E_Rv[0];
DI_v[0] = - (theta[0] + mu[0]) * I_v[0];
//Account for burn in for first age group - replace difference between leaving and staying
if (burn_in_logical == 1) {
replacement_people = (theta[0] + mu[0]) * N_age[0] - omega;
DS[0] += (1-gamma[0]) * replacement_people;
DS_v[0] += gamma[0] * replacement_people;
}
//Dynamic Equations Ageing and Vaccination all other age groups
for (i = 1; i < nstageA; i++) {
DS[i] = (1 - gamma[i]) * theta[i - 1] * S[i - 1] - theta[i] * S[i] - mu[i] * S[i];
DE_1[i] = theta[i - 1] * E_1[i - 1] - theta[i] * E_1[i] - mu[i] * E_1[i];
DE_2[i] = theta[i - 1] * E_2[i - 1] - theta[i] * E_2[i] - mu[i] * E_2[i];
DE_R[i] = theta[i - 1] * E_R[i - 1] - theta[i] * E_R[i] - mu[i] * E_R[i];
DI[i] = theta[i - 1] * I[i - 1] - theta[i] * I[i] - mu[i] * I[i];
DS_v[i] = gamma[i] * theta[i - 1] * S[i - 1] + theta[i - 1] * S_v[i - 1] - theta[i] * S_v[i] - mu[i] * S_v[i];
DE_1v[i] = theta[i - 1] * E_1v[i - 1] - theta[i] * E_1v[i] - mu[i] * E_1v[i];
DE_2v[i] = theta[i - 1] * E_2v[i - 1] - theta[i] * E_2v[i] - mu[i] * E_2v[i];
DE_Rv[i] = theta[i - 1] * E_Rv[i - 1] - theta[i] * E_Rv[i] - mu[i] * E_Rv[i];
DI_v[i] = theta[i - 1] * I_v[i - 1] - theta[i] * I_v[i] - mu[i] * I_v[i];
//Account for burn in for first age group - replace difference between leaving and staying
if (burn_in_logical == 1) {
replacement_people = (theta[i] + mu[i]) * N_age[i] - (theta[i - 1]) * N_age[i - 1];
DS[i] += (1-gamma[i]) * replacement_people;
DS_v[i] += gamma[i] * replacement_people;
}
}
// Dynamic Equations Disease
for (i=0; i < nstageA; i++) {
// Determine age group for TB mortality rate
if (i < 3) {
mu_t = mu_t_child;
}else if (3 >= i < 12) {
mu_t = mu_t_adult;
}else{
mu_t = mu_t_older_adult
}
//Determine age group for latent parameters
if (i < 1) {
epsilon_1 = epsilon_1_under_5;
kappa = kappa_under_5;
epsilon_2 = epsilon_2_under_5;
}else if (1 >= i < 3) {
epsilon_1 = epsilon_1_under_15;
kappa = kappa_under_15;
epsilon_2 = epsilon_2_under_15;
}else{
epsilon_1 = epsilon_1_adult;
kappa = kappa_adult;
epsilon_2 = epsilon_2_adult;
}
// Force of infection
foi = 0;
for (j = 0; j < nstageA; j++) {
foi += lambda * tau[j] * C[i * nstageA + j] * (I[j] + I_v[j] + M * iota[j] / nu) / N;
}
// Unvaccinated transitions
trans[0] = foi * S[i];
trans[1] = epsilon_1 * E_1[i];
trans[2] = kappa * E_1[i];
trans[3] = epsilon_2 * E_2[i];
trans[4] = foi * E_2[i];
trans[5] = (1 - delta[i]) * epsilon_1 * E_R[i];
trans[6] = kappa * E_R[i];
trans[7] = nu * I[i];
trans[16] = mu_t * I[i];
// Vaccinated transitions
trans[8] = (1 - chi[i]) * foi * S_v[i];
trans[9] = (1 - alpha[i]) * epsilon_1 * E_1v[i];
trans[10] = kappa * E_1v[i];
trans[11] = (1- alpha[i]) * epsilon_2 * E_2v[i];
trans[12] = foi * E_2v[i];
trans[13] = (1 - delta[i]) * epsilon_1 * E_Rv[i];
trans[14] = kappa * E_Rv[i];
trans[15] = nu * I_v[i];
trans[17] = mu_t * I_v[i];
//Update DS
DS[i] += - trans[0] + trans[7];
DE_1[i] += trans[0] - trans[1] - trans[2];
DE_2[i] += trans[2] + trans[6] - trans[3] - trans[4];
DE_R[i] += trans[4] - trans[5] - trans[6];
DI[i] += trans[1] + trans[3] + trans[5] - trans[7] - trans[16];
DS_v[i] += - trans[8] + trans[15];
DE_1v[i] += trans[8] - trans[9] - trans[10];
DE_2v[i] += trans[10] + trans[14] - trans[11] - trans[12];
DE_Rv[i] += trans[12] - trans[13] - trans[14];
DI_v[i] += trans[9] + trans[11] + trans[13] - trans[15] - trans[17];
// Transition rates
DN_E_1I[i] = trans[1];
DN_E_2I[i] = trans[3];
DN_E_RI[i] = trans[5];
DN_E_1I_v[i] = trans[9];
DN_E_2I_v[i] = trans[11];
DN_E_RI_v[i] = trans[13];
DDeaths[i] = trans[16] + trans[17];
//Account for burn in - replace those who die from TB
if (burn_in_logical == 1) {
DS[i] += (1-gamma[i]) * (trans[16] + trans[17]);
DS_v[i] += gamma[i] * (trans[16] + trans[17]);
}
}
")
## Generate observations
rObs <- Csnippet("
const double *N_E_1I = &N_E_1I1; const double *N_E_2I = &N_E_2I1;
const double *N_E_RI = &N_E_RI1;
const double *N_E_1I_v = &N_E_1I_v1; const double *N_E_2I_v = &N_E_2I_v1;
const double *N_E_RI_v = &N_E_RI_v1;
double *incidence = &incidence1;
double D_inc;
int i;
for(i=0;i < nstageA; i++) {
D_inc = N_E_1I[i] + N_E_2I[i] + N_E_RI[i] + N_E_1I_v[i] + N_E_2I_v[i] + N_E_RI_v[i];
if (k > 0) {
incidence[i] = rnbinom_mu(1.0 / k, rho * D_inc > 0 ? rho * D_inc : 0);
} else {
incidence[i] = rpois(rho * D_inc > 0 ? rho * D_inc : 0);
}
}
")
dObs <- Csnippet("
int i;
double f = 0;
double D_inc;
const double *N_E_1I = &N_E_1I1; const double *N_E_2I = &N_E_2I1;
const double *N_E_RI = &N_E_RI1;
const double *N_E_1I_v = &N_E_1I_v1; const double *N_E_2I_v = &N_E_2I_v1;
const double *N_E_RI_v = &N_E_RI_v1;
const double *incidence = &incidence1;
if (ISNA(incidence1)) {
lik = (give_log) ? 0 : 1;
} else {
for(i=0;i < nstageA; i++) {
D_inc= N_E_1I[i] + N_E_2I[i] + N_E_RI[i] + N_E_1I_v[i] + N_E_2I_v[i] + N_E_RI_v[i];
if (k > 0.0) {
f += dnbinom_mu(nearbyint(incidence[i]),1.0/k,rho * D_inc > 0 ? rho * D_inc : 0,1);
} else {
f += dpois(nearbyint(incidence[i]),rho * D_inc > 0 ? rho * D_inc : 0,1);
}
}
lik = (give_log) ? f : exp(f);
}")
init <- Csnippet("
double *S = &S1; double *E_1 = &E_11;
double *E_2 = &E_21; double *I = &I1;
double *E_R = &E_R1;
double *S_v = &S_v1; double *E_1v = &E_1v1;
double *E_2v = &E_2v1; double *I_v = &I_v1;
double *E_Rv = &E_Rv1;
const double *age_weight = &age_weight1;
double *gamma = &gamma1;
int j;
// Assume school age vaccination is in place when the model is initialised (for those at school age i.e 15-19 years old)
gamma[3] = gamma_school;
for (j = 0; j < nstageA; j++) {
S[j] = nearbyint(N_init * (1 - gamma[j]) * age_weight[j]) - nearbyint((init_I + init_E_1 + init_E_2 + init_E_R) * (1 - gamma[j]) * age_weight[j]);
S_v[j] = nearbyint(N_init * gamma[j] * age_weight[j]) - nearbyint((init_I + init_E_1 + init_E_2 + init_E_R) * gamma[j] * age_weight[j]);
E_1[j] = nearbyint(init_E_1 * (1 - gamma[j]) * age_weight[j]);
E_2[j] = nearbyint(init_E_2 * (1 - gamma[j]) * age_weight[j]);
E_R[j] = nearbyint(init_E_R * (1 - gamma[j]) * age_weight[j]);
E_1v[j] = nearbyint(init_E_1 * gamma[j] * age_weight[j]);
E_2v[j] = nearbyint(init_E_2 * gamma[j] * age_weight[j]);
E_Rv[j] = nearbyint(init_E_R * gamma[j] * age_weight[j]);
I[j] = nearbyint(init_I * ( 1 - gamma[j]) * age_weight[j]);
I_v[j] = nearbyint(init_I * gamma[j] * age_weight[j]);
}
")
## Specify Priors
##sample from prior
rprior <- Csnippet( '
int j;
// Disease parameters
lambda = runif(0, 1);
nu = 1 /(rgamma(1.056, 3.26));
rho = 0;
// Latent parameters
// Transition to active from high risk
epsilon_1_under_5 = 365.25 * rnorm(6.95e-3, 1.30e-3);
epsilon_1_under_15 = 365.25 * rnorm(2.8e-3, 5.61e-4);
epsilon_1_adult = 365.25 * rnorm(3.35e-4, 8.93e-5);
// Transition from high to low risk latent
kappa_under_5 = 365.25 * rnorm(1.33e-2, 2.42e-3);
kappa_under_15 = 365.25 * rnorm(1.20e-2, 2.07e-3);
kappa_adult = 365.25 * rnorm(7.25e-3, 1.91e-3);
// Transition from low risk latent to active disease
epsilon_2_under_5 = 365.25 * rnorm(8.00e-6 , 4.08e-6);
epsilon_2_under_15 = 365.25 * rnorm(9.84e-6, 4.67e-6);
epsilon_2_adult = 365.25 * rnorm(5.95e-6 , 2.07e-6);
//TB mortality
mu_t_child = 0;
mu_t_adult = 0;
mu_t_older_adult = 0;
//Reduction in activation due to previous infection
delta_all_age = rnorm(0.784, 0.0286);
//Proportion of cases that are smear positive
tau_child = rnorm(3.02e-1, 1.89e-2);
tau_adult = rnorm(6.52e-1, 5.18e-3);
tau_older_adult = rnorm(5.36e-1, 8.45e-3);
// Protection from infection in the BCG vaccinated
chi_school = rnorm(0.19, 0.1);
chi_neo = rnorm(0.19, 0.1);
//Vaccine coverage estimates
gamma_school = rnorm(0.75, 0.0255);
gamma_neo = rnorm(0.75, 0.0255);
//Vaccine effectiveness estimates
// Vaccine effectiveness in school aged
alpha_school[0] = rnorm(0.784, 0.0286);
alpha_school[1] = rnorm(0.784, 0.0286);
alpha_school[2] = rnorm(0.784, 0.0286);
alpha_school[3] = rnorm(0.784, 0.0286);
// Vaccine effectiveness at birth
alpha_neo[0] = rnorm(0.784, 0.0286);
alpha_neo[1] = rnorm(0.784, 0.0286);
alpha_neo[2] = rnorm(0.784, 0.0286);
alpha_neo[3] = rnorm(0.784, 0.0286);
')
## Parameter transforms
toEst <- Csnippet("
int i; int j;
const double *theta = &theta1;
const double *age_weight = &age_weight1;
const double *gamma = &gamma1;
const double *delta = &delta1;
const double *alpha = &alpha1;
const double *alpha_school = &alpha_school1;
const double *alpha_neo = &alpha_neo1;
const double *C = &C1;
double *Ttheta = &Ttheta1;
double *Tage_weight = &Tage_weight1;
double *Tgamma = &Tgamma1;
double *Tdelta = &Tdelta1;
double *Talpha = &Talpha1;
double *Talpha_school = &Talpha_school1;
double *Talpha_neo = &Talpha_neo1;
double *TC = &TC1;
// Single value parameters
Tlambda = log(lambda);
Tnu = log(nu);
Tdelta_all_age = logit(delta_all_age);
Tgamma_school = logit(gamma_school);
Tgamma_neo = logit(gamma_neo);
Tchi_school = logit(chi_school);
Tchi_neo= logit(chi_neo);
Tmu_t_child = log(mu_t_child);
Tmu_t_adult = log(mu_t_adult);
Tmu_t_older_adult = log(mu_t_older_adult);
Trho = logit(rho);
Tk = log(k);
TM = log(M);
Ttau_child = logit(tau_child);
Ttau_adult = logit(tau_adult);
Ttau_older_adult = logit(tau_older_adult);
//Latent parameters
Tepsilon_1_under_5 = log(epsilon_1_under_5);
Tepsilon_1_under_15 = log(epsilon_1_under_15);
Tepsilon_1_adult = log(epsilon_1_adult);
Tkappa_under_5 = log(kappa_under_5);
Tkappa_under_15 = log(kappa_under_15);
Tkappa_adult = log(kappa_adult);
Tepsilon_2_under_5 = log(epsilon_2_under_5);
Tepsilon_2_under_15 = log(epsilon_2_under_15);
Tepsilon_2_adult = log(epsilon_2_adult);
// Initial populations
to_log_barycentric(&TN_init, &N_init, 1);
to_log_barycentric(&Tinit_I, &init_I, 1);
to_log_barycentric(&Tinit_E_1, &init_E_1, 1);
to_log_barycentric(&Tinit_E_2, &init_E_2, 1);
to_log_barycentric(&Tinit_E_R, &init_E_R, 1);
// Vectorised parameters
for (i = 0; i < nstageA; i++) {
Ttheta[i] = log(theta[i]);
Tage_weight[i] = logit(age_weight[i]);
Tgamma[i] = logit(gamma[i]);
Talpha[i] = logit(alpha[i]);
Tdelta[i] = logit(delta[i]);
for (j = 0; j < nstageA; j++) {
TC[i * nstageA + j] = log(C[i * nstageA + j]);
}
}
// Vaccine effectiveness
for (i = 0; i < 3; i++) {
Talpha_school[i] = logit(alpha_school[i]);
Talpha_neo[i] = logit(alpha_neo[i]);
}
")
fromEst <- Csnippet("
int i; int j;
const double *theta = &theta1;
const double *age_weight = &age_weight1;
const double *gamma = &gamma1;
const double *alpha = &alpha1;
const double *alpha_school = &alpha_school1;
const double *alpha_neo = &alpha_neo1;
const double *delta = &delta1;
const double *C = &C1;
double *Ttheta = &Ttheta1;
double *Tage_weight = &Tage_weight1;
double *Talpha = &Talpha1;
double *Talpha_school = &Talpha_school1;
double *Talpha_neo = &Talpha_neo1;
double *Tdelta = &Tdelta1;
double *TC = &TC1;
// Single value parameters
Tlambda = exp(lambda);
Tnu = exp(nu);
Trho = expit(rho);
Tk = exp(k);
Tgamma_school = expit(gamma_school);
Tgamma_neo = expit(gamma_neo);
Tdelta_all_age = expit(delta_all_age);
Tchi_school = expit(chi_school);
Tchi_neo = expit(chi_neo);
mu_t_child = exp(Tmu_t_child);
mu_t_adult = exp(Tmu_t_adult);
mu_t_older_adult = exp(Tmu_t_older_adult);
M = exp(TM);
Ttau_child = expit(tau_child);
Ttau_adult = expit(tau_adult);
Ttau_older_adult = expit(tau_older_adult);
//Latent parameters
epsilon_1_under_5 = exp(Tepsilon_1_under_5);
epsilon_1_under_15 = exp(Tepsilon_1_under_15);
epsilon_1_adult = exp(Tepsilon_1_adult);
kappa_under_5 = exp(Tkappa_under_5);
kappa_under_15 = exp(Tkappa_under_15);
kappa_adult = exp(Tkappa_adult);
epsilon_2_under_5 = exp(Tepsilon_2_under_5);
epsilon_2_under_15 = exp(Tepsilon_2_under_15);
epsilon_2_adult = exp(Tepsilon_2_adult);
// Initial populations
from_log_barycentric(&TN_init, &N_init, 1);
from_log_barycentric(&Tinit_I, &init_I, 1);
from_log_barycentric(&Tinit_E_1, &init_E_1, 1);
from_log_barycentric(&Tinit_E_2, &init_E_2, 1);
from_log_barycentric(&Tinit_E_R, &init_E_R, 1);
// Vectorised parameters
for (i = 0; i < nstageA; i++) {
Ttheta[i] = exp(theta[i]);
Tage_weight[i] = expit(age_weight[i]);
Tgamma[i] = expit(gamma[i]);
Tdelta[i] = expit(delta[i]);
Talpha[i] = expit(alpha[i]);
for (j = 0; j < nstageA; j++) {
TC[i * nstageA + j] = exp(C[i * nstageA + j]);
}
}
// Vaccine effectiveness
for (i = 0; i < 3; i++) {
Talpha_neo[i] = expit(alpha_neo[i]);
Talpha_school[i] = expit(alpha_school[i]);
}
")
## Define global parmeters
globs <- paste0("static int nstageA = ", nstageA, ";")
## Define statenames
statenames <- c("S", "E_1", "E_2", "E_R", "I", "S_v", "E_1v", "E_2v", "E_Rv", "I_v",
"N_E_1I", "N_E_2I", "N_E_RI", "N_E_1I_v", "N_E_2I_v", "N_E_RI_v", "Deaths") %>%
add_pointer_struct(nstageA)
zeronames <- c("N_E_1I", "N_E_2I", "N_E_RI", "N_E_1I_v", "N_E_2I_v", "N_E_RI_v", "Deaths") %>%
add_pointer_struct(nstageA)
##Define vectorised params
params <- c(lambda = 1, tau_child = 1, tau_adult = 1, tau_older_adult = 1,nu = 1,
epsilon_1_under_5 = 1, epsilon_1_under_15 = 1, epsilon_1_adult = 1,
kappa_under_5 = 1, kappa_under_15 = 1, kappa_adult = 1,
epsilon_2_under_5 = 1, epsilon_2_under_15 = 1, epsilon_2_adult = 1,
N_init = 1, init_I = 1, init_E_1 = 1, init_E_2 = 1, init_E_R = 1,
k = 1, rho = 1, theta = 1:nstageA,
age_weight = 1:nstageA, gamma = 1:nstageA, gamma_school = 1, gamma_neo = 1,
alpha = 1:nstageA, alpha_school = 1:4, alpha_neo = 1:4, chi = 1:nstageA, chi_school = 1, chi_neo = 1,
delta = 1:nstageA, delta_all_age = 1, C = 1:nstageA^2 , mu_t_child = 1, mu_t_adult = 1,
scenario_no = 1, M = 1) %>%
names
## Define model
model <- df %>%
pomp(
skeleton = vectorfield(rskel),
#rprocess = euler.sim(step.fun = rsim,
# delta.t = timestep),
rmeasure = rObs,
dmeasure = dObs,
covar = covars,
tcovar = tcovar,
initializer = init,
rprior = rprior,
toEstimationScale = toEst,
fromEstimationScale = fromEst,
times = time,
t0 = t0,
globals = globs,
zeronames = zeronames,
paramnames = params,
statenames = statenames,
obsnames = obsnames,
...
)
return(model)
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment