Causal mediation analysis

This page picks up where the warning stopped

Confounding, DAGs and what to adjust for already covers what not to do: a variable measured after treatment, dropped into the regression, adjusts away part of the very path you wanted to estimate, and what comes out is neither the total effect nor a clean direct effect.

So what should you do instead? The answer is not to leave mediators alone. It is this: if you want to ask how much of a treatment effect travels through X, you first have to define the part that travels through X as a quantity, and only then estimate it. Putting the mediator into the regression fails as a mediation analysis because it skips the defining step — a coefficient does not become a direct effect just because that is what you had in mind.

The moment you are most likely to land here is a Methods section that says mediation analysis was performed using the Baron and Kenny approach, or causal mediation analysis (VanderWeele). After this page you should be able to ask four questions:

1. Was the mediator measured after treatment? If it was measured at baseline it is not a mediator.
2. Was mediator-outcome confounding dealt with? Randomising the exposure does not settle this one.
3. Was exposure-mediator interaction allowed for? Leaving the interaction term out is a stronger assumption, not a safer choice.
4. Is the denominator of the proportion mediated stable? When the total effect sits near zero, that proportion can be any number at all.

The data: a negative trial, and that is the point

The worked example is medicaldata::opt, the Obstetrics and Periodontal Therapy trial — does treating periodontal disease in pregnancy reduce preterm birth. The released dataset holds 823 women.

• Exposure Group: T (treatment) versus C (control), randomised
• Mediator V3.PD.avg: mean periodontal probing depth after treatment, adjusted for the baseline value of the same measure, BL.PD.avg
• Outcome Preg.ended…37.wk: whether the pregnancy ended before 37 weeks

More people drop out than you expect, and the mediator is why

The driver is not the outcome, it is the mediator: only 9 women are missing the outcome, while 139 are missing the post-treatment probing depth.

This is a fragility specific to mediation analysis. An ordinary efficacy analysis needs exposure and outcome; a mediation analysis needs exposure, mediator and outcome all three, so a complete-case analysis necessarily loses more. And the mediator is by construction measured after treatment, so its missingness often tracks the treatment experience itself — nobody measures probing depth for a woman who did not come back — which makes missing-at-random harder to defend. See missing data.

The two arms of the analysis set:

That is 66 preterm births in all, an event rate of 9.7%. The unadjusted risk difference is −0.44 percentage points.

NDE and NIE: turn the question into a quantity first

In counterfactual notation, $Y(a, m)$ is a person’s outcome when treatment is set to $a$ and the mediator is set to $m$, and $M(a)$ is the mediator value that person would have when treatment is set to $a$.

The natural direct effect (NDE) asks: switch treatment on, but hold the mediator at the value it would have taken without treatment — how much does the outcome move?

The natural indirect effect (NIE) asks: keep treatment switched on throughout, and move only the mediator from its untreated value to its treated value.

Add them and the middle term cancels:

This is not an approximation, and it is not a property to be tested. On the risk-difference scale it is an identity that falls out of the definitions. The residual actually produced by this page’s decomposition is 0, which is floating point and nothing else.

The three counterfactual risks on this dataset:

What the decomposition looks like on this trial

The main analysis is g-computation with the interaction term, on the risk-difference scale (units: percentage points):

The intervals are percentile intervals from 2000 resamples, with the mediator model and the outcome model refitted together in every replicate.

So the honest write-up for this trial is: in this analysis set, no indirect effect of periodontal treatment through post-treatment probing depth was detected, and no direct effect was detected either; the total effect likewise did not reach statistical significance.

What randomisation buys, and what it does not

Identifying NDE and NIE from observed data takes these six assumptions:

• consistency
• no exposure-outcome confounding (met by randomisation)
• no exposure-mediator confounding (met by randomisation)
• no mediator-outcome confounding (NOT met by randomisation; the assumption this page is about)
• no exposure-induced mediator-outcome confounder
• correctly specified mediator and outcome models

Randomisation turns two of them into facts of the design: no exposure-outcome confounding and no exposure-mediator confounding. That is a real purchase, and a clean one — treatment was decided by a coin, so nothing can determine both who was treated and whose probing depth ends up deeper.

The fourth one is untouched. The mediator was not randomised: how deep a woman’s probing depth is after treatment happened naturally. Anything that drives both periodontal status and preterm birth — smoking, socioeconomic position, systemic infectious burden — leaves a non-causal component inside the association between probing depth and preterm birth, and that association is precisely path $b$. The coin allocated the treatment, not the periodontal status.

The fifth assumption deserves its own sentence: a confounder of the mediator-outcome relationship that is itself affected by the exposure. You can neither adjust nor not adjust for such a variable — leave it out and path $b$ stays confounded; put it in and, being a descendant of the exposure, it blocks part of the NDE. Under that structure the natural effects stop being identifiable at all, and you have to move to interventional effects, which are beyond this page.

