run.jags running 20000 iterations regardless of samples statement

[Matthew McKinney]

Hi.
My version of run.jags seems to run 20000 iterations regardsless of what I specify in the "samples" statement.

> cjagsout<-run.jags(model="my model",
+                    monitor=parms2,data=forjags,n.chains=2,inits=j.inits,burnin=10000,adapt=1000,
+                    sample=10000,thin=2,jags.refresh=60,method="rjags",modules = "glm on",summarise = TRUE)
Compiling rjags model...
Calling the simulation using the rjags method...
Adapting the model for 1000 iterations...
  |++++++++++++++++++++++++++++++++++++++++++++++++++| 100%
Burning in the model for 10000 iterations...
  |**************************************************| 100%
Running the model for 20000 iterations...
  |**************************************************| 100%
Simulation complete
Calculating summary statistics...
Calculating the Gelman-Rubin statistic for 10 variables....
Finished running the simulation
> 

Any idea what I'm doing wrong?
Many thanks,
Matt

[Matt Denwood]

You have specified 10000 samples with a thin of 2, which means that JAGS is only storing every second iteration, so you need 20000 iterations to get 10000 samples, which is what happens. But if you look at the number of iterations stored (per chain) it should be 10000. As far as I can see this is the case e.g.:

library('coda')
library('runjags')

X <- 1:100
Y <- rnorm(length(X), 2*X + 10, 1)

model <- "model { 
for(i in 1 : length(Y)){ 
    Y[i] ~ dnorm(true.y[i], precision);
    true.y[i] <- (m * X[i]) + c
} 
m ~ dunif(-1000,1000)
c ~ dunif(-1000,1000) 
precision ~ dexp(1)
#data# Y, X
#monitor# m, c, precision
}"

(nc <- sample(1:4,1))
(sa <- sample(1:5,1)*1000)
(th <- sample(1:5,1))
results <- run.jags(model=model, n.chains=nc, sample=sa, thin=th, method="rjags")

# These should match:
cat('... ', sa*nc, ' samples (thin = ', th, '; chains = ', nc, '; ...\n', sep='')
results

# Convert to MCMC object:
mcmc <- as.mcmc.list(results)

# Reports the expected samples, chains and thin:
summary(mcmc)

# Total samples is as expected:
stopifnot(niter(mcmc) == sa)

If you find a combination of parameters that doesn't work, could you post a reproducible example?

Matt

[Matthew McKinney]

Thanks, Mr. Denwood. I just misunderstood the way the function works. I thought it would sample 10000, thinning by 2 along the way, to keep a remaining 5000.

Matt