Reputation: 4184
I'm running a JAGS (JAGS 4.3.2) model in R using the rjags
package (version 4.15), but I'm encountering an issue with the n.adapt
parameter. Specifically, when I set n.adapt
to any number like 1000, it seems like the adaptation phase is not occurring as expected because jags_model$iter() returns zero, and no error message is printed. The final results should show 1001:x, but they do not, which is another confirmation that the adaptation phase was not run. BTW. Even if I set a string for n.adapt, it will not print an error and will run successfully with 0 adaptation iterations.
The final model runs properly, and the manual pre-run works fine (update(..., n.iter = )
). Thus, the model itself seems healthy.
I look specifically for why n.adapt
in jags.model
is not working and why there is no JAGS error message.
library(rjags)
library(coda)
# Generate synthetic data
set.seed(123)
# Write the JAGS model to a text file
model_string <- "
model {
# Likelihood for arm-based data
for (i in study_ids) {
for (k in 1:num_arms[i]) {
theta[i, k] <- mu[i] + delta[i, k]
measurement[i, k] ~ dnorm(theta[i, k], prec[i, k])
prec[i, k] <- 1 / (error[i, k] * error[i, k])
deviance[i, k] <- pow(measurement[i, k] - theta[i, k], 2) * prec[i, k]
}
}
# Fixed effect model
for (i in study_ids) {
delta[i, 1] <- 0
for (k in 2:num_arms[i]) {
delta[i, k] <- d[treatment[i, 1], treatment[i, k]]
}
}
# Relative effect matrix
d[1, 1] <- 0
d[1, 2] <- d12
d[1, 3] <- d13
d[1, 4] <- d14
d[1, 5] <- d15
d[1, 6] <- d16
for (i in 2:num_treatments) {
for (j in 1:num_treatments) {
d[i, j] <- d[1, j] - d[1, i]
}
}
# Priors
prior_prec <- pow(re_prior_sd, -2)
for (i in study_ids) {
mu[i] ~ dnorm(0, prior_prec)
}
# Effect parameter priors
d12 ~ dnorm(0, prior_prec)
d13 ~ dnorm(0, prior_prec)
d14 ~ dnorm(0, prior_prec)
d15 ~ dnorm(0, prior_prec)
d16 ~ dnorm(0, prior_prec)
}
"
writeLines(model_string, con = "model.txt")
# Prepare data for JAGS
jags_data <- list(
treatment = matrix(sample(1:6, 96, replace = TRUE), nrow = 24, ncol = 4),
measurement = matrix(rnorm(96, -1.5, 1), nrow = 24, ncol = 4),
error = matrix(runif(96, 0.1, 1), nrow = 24, ncol = 4),
num_arms = rep(2:3, length.out = 24),
num_treatments = 6,
re_prior_sd = 1,
study_ids = 1:24
)
# Set up initial values (optional but recommended)
inits <- list(
list(
d12 = rnorm(1, 0, 1), d13 = rnorm(1, 0, 1), d14 = rnorm(1, 0, 1), d15 = rnorm(1, 0, 1), d16 = rnorm(1, 0, 1),
mu = rnorm(24, -1.5, 1)
),
list(
d12 = rnorm(1, 0, 1), d13 = rnorm(1, 0, 1), d14 = rnorm(1, 0, 1), d15 = rnorm(1, 0, 1), d16 = rnorm(1, 0, 1),
mu = rnorm(24, -1.5, 1)
),
list(
d12 = rnorm(1, 0, 1), d13 = rnorm(1, 0, 1), d14 = rnorm(1, 0, 1), d15 = rnorm(1, 0, 1), d16 = rnorm(1, 0, 1),
mu = rnorm(24, -1.5, 1)
),
list(
d12 = rnorm(1, 0, 1), d13 = rnorm(1, 0, 1), d14 = rnorm(1, 0, 1), d15 = rnorm(1, 0, 1), d16 = rnorm(1, 0, 1),
mu = rnorm(24, -1.5, 1)
)
)
# Define parameters to monitor
params <- c("d12", "d13", "d14", "d15", "d16", "mu", "theta", "deviance")
# Create and compile the JAGS model
jags_model <- jags.model("model.txt", data = jags_data, inits = inits, n.chains = 4, n.adapt = 1000)
# Check the iteration number
print(jags_model$iter())
# We can run the adapt manually but I am looking for why the n.adapt is not working
# Run properly
# update(jags_model, n.iter = 1000)
# samples <- coda.samples(jags_model, variable.names = params, n.iter = 5000)
Upvotes: 0
Views: 187
Reputation: 26
Expanding on the responses by @DaveArmstrong and @polkas, it is important to recognize that burn-in and adaptation refer to two different things and are implemented independently in JAGS/rjags
.
Here is an explanation of the differences between burn-in and adaptation: r2jags - How to interpret some syntax (n.adapt, update..) in jags? - Stack Overflow. However, this explanation does not explain when adaptation is run or not.
Burn-in is implemented in rjags using the update
function. If specified, the model will always run with the burn-in iterations.
Adaptation is implemented using the n.adapt
argument in the jags.model
function. However, the nature of the model itself determines whether the adaptation process will run. A deeper dive into JAGS C++ code shows that adaptation is based on an evaluation of isAdaptive
, which is only evaluated as TRUE (and thus adaptation is run) if at least one model is non-conjugate. Please refer to this link, where it is evident that the presence of an adaptive sampler is needed for adaptive mode.
The original example contains only conjugates and thus no adaptation is performed. This can be seen by running the below code:
rjags::list.samplers(jags_model)
$`bugs::ConjugateNormal`
[1] "mu[24]"
$`bugs::ConjugateNormal`
[1] "mu[23]"
$`bugs::ConjugateNormal`
[1] "mu[22]"
…
Thus, when running a JAGS model, it is important to consider how you want to implement burn-in and adaptation and be aware that the latter is model-driven and not always performed even when specified.
Upvotes: 1