Search code examples
rbayesianrjags

Regarding a warning message in JAGS for R

I am experimenting with rjags r package while learning Bayesian. I fitted following multivariate model using tobacco data set in rrr package.

enter image description here

where \gamma is a random intercept which measures the correlation between the outcomes. Also I purposefully make some values of the response as missing to explore the results with a unbalanced structure of the data. At the same time, I added AID variable which denotes the subject ids.

I fitted this model using rjags as follows.

library(rrr)
require(dplyr)
library(rjags)
data("tobacco")
tobacco <- as_data_frame(tobacco)

tobacco$AID=seq(1:25)
tobacco[4,1]=NA
tobacco[14,1]=NA
tobacco[8,1]=NA
tobacco[6,2]=NA
tobacco[1,2]=NA
tobacco[19,2]=NA
tobacco[21,2]=NA


N1=length(tobacco$Y1.BurnRate)
bayes_model_mul="model {
for(i in 1:N1){
Y1.BurnRate[i]~dnorm(mu1[i],tau1)
Y2.PercentSugar[i]~dnorm(mu2[i],tau2)
mu1[i]=beta1[1] + beta1[2]*X2.PercentChlorine[i] + beta1[3]*X3.PercentPotassium[i] + gamma[AID[i]]
mu2[i]=beta2[1] + beta2[3]*X2.PercentChlorine[i] + beta2[2]*X1.PercentNitrogen[i]+
beta2[4]*X3.PercentPotassium[i]+gamma[AID[i]]
gamma[i] ~ dmnorm(0,tau_u)

}
for (l in 1:3) { beta1[l] ~dnorm(0, 0.001) }
for (l in 1:4) { beta2[l] ~dnorm(0, 0.001) }
tau1 ~ dgamma(.01,.01)
sigma_tau1 = 1/tau1 

tau2 ~ dgamma(.01,.01)
sigma_tau2 = 1/tau2

tau_u ~ dgamma(.01,.01)
sigma_tau_u = 1/tau_u 
}"


model3 <- jags.model(textConnection(bayes_model_mul), 
                     data = list(Y1.BurnRate=tobacco$Y1.BurnRate,
                                 Y2.PercentSugar=tobacco$Y2.PercentSugar
                                 ,X1.PercentNitrogen=tobacco$X1.PercentNitrogen,
                                 N1=N1,X2.PercentChlorine=tobacco$X2.PercentChlorine,
                                 X3.PercentPotassium=tobacco$X3.PercentPotassium,AID=tobacco$AID),
                     n.chains=1)
params <- c('beta1','sigma_tau1','sigma_tau2','beta2','sigma_tau_u','gamma')
samps.1 <- coda.samples(model3, params, n.iter = 10000)
burn.in=1000
summary.model.1=summary(window(samps.1, start = burn.in))

I didnt get any error. But I got following warning message.

Warning message:
In FUN(X[[i]], ...) : start value not changed

Can anyone help me to figure what is this error message is about ?

Thank you.

Solution

  • The window function is used to select a range of samples. It takes start and end parameters. The value that start takes should be greater than the the sum of the number of iterations passed to n.adapt in jags.models and n.iter in update, if explicitly setting a number of burnin samples.

    You can have a look at the object retuned by jags.model to see how the number of iterations changes.

    # Define the model
model3 <- jags.model(textConnection(bayes_model_mul), n.adapt=1234,
                     data = list(Y1.BurnRate=tobacco$Y1.BurnRate,Y2.PercentSugar=tobacco$Y2.PercentSugar,X1.PercentNitrogen=tobacco$X1.PercentNitrogen,N1=N1,X2.PercentChlorine=tobacco$X2.PercentChlorine,X3.PercentPotassium=tobacco$X3.PercentPotassium,AID=tobacco$AID),
                     n.chains=1)
# Look at the number of iterations (== n.adapt)
model3$iter()
#1234

# Add some burnin iterations: iterations are updated to n.adapt + burnin
update(model3, n.iter=1111)
model3$iter()
# 2345

# Draw more samples -> iterations updated to n.adapt + burnin + n.iter
samps.1 <- coda.samples(model3, n.iter = 10000, variable.names=c('beta1','sigma_tau1','sigma_tau2','beta2','sigma_tau_u','gamma'))
model3$iter()
# 12345

    When we use window the start value is expected to be greater than n.adapt + burnin as these iterations have been discarded. If it is not then a start value is extracted from the coda.samples object (see mcpar(samps.1[[1]])) and your manual start value is not used.

    So using w=window(samps.1, start=2345) gives the warning that you see

    Warning message:
    In FUN(X[[i]], ...) : start value not changed

    But the following is okay as the start is greater than n.adapt + burnin

    w=window(samps.1, start=2346)

    Remember these are just warnings and in this case are not that important. But if you rerun coda.samples without rerunning jags.model this will again update the iterations which makes it a bit harder to track what values to pass to window (and can then throw errors)