Simpson Reversals In Marginal GLMM Comparisons

Why a Population-Averaged Comparison Can Reverse

An exposure can be associated with a lower outcome probability in every covariate stratum, yet with a higher outcome probability when each exposure group is averaged over its own covariate distribution. This is a Simpson reversal: the within-stratum comparisons and the unadjusted population averages point in opposite directions because the exposure groups contain different mixes of strata.

In a GLMM, the distinction between conditional and marginal quantities matters as well. A conditional prediction fixes the random effect; a marginal probability averages over its distribution. For the reversal studied here, compare exposure groups using stratum-specific probabilities that have already been marginalized over the random effect. The remaining source of reversal is the different distribution of the covariate strata across exposure groups.

Definition: Let \(A\in\{0,1\}\) denote exposure, \(Z\in\{1,\ldots,K\}\) a categorical covariate defining strata, and \(Y\in\{0,1\}\) a binary outcome. Define \(m_{ak}=P(Y=1\mid A=a,Z=k)\), where the probability is marginal over the GLMM random effects. Let \(w_{ak}=P(Z=k\mid A=a)\). The group-specific marginal probability is \(M_a=\sum_{k=1}^K w_{ak}m_{ak}\); the crude marginal risk difference is \(M_1-M_0\). A standardized probability instead uses the same target weights \(q_k\) for both exposure levels: \(M_a^{(q)}=\sum_{k=1}^K q_km_{ak}\), where \(q_k\geq 0\) and \(\sum_k q_k=1\).

The quantities answer different questions. \(M_1-M_0\) compares the two groups as they are composed in the population being summarized. \(M_1^{(q)}-M_0^{(q)}\) compares exposure levels after placing both on the same specified covariate distribution. Neither automatically has a causal interpretation: that requires an appropriate causal design or assumptions, including adequate control of confounding.

GLMM Setup and the Key Decomposition

For example, a binary-outcome random-intercept model with a stratum-by-exposure interaction can be written

Y_{ij}\mid b_j,A_{ij},Z_{ij}\sim \operatorname{Bernoulli}(\pi_{ij}),\qquad \operatorname{logit}(\pi_{ij})=\beta_0+\beta_A A_{ij}+\beta_Z Z_{ij}+\beta_{AZ}A_{ij}Z_{ij}+b_j,\qquad b_j\sim N(0,\tau^2).

Here \(i\) indexes individuals and \(j\) indexes clusters. The Bernoulli outcomes are assumed conditionally independent given the covariates and random effects, and the random intercepts are assumed to follow the stated distribution. For categorical \(Z\) with more than two levels, the covariate and interaction terms are represented by the corresponding indicator variables. To obtain \(m_{ak}\), hold \(A=a\) and \(Z=k\) fixed and integrate the response probability over the random-intercept distribution. The stratum-specific marginal probabilities need not be calculated from a single fixed-effect coefficient; they are response-scale quantities.

The decomposition below explains precisely how the reversal can occur. Add and subtract \(\sum_k w_{1k}m_{0k}\) in the crude contrast:

\begin{aligned} M_1-M_0 &=\sum_k w_{1k}m_{1k}-\sum_k w_{0k}m_{0k}\\ &=\sum_k w_{1k}(m_{1k}-m_{0k}) +\sum_k(w_{1k}-w_{0k})m_{0k}. \end{aligned}

The first term averages the within-stratum risk differences using the exposed group's stratum weights. The second is a composition term: it measures how the exposure groups' stratum weights differ, evaluated against the unexposed probabilities. A sufficiently large positive composition term can outweigh negative within-stratum differences. This decomposition is an identity, not an assumption that composition and exposure have separate causal effects; the numerical split depends on which group's probabilities are used as the reference.

By contrast, for a common target population,

M_1^{(q)}-M_0^{(q)} =\sum_k q_k(m_{1k}-m_{0k}).

