How to fit a model with and without an interaction in a JAGS regression model

Question

I'm using this tutorial to wrap my head around JAGS code. In the section 'Same model with an additional categorical predictor' it states that "This model includes an interaction between sex and body length". How can I remove this so that there's no interaction?

Here's the full setup and model in R and JAGS.

First the data:

set.seed(42)

samplesize <- 50 # Larger sample size because we're fitting a more complex model
b_length <- sort(rnorm(samplesize)) # Body length
sex <- sample(c(0, 1), size = samplesize, replace = T) # Sex (0: female, 1: male)

int_true_f <- 30 # Intercept of females
int_true_m_diff <- 5 # Difference between intercepts of males and females
slope_true_f <- 10 # Slope of females
slope_true_m_diff <- -3 # Difference between slopes of males and females

mu <- int_true_f + sex * int_true_m_diff + (slope_true_f + sex * slope_true_m_diff) * b_length # True means
sigma <- 5 # True standard deviation of normal distributions

b_mass <- rnorm(samplesize, mean = mu, sd = sigma) # Body mass (response variable)

# Combine into a data frame:
snakes2 <- data.frame(b_length = b_length, b_mass = b_mass, sex = sex)
head(snakes2)

jagsdata_s2 <- with(snakes2, list(b_mass = b_mass, b_length = b_length, sex = sex, N = length(b_mass)))

JAGS code:

lm2_jags <- function(){
    # Likelihood:
    for (i in 1:N){
        b_mass[i] ~ dnorm(mu[i], tau) # tau is precision (1 / variance)
        mu[i] <- alpha[1] + sex[i] * alpha[2] + (beta[1] + beta[2] * sex[i]) * b_length[i]
    }
    # Priors:
    for (i in 1:2){
        alpha[i] ~ dnorm(0, 0.01)
        beta[i] ~ dnorm(0, 0.01)
    }
    sigma ~ dunif(0, 100)
    tau <- 1 / (sigma * sigma)
}

Initial values and run:

init_values <- function(){
    list(alpha = rnorm(2), beta = rnorm(2), sigma = runif(1))
}

params <- c("alpha", "beta", "sigma")

fit_lm2 <- jags(data = jagsdata_s2, inits = init_values, parameters.to.save = params, model.file = lm2_jags,
               n.chains = 3, n.iter = 12000, n.burnin = 2000, n.thin = 10, DIC = F)

score 1 · Accepted Answer · answered Feb 22 '18 at 17:59

The interaction term is contained in your calculation of mu. The sex changes how the formula between body length and body mass is defined, via the slope terms. To build a model where sex and body length are treated as independent with respect to how they affect body mass, you could do something like this:

mu <- int_true_f + (sex * int_true_m_diff) + b_length

The JAGS code would then become

lm2_jags <- function(){
  # Likelihood:
  for (i in 1:N){
    b_mass[i] ~ dnorm(mu[i], tau) # tau is precision (1 / variance)
    mu[i] <- alpha[1] + (sex[i] * alpha[2]) + (b_length[i] * alpha[3])
  }
  # Priors:
  for (i in 1:3){
    alpha[i] ~ dnorm(0, 0.01)
  }
  sigma ~ dunif(0, 100)
  tau <- 1 / (sigma * sigma)
}

Thank you, to clarify why are there 3 alpha values if we're calculating two intercepts? Or is that for the 2 intercepts and 1 common slope? — adkane, Feb 22 '18 at 20:50
The latter, to my understanding. It is for he two coefficients and the common one. — data princess, Feb 22 '18 at 23:51

How to fit a model with and without an interaction in a JAGS regression model

1 Answers1