The sixth, correct model specification, is more expensive here than in ordinary regression: both models have to be right, and misspecification of the mediator model propagates through the integral into both effects. Measurement error in the mediator pushes the NIE towards zero in the same way — see measurement error.

Why Baron-Kenny is not enough

The Baron and Kenny steps go like this: fit $M \sim A$ for $a$, fit $Y \sim A + M$ for $b$ and $c'$, then declare the indirect effect to be $a \times b$, the direct effect $c'$, and the total $c' + ab$. In a world of linear models with no interaction, that is correct. Clinical papers usually satisfy neither condition.

First gap: it cannot represent exposure-mediator interaction

The single $b$ of the stepwise approach means one unit more mediator moves the outcome by this much — the same amount in the treated arm and the control arm. But treatment may change how strongly the mediator acts, and that is an $A \times M$ interaction term.

In these data the interaction term is +0.837 on the log odds scale (95% CI −0.442 to +2.098). The interval crosses zero, so it did not reach statistical significance. That is not a licence to treat it as zero, particularly when the interval is wide enough to have opposite signs at its two ends. Keeping the term versus dropping it moves the decomposition like this:

All in percentage points. The total barely moves (+0.10), but the direct and indirect pieces each shift by more than two percentage points, in opposite directions. That is what an interaction term does: it does not change the sum, it changes the split.

So omitting the interaction term is not the conservative option. It is an extra assertion — that treatment does not change how strongly the mediator acts. Either you have a reason for it, or you report both versions.

Second gap: on a non-linear model, coefficients do not multiply

The derivation of $c = c' + ab$ leans on coefficients adding and multiplying, which is a property of linear models. Here the outcome is binary, so the outcome model is logistic regression, and the log odds ratio is non-collapsible: the coefficient from $Y \sim A$ and the coefficient from $Y \sim A + M$ do not agree even with no confounding whatsoever, because they answer different questions — marginal versus mediator-conditional. The full account is in marginal effects and adjusted risk ratios.

Four routines side by side:

The first two rows and the last two rows are not comparable in size; the scales differ. To see how far the shortcut is off you have to compare on one scale — the last two rows. Both are log odds ratios, both target the same quantity, and the only difference is that the VanderWeele closed form carries the interaction term. Baron-Kenny puts the NIE at −0.147; the closed form puts it at −0.352, which is 2.39 times as large. Same data, same scale: the shortcut shrank the indirect effect to less than half.

Proportion mediated: the trouble is the denominator

The definition is simple enough:

On this page the point estimate is 9.89, with a percentile interval of -24.67 to 23.33. A quantity called a proportion, estimated well above 1, with an interval running through negative values — that is not an arithmetic error, it is what this quantity does under these conditions.

The reason is the denominator. The total effect in this trial is −0.35 percentage points, itself short of statistical significance and numerically very close to zero. When the denominator approaches zero, the ratio explodes under any small perturbation of either part, and if the denominator changes sign the ratio flips wholesale.

This page quantifies that rather than describing it. Across 2000 resamples:

Nearly half the resamples flip the sign of the denominator, and more than four fifths produce a proportion that is not between 0 and 1. The interval above therefore cannot be read as a range for a mediated share — it is not a proportion at all.

How to compute it

library(medicaldata)
data(opt, package = "medicaldata")

# The column is Preg.ended...37.wk, with THREE dots. A typo does not error:
# opt$Preg.ends returns NULL and sum(is.na(NULL)) is 0, which reads as "no
# missing values". This line is load-bearing.
stopifnot("Preg.ended...37.wk" %in% names(opt))

# The factor levels are padded: "   ", "No ", "Yes". A `== "No"` without
# trimws() matches nothing at all, silently and with no warning.
yraw <- trimws(as.character(opt[["Preg.ended...37.wk"]]))

d <- data.frame(
  a  = as.integer(opt$Group == "T"),                            # exposure, randomised
  y  = ifelse(yraw == "Yes", 1L, ifelse(yraw == "No", 0L, NA)), # outcome
  m  = opt$V3.PD.avg,                                           # post-treatment mediator
  bl = opt$BL.PD.avg                                            # baseline of the same measure
)
d <- d[complete.cases(d), ]

# Two models. The outcome model carries the a:m term -- dropping it is an extra
# assumption, not a safer default.
mfit <- lm(m ~ a + bl, data = d)
ofit <- glm(y ~ a * m + bl, data = d, family = binomial())

# g-computation: integrate the mediator distribution away to get counterfactual
# risks directly. `astar` picks which arm's mediator distribution to draw from,
# `a` picks which arm the outcome model is evaluated under. They are allowed to
# differ, and that mismatch IS the cross-world quantity E{Y(1, M(0))}.
set.seed(20260823)
R <- 500
ey <- function(a, astar) {
  mu <- predict(mfit, newdata = transform(d, a = astar))
  md <- rnorm(R * nrow(d), rep(mu, R), sigma(mfit))
  mean(predict(ofit, type = "response",
               newdata = data.frame(a = a, m = md, bl = rep(d$bl, R))))
}