If \(m_{1k}-m_{0k}<0\) in every stratum with positive target weight, then every term in this sum is nonpositive and at least one is negative. Therefore the standardized risk difference is negative. A common target distribution cannot reverse the direction of stratum-specific risk differences when they all have the same sign. This conclusion concerns risk differences; other effect scales require interpreting their own contrasts.

Worked Example: A Crude Risk Difference Reverses Direction

Suppose a hypothetical GLMM for a binary clinical outcome yields the following fitted probabilities, marginalized over the random-intercept distribution. The covariate \(Z\) has two levels, lower-risk and higher-risk. The table also gives the numbers of individuals in each exposure-by-stratum group in the population being summarized.

Stratum	Unexposed: number	Unexposed probability \(m_{0k}\)	Exposed: number	Exposed probability \(m_{1k}\)
Lower-risk	450	0.10	100	0.08
Higher-risk	50	0.60	400	0.55

Within the lower-risk stratum, the exposed probability is \(0.08-0.10=-0.02\), or 2 percentage points lower. Within the higher-risk stratum, it is \(0.55-0.60=-0.05\), or 5 percentage points lower. Both stratum-specific comparisons favor lower outcome probability among the exposed.

But 400 of 500 exposed individuals, or \(400/500=0.80\), are in the higher-risk stratum. Only \(50/500=0.10\) of unexposed individuals are in that stratum. Thus

M_1=(0.20)(0.08)+(0.80)(0.55)=0.016+0.440=0.456, \qquad M_0=(0.90)(0.10)+(0.10)(0.60)=0.090+0.060=0.150.

The crude difference is \(0.456-0.150=0.306\): the exposed group's average fitted outcome probability is 30.6 percentage points higher. This does not contradict the lower exposed probability in each stratum. It reflects the much larger share of exposed individuals in the higher-risk stratum. The group-specific averages describe these different mixes, not a comparison standardized to one shared covariate distribution.

Worked Example: Standardizing Both Groups to One Target

Use the fitted probabilities in the first example, but choose a hypothetical target population with equal proportions in the two strata: \(q_{\text{lower}}=q_{\text{higher}}=0.50\). The standardized exposed probability is

M_1^{(q)}=(0.50)(0.08)+(0.50)(0.55)=0.040+0.275=0.315.

The standardized unexposed probability is

M_0^{(q)}=(0.50)(0.10)+(0.50)(0.60)=0.050+0.300=0.350.

Therefore the standardized risk difference is \(0.315-0.350=-0.035\). For this equal-stratum target, the exposed probability is 3.5 percentage points lower. This direction agrees with both within-stratum risk differences. The answer depends on the chosen target weights: equal weights are used here to demonstrate the calculation, not because they must be the scientifically appropriate population.

In applied work, the target weights should correspond to the population addressed by the research question. For instance, one might use the proportions in a defined eligible population. The same fitted stratum probabilities can yield different standardized averages for different target populations, even though the within-stratum comparisons stay fixed.

Checking the Arithmetic and the Interpretation

Worked Example: Separating Within-Stratum and Composition Terms

For the first example, the exposed weights are \(w_{1,\text{lower}}=0.20\) and \(w_{1,\text{higher}}=0.80\); the unexposed weights are \(w_{0,\text{lower}}=0.90\) and \(w_{0,\text{higher}}=0.10\). The within-stratum component of the decomposition is

$$ (0.20)(0.08-0.10)+(0.80)(0.55-0.60) =(0.20)(-0.02)+(0.80)(-0.05)=-0.004-0.040=-0.044. $$

The composition component, using the unexposed probabilities as the reference, is

$$ (0.20-0.90)(0.10)+(0.80-0.10)(0.60) =(-0.70)(0.10)+(0.70)(0.60)=-0.070+0.420=0.350. $$

Adding the terms gives \(-0.044+0.350=0.306\), exactly the crude difference calculated earlier. The within-stratum component is negative, but the composition component is larger and positive, producing a positive crude contrast. The decomposition makes the reversal transparent; it does not establish that the exposure caused either component.

Several checks help determine whether an apparent reversal is meaningful or is instead a warning about the analysis:

