I an running a stan programme in r but getting error in stanc(file = file, model_code = model_code, model_name = model_name, : parser failed badly
data {
int <lower=1> M; // number of states
int <lower=1> P; // number of covariates
int <lower=1> N0; // number of days for which to impute infections
int <lower=1> N[M]; // days of observed data for country m. each entry must be <= N2
int <lower=1> N2; // days of observed data + # of days to forecast
int cases[N2,M]; // reported cases
int deaths[N2, M]; // reported deaths -- the rows with i > N contain -1 and should be ignored
matrix[N2, M] f; // h * s
matrix[N2, P] X[M]; // features matrix
int EpidemicStart[M];
real pop[M];
real SI[N2]; // fixed pre-calculated SI using emprical data from Neil
}
transformed data {
vector[N2] SI_rev; // SI in reverse order
vector[N2] f_rev[M]; // f in reversed order
for(i in 1:N2)
SI_rev[i] = SI[N2-i+1];
for(m in 1:M){
for(i in 1:N2) {
f_rev[m, i] = f[N2-i+1,m];
}
}
}
parameters {
real<lower=0> mu[M]; // intercept for Rt
real<lower=0> alpha_hier[P]; // sudo parameter for the hier term for alpha
real<lower=0> kappa;
real<lower=0> y[M];
real<lower=0> phi;
real<lower=0> tau;
real <lower=0> ifr_noise[M];
}
transformed parameters {
vector[P] alpha;
matrix[N2, M] prediction = rep_matrix(0,N2,M);
matrix[N2, M] E_deaths = rep_matrix(0,N2,M);
matrix[N2, M] Rt = rep_matrix(0,N2,M);
matrix[N2, M] Rt_adj = Rt;
{
matrix[N2,M] cumm_sum = rep_matrix(0,N2,M);
for(i in 1:P){
alpha[i] = alpha_hier[i] - ( log(1.05) / 6.0 );
}
for (m in 1:M){
prediction[1:N0,m] = rep_vector(y[m],N0); // learn the number of cases in the first N0 days
cumm_sum[2:N0,m] = cumulative_sum(prediction[2:N0,m]);
Rt[,m] = mu[m] * exp(-X[m] * alpha);
Rt_adj[1:N0,m] = Rt[1:N0,m];
for (i in (N0+1):N2) {
real convolution = dot_product(sub_col(prediction, 1, m, i-1), tail(SI_rev, i-1));
cumm_sum[i,m] = cumm_sum[i-1,m] + prediction[i-1,m];
Rt_adj[i,m] = ((pop[m]-cumm_sum[i,m]) / pop[m]) * Rt[i,m];
prediction[i, m] = Rt_adj[i,m] * convolution;
}
E_deaths[1, m]= 1e-15 * prediction[1,m];
for (i in 2:N2){
E_deaths[i,m] = ifr_noise[m] * dot_product(sub_col(prediction, 1, m, i-1), tail(f_rev[m], i-1));
}
}
}
}
model {
tau ~ exponential(0.03);
for (m in 1:M){
y[m] ~ exponential(1/tau);
}
phi ~ normal(0,5);
kappa ~ normal(0,0.5);
mu ~ normal(3.28, kappa); // citation: https://academic.oup.com/jtm/article/27/2/taaa021/5735319
alpha_hier ~ gamma(.1667,1);
ifr_noise ~ normal(1,0.1);
for(m in 1:M){
deaths[EpidemicStart[m]:N[m], m] ~ neg_binomial_2(E_deaths[EpidemicStart[m]:N[m], m], phi);
}
}
generated quantities {
matrix[N2, M] prediction0 = rep_matrix(0,N2,M);
matrix[N2, M] E_deaths0 = rep_matrix(0,N2,M);
{
matrix[N2,M] cumm_sum0 = rep_matrix(0,N2,M);
for (m in 1:M){
for (i in 2:N0){
cumm_sum0[i,m] = cumm_sum0[i-1,m] + y[m];
}
prediction0[1:N0,m] = rep_vector(y[m],N0);
for (i in (N0+1):N2) {
real convolution0 = 0;
for (j in 1:(i-1)) {
convolution0 += prediction0[j, m] * SI[i-j];
}
cumm_sum0[i,m] = cumm_sum0[i-1,m] + prediction0[i-1,m];
prediction0[i, m] = ((pop[m]-cumm_sum0[i,m]) / pop[m]) * mu[m] * convolution0;
}
E_deaths0[1, m]= 1e-15 * prediction0[1,m];
for (i in 2:N2){
for (j in 1:(i-1)) {
E_deaths[i,m] += prediction0[j, m] * f[i-j,m] * ifr_noise[m];
}
}
}
}
}
[base.stan|attachment](upload://mlB6qzA1GKRB2odpEiewSlD0e6T.stan) (3.6 KB)
If you need to use math formula, use Latex syntax:
Y \sim N(\mu, \sigma)
Don’t forget to attach tags (top right of this form) for application area/class of models or other general subject areas your topic touches on.