Skip to content
Statistical RethinkingStatistical Rethinking

12.1 Example: Multilevel tadpoles

Multilevel models are models within models. The data are grouped, and each group has its own parameters. The group parameters are drawn from a common distribution, which links the groups together. This partial pooling shares information across groups.

12.1.1 The tadpole example

Reed frog tadpoles are raised in tanks with varying densities. Survival varies by tank, but all tanks share a common environment. We model each tank's survival as a draw from a common distribution:

logit(pj)=α+γj\text{logit}(p_j) = \alpha + \gamma_j γj∼Normal(0,σγ)\gamma_j \sim \text{Normal}(0, \sigma_\gamma)

where jj indexes tanks, α\alpha is the population average, and γj\gamma_j is the tank-specific deviation.

# Multilevel model: tadpole survival
# Simulate data inspired by the reed frog example

set.seed(42)
n_tanks <- 48
density <- sample(c(10, 25, 35), n_tanks, replace = TRUE)

# Population-level survival

a_pop <- 0.5 # log-odds
sigma_tank <- 1.5 # between-tank variation

# Tank-specific deviations

gamma_tank <- rnorm(n_tanks, 0, sigma_tank)
p_tank <- 1 / (1 + exp(-(a_pop + gamma_tank)))
survival <- rbinom(n_tanks, density, p_tank)

d <- data.frame(tank = 1:n_tanks, density = density, survival = survival)

cat("Tadpole survival data:\n")
cat(" Tanks:", n_tanks, "\n")
cat(" Densities:", paste(sort(unique(density)), collapse = ", "), "\n")
cat(" Survival proportions range:",
round(min(survival/density), 2), "to", round(max(survival/density), 2), "\n")

# No-pooling estimates (separate per tank)

p_nopool <- survival / density

# Complete-pooling estimate (one common p)

p_pool <- sum(survival) / sum(density)

cat("\nComplete pooling p:", round(p_pool, 3), "\n")
cat("No pooling range:", round(min(p_nopool), 3), "to", round(max(p_nopool), 3), "\n")

# Plot

plot(1:n_tanks, p_nopool, pch = 16, col = "steelblue",
ylim = c(0, 1), xlab = "Tank", ylab = "Survival proportion",
main = "Tadpole survival: no-pooling estimates")
abline(h = p_pool, lwd = 2, col = "firebrick", lty = 2)
legend("topright", legend = c("No pooling", "Complete pooling"),
col = c("steelblue", "firebrick"), pch = c(16, NA), lty = c(NA, 2), lwd = c(NA, 2))
Click Run to execute this code in your browser with WebR

12.1.2 Fitting the multilevel model

The multilevel model estimates both the population average α\alpha and the between-tank variation σγ\sigma_\gamma, while producing tank-specific estimates that are shrunk toward the population mean.

# Fit multilevel model using lme4
library(lme4)

# Binomial multilevel model with tank-level random intercepts

m <- glmer(cbind(survival, density - survival) ~ (1 | tank),
family = binomial, data = d)

cat("Multilevel model:\n")
cat(" Fixed intercept (log-odds):", round(fixef(m)[1], 3), "\n")
cat(" Population p:", round(plogis(fixef(m)[1]), 3), "\n")
cat(" Between-tank SD (log-odds):", round(attr(VarCorr(m)$tank, "stddev"), 3), "\n")

# Extract tank-specific estimates (partial pooling)

p_partial <- plogis(fixef(m)[1] + ranef(m)$tank[, 1])
p_nopool <- d$survival / d$density
p_pool <- sum(d$survival) / sum(d$density)

cat("\nComparison (first 10 tanks):\n")
cat(" Tank No-pool Partial Complete\n")
for (i in 1:10) {
cat(sprintf(" %4d %7.3f %7.3f %7.3f\n",
i, p_nopool[i], p_partial[i], p_pool))
}

# Plot: shrinkage

plot(p_nopool, p_partial, pch = 16, col = "steelblue",
xlim = c(0, 1), ylim = c(0, 1),
xlab = "No-pooling estimate", ylab = "Partial pooling estimate",
main = "Shrinkage: no-pooling vs partial pooling")
abline(0, 1, lty = 2, col = "gray50")
abline(h = p_pool, lty = 3, col = "firebrick")
abline(v = p_pool, lty = 3, col = "firebrick")

cat("\nShrinkage: partial pooling estimates are pulled toward\n")
cat("the population mean, especially for tanks with few tadpoles.\n")
Click Run to execute this code in your browser with WebR
Rethinking

Partial pooling is adaptive regularization. The multilevel model automatically regularizes: tanks with few observations are shrunk more toward the population mean, while tanks with many observations are shrunk less. The amount of shrinkage is determined by the data — specifically, by the ratio of within-tank to between-tank variation. If tanks are very different (large σγ\sigma_\gamma), there is little shrinkage. If tanks are similar (small σγ\sigma_\gamma), there is a lot of shrinkage. This is adaptive regularization: the model learns how much to regularize from the data itself.

Summary of Section 12.1

  • Multilevel models have parameters within parameters: group-level effects drawn from a population distribution.
  • Partial pooling shares information across groups — better estimates for small groups.
  • The tadpole example: logit(pj)=α+γj\text{logit}(p_j) = \alpha + \gamma_j, γj∼Normal(0,σγ)\gamma_j \sim \text{Normal}(0, \sigma_\gamma).
  • Use lme4::glmer() for binomial multilevel models in R.
  • Partial pooling = adaptive regularization — shrinkage amount is learned from data.