Check stratum definitions and cell estimates. The strata should represent the covariate comparisons intended, and each exposure group should have enough information in each relevant stratum. Sparse or empty cells can make fitted probabilities unstable or require model-based extrapolation.
Check the weighting population. State whether the reported average uses each group's observed covariate distribution or a common target distribution. A reversal can be entirely due to using different weights.
Check model adequacy. The GLMM's mean structure, including relevant interactions, should be defensible. Inspect fit and residual diagnostics appropriate to the outcome and model, and assess whether predicted probabilities are plausible across exposure-stratum cells. A misspecified interaction structure can obscure genuine stratum differences.
Check the random-effects assumptions and target. Integration uses a specified random-effects distribution. If cluster effects are associated with exposure or covariates, assuming a common independent random-effects distribution may not represent the intended population average.
Check overlap and interpretation. When exposure groups have little covariate overlap, a common-target comparison may rely on extrapolation. Report the target population and avoid presenting an associational contrast as a causal effect without the required design and assumptions.

Comparison	Weights applied	Question answered	Key caution
Crude group-specific average	Each group uses its own \(w_{ak}\)	What is the average fitted probability in each observed exposure group?	Different covariate mixes can reverse the direction.
Stratum-specific comparison	Compare \(m_{1k}\) with \(m_{0k}\) within each stratum	How do fitted probabilities differ at each covariate level?	Does not by itself summarize one target population.
Standardized comparison	Both exposure levels use the same \(q_k\)	What are the probabilities and contrast in a specified target population?	Depends on target choice and adequate support for the comparisons.

A Simpson reversal is not the same as odds-ratio noncollapsibility. Noncollapsibility can make a marginal odds ratio differ from a conditional odds ratio even when there is no confounding or covariate imbalance. The examples here demonstrate a different mechanism: unequal stratum weights make the two group-specific averages refer to different mixtures. On the probability scale, the mixture calculation is explicit.

Common Mistake

A common error is to see a crude exposed-group probability above the unexposed probability and conclude that exposure raises risk in every type of individual—or, conversely, to see lower probabilities in every stratum and claim that the crude difference must also be negative. Neither follows. The crude averages use different weights, so their difference combines within-stratum contrasts and composition. State which comparison is being reported, calculate its weights, and use a common target distribution when the intended question requires one.

Key takeaway: In a GLMM, first define the response-scale probabilities and the population over which random effects are marginalized. Then distinguish averaging each exposure group over its own covariate distribution from standardizing both groups to a common target. Unequal stratum weights can reverse a crude contrast even when every stratum-specific risk difference has the same sign.

Check Your Understanding

Work through these questions using risk differences and probabilities.

In the first example, what are the lower-risk and higher-risk stratum-specific risk differences?
If a target population has 75% lower-risk and 25% higher-risk individuals, calculate the standardized probability for each exposure group and the standardized risk difference.
In the first example's decomposition, why can the crude contrast be positive even though the within-stratum component is negative?
What is the key difference between a group-specific crude average and a standardized average?

Solutions.

1. The lower-risk difference is \(0.08-0.10=-0.02\), or 2 percentage points lower among the exposed. The higher-risk difference is \(0.55-0.60=-0.05\), or 5 percentage points lower.

2. With weights \(0.75\) and \(0.25\), the exposed probability is \((0.75)(0.08)+(0.25)(0.55)=0.1975\). The unexposed probability is \((0.75)(0.10)+(0.25)(0.60)=0.2250\). The standardized risk difference is \(0.1975-0.2250=-0.0275\), or 2.75 percentage points lower in the exposed group for this target.

3. The within-stratum component is \(-0.044\), while the composition component is \(0.350\). Their sum is \(0.306\), so the positive composition term more than offsets the negative within-stratum term.

4. A crude average weights each exposure group's stratum probabilities by that group's own stratum proportions. A standardized average uses the same chosen target weights for both exposure levels, so it compares them in a common covariate mix.

Simpson Reversals In Marginal GLMM Comparisons

What You'll Learn