Skip to content
Statistical RethinkingStatistical Rethinking

13.3 Example: Cross-classified chimpanzees with varying slopes

The chimpanzee data contain multiple crossed random effects: actor, block, and treatment. A multilevel model with varying slopes handles this complex structure.

13.3.1 The crossed varying slopes model

Each actor may respond differently to treatments, and each block may have its own baseline. The model includes:

  • Fixed treatment effects: Average effect of each treatment.
  • Varying actor intercepts: Some chimps pull more than others.
  • Varying block intercepts: Some experimental sessions have different baselines.

logit(pi)=α+βtreat[i]+αactor[i]+αblock[i]\text{logit}(p_i) = \alpha + \beta_{\text{treat}[i]} + \alpha_{\text{actor}[i]} + \alpha_{\text{block}[i]}

# Cross-classified chimpanzee model with varying slopes
set.seed(42)
n_actors <- 7
n_blocks <- 6
n_treatments <- 4
n_trials <- 200

# Simulate data

actor <- sample(1:n_actors, n_trials, replace = TRUE)
block <- sample(1:n_blocks, n_trials, replace = TRUE)
treatment <- sample(1:n_treatments, n_trials, replace = TRUE)

# True effects

a_bar <- 0
a_actor <- rnorm(n_actors, 0, 1.5)
a_block <- rnorm(n_blocks, 0, 0.5)
b_treatment <- c(0, 0.5, 0.3, 0.7) # treatment effects

# Generate outcome

logit_p <- a_bar + a_actor[actor] + a_block[block] + b_treatment[treatment]
p <- 1 / (1 + exp(-logit_p))
y <- rbinom(n_trials, 1, p)

d <- data.frame(
actor = factor(actor),
block = factor(block),
treatment = factor(treatment),
y = y
)

# Fit cross-classified model

library(lme4)
m <- glmer(y ~ treatment + (1 | actor) + (1 | block),
family = binomial, data = d)

cat("Cross-classified chimpanzee model:\n")
cat(" Fixed effects (treatment baseline = 1):\n")
fe <- fixef(m)
for (i in seq_along(fe)) {
cat(sprintf(" %s: %.3f\n", names(fe)[i], fe[i]))
}
cat(" Random SD: actor =", round(attr(VarCorr(m)$actor, "stddev"), 3),
  "  block =", round(attr(VarCorr(m)$block, "stddev"), 3), "\n")
Click Run to execute this code in your browser with WebR

13.3.2 Adding varying slopes by actor

We can also let each actor respond differently to the treatment:

logit(pi)=α+βtreat[i]+αactor[i]+βactor[i]⋅treati+αblock[i]\text{logit}(p_i) = \alpha + \beta_{\text{treat}[i]} + \alpha_{\text{actor}[i]} + \beta_{\text{actor}[i]} \cdot \text{treat}_i + \alpha_{\text{block}[i]}

# Cross-classified with varying slopes by actor
set.seed(42)
n_actors <- 7
n_blocks <- 6
n_trials <- 200

actor <- sample(1:n_actors, n_trials, replace = TRUE)
block <- sample(1:n_blocks, n_trials, replace = TRUE)
treatment <- sample(1:n_treatments, n_trials, replace = TRUE)

a_actor <- rnorm(n_actors, 0, 1.5)
a_block <- rnorm(n_blocks, 0, 0.5)
b_treatment <- c(0, 0.5, 0.3, 0.7)

# Actor-specific treatment deviations

b_actor_treat <- matrix(rnorm(n_actors * n_treatments, 0, 0.3), nrow = n_actors)

logit_p <- rep(0, n_trials)
for (i in 1:n_trials) {
logit_p[i] <- a_actor[actor[i]] + a_block[block[i]] +
b_treatment[treatment[i]] + b_actor_treat[actor[i], treatment[i]]
}
p <- 1 / (1 + exp(-logit_p))
y <- rbinom(n_trials, 1, p)

d <- data.frame(
actor = factor(actor),
block = factor(block),
treatment = factor(treatment),
y = y
)

library(lme4)

# Varying intercepts + slopes by actor, varying intercepts by block

m <- glmer(y ~ treatment + (1 + treatment | actor) + (1 | block),
family = binomial, data = d)

cat("Cross-classified with varying slopes:\n")
cat(" Fixed effects:\n")
fe <- fixef(m)
for (i in seq_along(fe)) {
cat(sprintf(" %s: %.3f\n", names(fe)[i], fe[i]))
}
cat(" Actor SD (intercept):", round(attr(VarCorr(m)$actor, "stddev")[1], 3), "\n")
cat("  Block SD (intercept):", round(attr(VarCorr(m)$block, "stddev"), 3), "\n")

cat("\n Varying slopes by actor allow each chimp\n")
cat(" to respond differently to each treatment.\n")
Click Run to execute this code in your browser with WebR
Rethinking

More structure, more shrinkage. Adding varying slopes by actor means the model estimates more parameters, but it does not overfit — because each slope is shrunk toward the population average. An actor with few trials in a treatment will have a slope close to the average. An actor with many trials can deviate more. This is the power of partial pooling: you can add structure without adding uncontrolled complexity. The model learns how much variation there is and adjusts the shrinkage accordingly.

Summary of Section 13.3

  • Cross-classified models handle multiple crossed random effects simultaneously.
  • Actor and block effects are crossed — each actor appears in each block.
  • Varying slopes by actor allow individual chimps to respond differently to treatments.
  • More structure does not mean overfitting — partial pooling controls complexity.
  • Use (1 + treatment | actor) + (1 | block) in lme4.