← The Analysis Clinic
Analysis Clinic · Case 004 · Diagnosed
A reviewer or editor challenged my analysis

Reviewer: “Why Didn’t You Use a Multilevel Model?”

Choose an analysis that respects the unit of assignment and dependence, then show what changes in the inference.

Clustered dataMixed modelsR

Symptoms

What this usually means

The concern is dependence, not the name of the model. In a cluster-randomized study, treatment is assigned to clusters; additional participants within those clusters do not create additional independent treatment assignments. A mixed model is one response, but an appropriate cluster-level analysis or suitably corrected marginal analysis may also answer the question.

Common causes

Run these checks

  1. Map assignment, measurement, and sampling units before fitting anything.
  2. Count independent clusters in each arm and describe their sizes; distinguish cluster-level from within-cluster predictors.
  3. Check whether cluster size is related to outcomes and whether the target averages participants or clusters.
  4. Compare effect estimates and uncertainty under an analysis that respects the design.

What not to do

Do not decide from an arbitrary ICC cutoff, or add a random intercept solely to satisfy a reviewer. A nonsignificant variance estimate does not undo cluster randomization. Avoid relying on uncorrected asymptotic inference with very few clusters.

Treatment options

For balanced parallel cluster randomization with a continuous outcome, compare cluster means. For participant covariates or more complex dependence, consider mixed models, GEE, or cluster-robust inference with suitable small-sample corrections. Unequal cluster sizes, informative size, binary outcomes, and observational exposures need additional decisions about weighting and the estimand.

Worked example

This synthetic study has 20 clusters, 25 participants per cluster, and equal allocation. The script compares an independent-participant t test with a t test on the cluster means. Both estimate the intervention-minus-control difference as 0.264.

AnalysisSE95% CIp-value
Independent participants0.1200.029 to 0.4990.0281
Cluster means; 18 df0.436−0.652 to 1.1790.5525

The point estimate agrees because the clusters are equally sized. The uncertainty changes sharply. These simulated results demonstrate the assignment-unit issue; they do not establish that cluster means are optimal for every clustered dataset.

What to tell the reviewers

Adapt the following only after running the corresponding checks:

Treatment was assigned to clusters. We therefore reanalyzed the continuous outcome using cluster means with equal weight in this balanced design. The estimated difference was unchanged, but its confidence interval widened. We revised the inference and reported the cluster counts, sizes, and analysis unit.

See it in R and Python

Python uses the same synthetic participant data as R and reproduces both the independent-participant and cluster-mean t tests.

Python dependencies: NumPy and SciPy. Install with python -m pip install numpy scipy.

Download the shared synthetic CSV beside the Python script before running it. Identical seeds in R and Python do not generally generate identical samples.

# Synthetic balanced cluster-randomized example; base R only.
# Treatment is allocated to clusters. Cluster-mean t test respects that unit.
set.seed(4004)
k <- 20; m <- 25
cluster <- rep(seq_len(k), each=m)
arm <- rep(rep(0:1, each=k/2), each=m)
u <- rnorm(k, sd=1)
y <- 0.5*arm + u[cluster] + rnorm(k*m)
d <- data.frame(cluster, arm, y)
naive <- t.test(y ~ arm, data=d, var.equal=TRUE)
means <- aggregate(y ~ cluster + arm, d, mean)
valid <- t.test(y ~ arm, data=means, var.equal=TRUE)
cat(sprintf('Participants=%d; clusters=%d; cluster size=%d\n', nrow(d),k,m))
cat(sprintf('Difference=%.3f; naive SE=%.3f; cluster SE=%.3f\n', diff(tapply(y,arm,mean)), naive$stderr, valid$stderr))
cat(sprintf('Naive CI=[%.3f, %.3f]; cluster CI=[%.3f, %.3f]\n', -naive$conf.int[2],-naive$conf.int[1],-valid$conf.int[2],-valid$conf.int[1]))
cat(sprintf('Naive p=%.4f; cluster p=%.4f; cluster df=%.0f\n',naive$p.value,valid$p.value,valid$parameter))
stopifnot(valid$stderr > naive$stderr, nrow(means)==k)

# Maintainer-only export; run from repository root. Synthetic data only.
if (Sys.getenv("CLINIC_EXPORT_DATA") == "1") {
write.csv(d,'analysis-clinic/examples/case-004-shared-data.csv',row.names=FALSE,quote=FALSE)
}
# DASS Analysis Clinic Case 004; Python dependencies: numpy, scipy.
from pathlib import Path
import numpy as np
from scipy import stats, optimize
HERE = Path(__file__).resolve().parent

def ols(X, y):
    b = np.linalg.lstsq(X, y, rcond=None)[0]
    residual = y - X @ b
    df = len(y) - X.shape[1]
    cov = (residual @ residual / df) * np.linalg.inv(X.T @ X)
    se = np.sqrt(np.diag(cov))
    ci = np.column_stack((b-stats.t.ppf(.975, df)*se,b+stats.t.ppf(.975,df)*se))
    return b, cov, ci

def power(n, d):
    # Matches R power.t.test(strict=FALSE): rejection tail in effect direction.
    return stats.nct.sf(stats.t.ppf(.975,2*n-2),2*n-2,d*np.sqrt(n/2))

d=np.genfromtxt(HERE/'case-004-shared-data.csv',delimiter=',',names=True)
def compare(a,b):
 effect=b.mean()-a.mean(); df=len(a)+len(b)-2
 pooled=((len(a)-1)*a.var(ddof=1)+(len(b)-1)*b.var(ddof=1))/df
 se=np.sqrt(pooled*(1/len(a)+1/len(b)))
 return effect,se,effect+np.array([-1,1])*stats.t.ppf(.975,df)*se,2*stats.t.sf(abs(effect/se),df)
naive=compare(d['y'][d['arm']==0],d['y'][d['arm']==1])
clusters=np.unique(d['cluster']); means=np.array([d['y'][d['cluster']==k].mean() for k in clusters]); arms=np.array([d['arm'][d['cluster']==k][0] for k in clusters])
valid=compare(means[arms==0],means[arms==1])
print('Independent participants:',naive); print('Cluster means:',valid)
assert np.isclose(naive[0],valid[0]) and valid[1]>naive[1] and np.isclose(valid[1],.436,atol=.0005)

Download Python script · Download shared synthetic CSV

Reproduce this Case

Every number above comes from one base-R script, with no packages to install.

Download case-004-multilevel-model.R →

Sources

← Case 003: Can I Power My Grant on My Pilot's Effect Size?
Case 005: Powering a Cluster-Randomized Grant When the ICC Is Unknown →

Still stuck after the first checks?

Some problems turn on the details of your design, data, or the exact reviewer comment. A free 30-minute consult can identify the next defensible step and what it would take.

Book a free consult