y11 <- ey(1, 1); y10 <- ey(1, 0); y00 <- ey(0, 0)
c(NDE = y10 - y00, NIE = y11 - y10, TE = y11 - y00)   # first two add to the third

# There is no closed form for the intervals, so bootstrap. Every replicate must
# refit BOTH models: refitting only the outcome model understates the variance,
# because uncertainty in the mediator model propagates into the NIE.
# The script behind this page integrates with 30-node Gauss-Hermite quadrature
# instead of the random draws above, so the bootstrap resamples the data only
# and no Monte Carlo noise is layered on top.

import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf

opt = sm.datasets.get_rdataset("opt", "medicaldata").data

# The outcome levels are blank-padded ("   ", "No ", "Yes"), so strip first --
# otherwise the map returns an all-NaN column, silently and with no warning.
y = opt["Preg.ended...37.wk"].str.strip()

d = pd.DataFrame({
    "a": (opt["Group"] == "T").astype(int),
    "y": y.map({"Yes": 1.0, "No": 0.0}),
    "m": opt["V3.PD.avg"],
    "bl": opt["BL.PD.avg"],
}).dropna()

mfit = smf.ols("m ~ a + bl", data=d).fit()
ofit = smf.glm("y ~ a * m + bl", data=d, family=sm.families.Binomial()).fit()
sigma = np.sqrt(mfit.scale)

rng = np.random.default_rng(20260823)
R = 500

def ey(a, astar):
    mu = mfit.predict(d.assign(a=astar))
    md = rng.normal(np.tile(mu, R), sigma)
    nd = pd.DataFrame({"a": a, "m": md, "bl": np.tile(d.bl, R)})
    return ofit.predict(nd).mean()

y11, y10, y00 = ey(1, 1), ey(1, 0), ey(0, 0)
print({"NDE": y10 - y00, "NIE": y11 - y10, "TE": y11 - y00})

Sensitivity analysis: the mediator at a different time point

The same trial measured probing depth again at a later visit, V5.PD.avg. Rerunning the identical analysis on it:

In percentage points. NDE and NIE point the same way as in the main analysis and are smaller (the indirect effect is about half the size), and all three again fall short of statistical significance — the indirect effect spans −6.96 to +2.27.

The total effect, however, changes sign: −0.35 percentage points in the main analysis against +0.86 here. These are not two contradictory conclusions. They are the previous section seen from another angle — when the total effect hugs zero, even its sign is unstable, and swapping in a mediator from a different visit is enough to turn it over. It also shows again why the proportion mediated cannot be reported on these data: the sign of the denominator is not yet settled.

How to read the report

A published mediation analysis usually gives you a path diagram plus one table holding NDE, NIE, total effect and proportion mediated. Each of the four opening questions maps onto a specific place in that report:

1. When the mediator was measured. The Methods should say which post-treatment visit. A mediator measured at baseline cannot be a mediator; a mediator measured over a window that overlaps the outcome may be reverse causation instead.
2. Mediator-outcome confounding. Look for a sentence about adjustment in the mediator-outcome relationship, and for a sensitivity analysis for unmeasured confounding. This study was randomised is not an answer, for the reason set out above.
3. Whether an interaction term is present. The modern approach assumes it by default; a paper reporting only $a \times b$ is asserting that the interaction is zero. Where both are given, look at how far apart they are.
4. The denominator behind the proportion mediated. If the total effect’s interval crosses the null, that percentage should not be quoted. When you read that some percentage of the effect was mediated by Y, go back to the total-effect row first.

One more thing that is easy to miss: the scale. NDE and NIE on the risk-difference, risk-ratio and log odds ratio scales are three different sets of numbers, and so is the proportion mediated. A table that does not name its scale cannot be reproduced and cannot be compared with anyone else’s.

Common misuses

How this page relates to the others

• The what-not-to-do lives next door — see confounding, DAGs and what to adjust for. That page covers what happens when the mediator goes into the regression; this one covers what to do instead.
• The interaction term itself — see interaction. The $A \times M$ term here is the same object, except that what it changes is how the effect splits.
• Why coefficients do not multiply — see logistic regression and marginal effects. Non-collapsibility is the technical reason Baron-Kenny is not the main analysis on this page.
• g-computation in general — see marginal effects and target trial emulation. This page only applies it to an estimand with a cross-world term in it.
• Missing data — see missing data. A mediation analysis needs three variables present at once, so complete cases always cost more here.
• Measurement error in the mediator — see measurement error. It pushes the NIE towards zero, so part of an undetected indirect effect may simply be imprecise measurement.
• The clinical setting of these data — see the design chapter on randomised controlled trials.

Rerunning every number on this page

/opt/homebrew/bin/Rscript figures/scripts/B6-09-mediation.R

Report a problem on this page