Search code examples
rjagsrjags

JAGS n.adapt Parameter Not Working as Expected


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)



Solution

  • 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.