0
votes

I am new to bayesian statistical analysis with JAGS/WinBUGS/STAN languages. I am trying to fit an extended SIR model, it seems the error is located in line 2 since I always get the warning: "Error parsing model file: syntax error on line 2 near "for""

But I can't get my head around this error, the data: y is in matrix form (T x K), N vector size K, L vector of size K.

bayesian.model <- "
                  for(r in 1:K){
                    for(t in 2:T){

                         S[t,r] <- S[t-1,r] - y[t,r] + b*R[t-1,r]
                         E[t,r] <- E[t-1,r] + y[t,r] - c*E[t-1,r]
                         I[t,r] <- I[t-1,r] - a*I[t-1,r] + c*E[t-1,r]
                         R[t,r] <- R[t-1,r] + a*I[t-1,r] - b*R[t-1,r]
 
                         y[t,r] ~ dbin(theta[t,r], S[t-1,r])
                         logit(theta[t,r]) <-  alpha*log(I[t-1,r]) +  rzero + beta*L[r]
                     }
                   }

                       # PARAMETERS OF THE MODEL
                              a <- 1/17.5
                              b <- 1/m
                              m ~ dnorm(300, 20)
                              c <- 1/11.17
                              for(r in 1:K){
                                    S[1,r] <- N[r] - 100 - 0.5*100 
                                    E[1,r] <- 0.5*100
                                    I[1,r] <- 100
                                    R[1,r] <- 0
                                    y[1,r] ~ dbin(theta[1,r], S[1,r])
                                    logit(theta[1,r]) <- alpha*log(I[1,r]+100) + rzero + 
                                                                    beta*L[r]                                                 
                              }
                                  rzero ~ dbeta(0.1, 0,3)
                                  beta ~ dbeta(0.1, 0.3)
                                  alpha <- 0.8
"
K <- ncol(y)
T <- nrow(y)
jags.data <- list(y = data, N = data$Population, L = data$Location, K = K, T = T)

parameters <- c("theta")
model <- jags(data = jags.data, parameters.to.save= parameters,
              n.iter = 10000,n.burnin = 0, n.thin = 1,
              n.chains = 2, model = textConnection(bayesian.model))
1
do you need model { ...} in front of your model-string? see here: stats.stackexchange.com/questions/85690/…holzben
Thank you!!! It was exactly that (and a comma in the beta distribution haha). Thanks :)Ajam
how bout formalizing your answer so it can be accepted @holzbenGmichael

1 Answers

0
votes

Like suggested in the comments adding { ...} and replacing rzero ~ dbeta(0.1, 0,3) with rzero ~ dbeta(0.1, 0.3) should work:

bayesian.model <- "
{
                  for(r in 1:K){
                    for(t in 2:T){

                         S[t,r] <- S[t-1,r] - y[t,r] + b*R[t-1,r]
                         E[t,r] <- E[t-1,r] + y[t,r] - c*E[t-1,r]
                         I[t,r] <- I[t-1,r] - a*I[t-1,r] + c*E[t-1,r]
                         R[t,r] <- R[t-1,r] + a*I[t-1,r] - b*R[t-1,r]
 
                         y[t,r] ~ dbin(theta[t,r], S[t-1,r])
                         logit(theta[t,r]) <-  alpha*log(I[t-1,r]) +  rzero + beta*L[r]
                     }
                   }

                       # PARAMETERS OF THE MODEL
                              a <- 1/17.5
                              b <- 1/m
                              m ~ dnorm(300, 20)
                              c <- 1/11.17
                              for(r in 1:K){
                                    S[1,r] <- N[r] - 100 - 0.5*100 
                                    E[1,r] <- 0.5*100
                                    I[1,r] <- 100
                                    R[1,r] <- 0
                                    y[1,r] ~ dbin(theta[1,r], S[1,r])
                                    logit(theta[1,r]) <- alpha*log(I[1,r]+100) + rzero + 
                                                                    beta*L[r]                                                 
                              }
                                  rzero ~ dbeta(0.1, 0.3)
                                  beta ~ dbeta(0.1, 0.3)
}                                  alpha <- 0.8 
"