Missing data and multiple imputation
Complete-case analysis is not a neutral default; it is a choice with an assumption attached. This page uses a dataset that is naturally incomplete to show that the people who get dropped really do differ from the people who stay, then takes multiple imputation apart to reveal what it is — a regression plus noise, followed by two lines of arithmetic.
“Complete-case analysis” is a choice, not a default
coxph(), lm() and glm() silently skip any row with a missing value. No warning, no error — the model’s n is simply smaller than you thought.
That behaviour is called complete-case analysis, and it carries a strong assumption: the people who were dropped do not differ systematically from the people who stayed, in any way relevant to the analysis.
The first thing this page does is check whether that assumption holds in real data.
What the missingness looks like in this dataset
survival::pbc has 418 patients and 20 variables, of which
12 have missing values. The worst offenders:
| Variable | Missing | Proportion missing |
|---|---|---|
trig | 136 | 32.5% |
chol | 134 | 32.1% |
copper | 108 | 25.8% |
trt | 106 | 25.4% |
ascites | 106 | 25.4% |
hepato | 106 | 25.4% |
spiders | 106 | 25.4% |
alk.phos | 106 | 25.4% |
figures/scripts/B8-01-missing-data.RThat block is the point. In this dataset 312 patients entered the trial and 106 declined randomisation and were followed in a registry instead. The trial laboratory tests were part of the protocol, so the registry group has none of those values by definition:
| Test | Missing in the trial group | Missing in the registry group |
|---|---|---|
chol | 28 / 312 (9%) | 106 / 106 (100%) |
copper | 2 / 312 (1%) | 106 / 106 (100%) |
trig | 30 / 312 (10%) | 106 / 106 (100%) |
platelet | 4 / 312 (1%) | 7 / 106 (7%) |
alk.phos | 0 / 312 (0%) | 106 / 106 (100%) |
ast | 0 / 312 (0%) | 106 / 106 (100%) |
platelet is the exception in that table: it has only scattered missing values in both groups (in the registry group, 7 / 106),
so it is ordinary sporadic missingness rather than the structural kind that comes from “the registry group does not have this by definition”. chol and trig also have a non-zero amount missing inside the trial group.
The two that fall entirely within the registry group are alk.phos and ast.
Missingness is not sprinkled at random across the data. It has structure, and that structure is tied to a substantive characteristic of the patient — whether or not they were willing to enter the trial.
Three missingness mechanisms
| Mechanism | What it means | What complete-case analysis does |
|---|---|---|
| MCAR (missing completely at random) | Missingness is unrelated to any variable, including the missing value itself | Unbiased, merely inefficient |
| MAR (missing at random) | Missingness can be explained by the variables you did observe | Not automatically biased: a regression coefficient stays unbiased as long as missingness is independent of the outcome given the model’s covariates, and only loses efficiency; it is biased when that fails, and a marginal quantity such as a mean can be biased even when it holds, because the complete cases are a covariate-selected subgroup. Multiple imputation covers both |
| MNAR (missing not at random) | Missingness depends on the unobserved value itself | Biased; imputation cannot save it, and you need sensitivity analyses |
The people who get dropped are not like the people who stay
The model uses five variables, which leaves 310 complete cases and drops 108 patients (25.8%).
So what do the dropped patients look like?
| Variable | Kept (n = 310) | Dropped (n = 108) | SMD |
|---|---|---|---|
| age | 49.95 ± 10.57 | 53.00 ± 9.78 | 0.299 |
| bili | 3.27 ± 4.54 | 3.08 ± 4.02 | -0.046 |
| albumin | 3.52 ± 0.42 | 3.43 ± 0.43 | -0.222 |
| Observed event rate | 40.0% | 34.3% | — |
The dropped patients are on average 3.0 years older (SMD 0.299, comfortably past the 0.1 balance threshold), their albumin is lower (SMD -0.222), and their observed event rate differs.
What multiple imputation is doing
The idea behind multiple imputation (MI) fits into three sentences:
- Predict the missing values from the other variables, and add random noise — the noise is the essential part, because it represents the fact that we do not know the true value
- Repeat M times to obtain M complete datasets, and run the analysis once in each
- Combine the M results into a single estimate and standard error using Rubin’s rules
The noise in step 1 is routinely misunderstood. Filling in the conditional mean from a regression and nothing else (single imputation) amounts to claiming you know what the value was — the standard error comes out too small and the confidence interval too narrow.
The arithmetic in step 3 is:
- Pooled estimate = the mean of the M estimates
- Total variance = mean within-imputation variance + (1 + 1/M) × between-imputation variance
Within-imputation variance is “how uncertain this estimate would be if the data were complete”; between-imputation variance is “how much extra uncertainty comes from having had to impute”.
library(survival)
data(pbc, package = "survival")
d <- pbc
d$log_bili <- log(d$bili)
d$log_copper <- log(d$copper) # 108 NAs; they are what this page is about
# One imputation: regression prediction + a random residual
impute_once <- function(data, target, predictors) {
obs <- data[!is.na(data[[target]]), ]
fit <- lm(reformulate(predictors, target), data = obs)
need <- is.na(data[[target]])
# The noise is what makes it an imputation; the conditional mean alone
# understates the standard error
data[[target]][need] <- predict(fit, newdata = data[need, ]) +
rnorm(sum(need), 0, sigma(fit))
data
}
M <- 20
est <- se <- numeric(M)
for (m in seq_len(M)) {
dm <- impute_once(d, "log_copper", c("age", "log_bili", "albumin"))
fm <- coxph(Surv(time, status == 2) ~ age + log(bili) + albumin +
log(exp(dm$log_copper)) + factor(stage), data = dm)
est[m] <- summary(fm)$coefficients["log(bili)", "coef"]
se[m] <- summary(fm)$coefficients["log(bili)", "se(coef)"]
}
# Rubin's rules
q_bar <- mean(est) # pooled estimate
u_bar <- mean(se^2) # mean within-imputation variance
b <- var(est) # between-imputation variance
total <- u_bar + (1 + 1/M) * b # total variance
fmi <- ((1 + 1/M) * b) / total # fraction of missing information
# In practice this is one line:
# library(mice); imp <- mice(d, m = 20); pool(with(imp, coxph(...)))Verified with R 4.6.0 and survival 3.8.6. This is deliberately hand-written in base R so that every step of MI is visible; in practice you would use the mice package.
import numpy as np
import statsmodels.api as sm
from statsmodels.imputation.mice import MICEData
pbc = sm.datasets.get_rdataset("pbc", "survival").data
predictors = ["age", "log_bili", "albumin", "log_copper"]
df = pbc.assign(log_bili=lambda x: np.log(x["bili"]),
log_copper=lambda x: np.log(x["copper"]),
status=lambda x: (x["status"] == 2).astype(int))[
predictors + ["time", "status"]]
# MICEData performs chained-equation imputation
imp = MICEData(df)
est, se = [], []
for _ in range(20):
imp.update_all()
fitted = sm.PHReg(imp.data["time"], imp.data[predictors],
status=imp.data["status"]).fit()
# PHReg's params is a numpy array, not a Series, so look the position up:
# a hard-coded params[0] follows the order of `predictors` and hands you
# age rather than log(bili). The R model also adjusts for factor(stage),
# so the two sides are close but not identical.
k = predictors.index("log_bili")
est.append(fitted.params[k])
se.append(fitted.bse[k])
q_bar = np.mean(est)
u_bar = np.mean(np.square(se))
b = np.var(est, ddof=1)
total = u_bar + (1 + 1/len(est)) * b
fmi = ((1 + 1/len(est)) * b) / totalOn the Python side, use statsmodels' MICEData, or sklearn's IterativeImputer (the latter performs a single imputation by default, so you have to repeat and pool it yourself).
The result: the two approaches side by side
| Approach | Patients analysed | HR for log(bilirubin) | 95% CI |
|---|---|---|---|
| Complete-case analysis | 310 | 2.357 | 1.897–2.929 |
| Multiple imputation (M = 20) | 418 | 2.289 | 1.893–2.767 |
The point estimates are close, but the confidence interval is narrower — because imputation brings back the other information those 108 patients carried (their age, bilirubin and albumin were all recorded).
The fraction of missing information (FMI) is 0.052, meaning that of the uncertainty in this estimate, only 5.2% comes from the missingness and the rest comes from the sample itself. When FMI is high, how well the imputation model is specified matters a great deal to the conclusion; when it is low, things are more robust.
There is a third model as well: simply leave copper out. The R script does run it — the same Cox model, but with copper removed entirely: 412 patients, HR 2.459 (2.083–2.903), the narrowest confidence interval of the three.
But it does not belong in the same comparison. The two rows above ask the same question (“adjusting for age, bilirubin, albumin, copper and stage, what is the HR for bilirubin?”) and differ only in how they handle the missingness; dropping copper changes the set of adjustment variables, so it estimates a different quantity. Its interval is narrower not because it is more efficient, but because it was bought by accepting the residual confounding copper would have absorbed — deleting a covariate will of course take n from 310 back up to 412, but you have given up adjusting for it.
Its role here is as a sensitivity analysis: if even changing the adjustment set does not flip the conclusion, then the conclusion is relatively robust to how copper is handled. It is not a third way of handling missing data.
Common misuses
| Misuse | Why it is wrong |
|---|---|
| Fitting the model and letting the package throw away incomplete rows | That is complete-case analysis, a choice with an assumption attached, and it does not warn you |
| Not reporting how many people were excluded, or what they were like | That comparison table is the only basis for judging the direction of the bias |
| Filling in the mean or median | Single imputation understates variance and distorts the correlation structure between variables |
| Imputing only once and analysing that | There is no between-imputation variance, so the standard error is too small |
| Leaving the outcome out of the imputation model | It systematically dilutes the association between exposure and outcome |
| Applying standard MI to MNAR data and calling it handled | MI assumes MAR; MNAR needs sensitivity analyses such as delta adjustment |
| Choosing an approach on the basis of a test for MCAR | MAR and MNAR are indistinguishable in the data; the judgement rests on how the data were generated |
| Treating imputed values as real observations in further analyses | Imputed values carry uncertainty, and only Rubin’s rules give the correct standard error |
Reproducing every number on this page
/opt/homebrew/bin/Rscript figures/scripts/B8-01-missing-data.RRead the figure
The answer comes from the same statistical output that produced this page's figures, not from a number typed in beside them.
The dataset holds 418 rows, and the Cox model prints an analysed count of 310. How should the gap between them be read?
Show the answer and why
Correct answer: Those 108 were dropped silently by the fitting function - complete-case analysis is a choice with an assumption behind it, not a default
The printed count is who stayed. The gap of 108 comes from the fitting function skipping any row with a missing covariate, with no warning and no error. 106 is the number who declined randomisation and entered the registry; the two are close but not the same group, because what was dropped is every row missing at least one of the five covariates this model uses, which merely overlaps heavily with the registry. 124 is the event count among the rows that stayed - a Cox model excludes nobody for not having an event, and censored rows still contribute to the risk sets. Whether dropping these rows matters is decided by how much they resemble the rows that stayed, not by how many went.
The 108 excluded rows and the 310 that stayed come with a comparison table carrying a standardised mean difference for each of three variables. Which statement is right?
Show the answer and why
Correct answer: The standardised mean difference for age is 0.299, well past the usual balance threshold, and those excluded are the older ones
Age at 0.299 is already far from balance: those excluded were about three years older on average. Bilirubin at -0.046 really is close to zero, but a table has to be read across every row - age is imbalanced, albumin is imbalanced, and the two groups differ in observed event rate as well, so letting one row speak for the whole table misses the other two facts. Albumin at -0.222 has its direction reversed: those excluded had lower albumin, not higher, and which group is the minuend decides what the sign means. That is the cell this kind of table is most often misread on. Complete-case analysis does not merely cost sample size, it shifts the cohort's observed event rate too. Almost no paper prints this table, yet it is the only basis for judging whether the exclusions matter.
Multiple imputation splits the uncertainty in two under Rubin's rules. The total variance here is 0.0094. Which statement is right?
Show the answer and why
Correct answer: The within-imputation variance is 0.008878, measuring the uncertainty this estimate would carry with complete data, and what imputation adds on top is comparatively small
The total is the average within-imputation variance plus a factor slightly above 1 times the between-imputation variance. The within part, 0.008878, is the uncertainty the estimate would carry with complete data; the between part, 0.000463, is what imputation adds because the missing values are not known. They differ by more than an order of magnitude, so the total of 0.0094 comes almost entirely from the sample itself. 0.051889 is not a variance at all but the fraction of missing information, the share of the total that the between part accounts for, and reading it as a variance gets the units wrong. One caution: the imputation on this page is deliberately simplified and systematically understates the between-imputation variance, so this fraction demonstrates a mechanism rather than reporting a result about this dataset.
After multiple imputation the confidence interval for the log-bilirubin hazard ratio is narrower than under complete-case analysis. Why?
Show the answer and why
Correct answer: Because imputation brings the analysis back to all 418 rows, so what the discarded rows do carry - age, bilirubin, albumin - returns to the model
The point is that the discarded rows are not blank: age, bilirubin and albumin are all recorded, and only copper or stage is missing. Complete-case analysis throws away a whole row over one missing cell; imputation fills that cell so the row rejoins the 418, and the extra information is where the narrowing comes from. Saying that 310 rows still enter the model is the single-imputation misreading, since the analysed count really does change. The 412 row is a different thing again: it drops copper entirely, which changes the set of covariates and therefore estimates a different quantity. Its interval is narrower not because it is more efficient but because it buys that width by accepting residual confounding, so it cannot be set alongside the other two.
The point estimates from complete-case analysis and from multiple imputation come out close. Some read that as: the missingness sits in a variable unrelated to the exposure, so it does not matter. Does this dataset support that reading?
Show the answer and why
Correct answer: No. Copper and bilirubin have a Spearman correlation of 0.63, so they are not going their separate ways; the closeness came out of the run and was not derivable in advance
0.63 is a moderately strong monotone correlation: when copper is missing, the information it shares with bilirubin goes missing with it, so unrelated to the exposure does not hold in this dataset. The bilirubin row of the comparison table is indeed near zero, but it says something else - the excluded and the retained have similar mean bilirubin, which is neither evidence that the missingness mechanism ignores the exposure nor a guarantee that the adjusted coefficient will not move. 0.26 is the missing fraction, and a fraction on its own settles nothing in either direction: it cannot show the estimate is unbiased, and it cannot guarantee bias the way that option does - the same fraction under a different mechanism gives very different answers, and here the two estimates did come out close. Change the structure - a higher missing rate, or missingness in a variable still strongly related to the outcome after adjustment - and the two approaches can diverge. Closeness is a result, not a premise.
A solid red block spans several laboratory columns near the top of the missingness map. The page says it is not noise. Which statement is right?
Show the answer and why
Correct answer: It is the 106 who declined randomisation and entered registry follow-up; the trial laboratories were never drawn for them by protocol, which makes this structural missingness
The block is 106 multiplied by several laboratory columns. These are the ones who declined randomisation and entered registry follow-up, and those laboratories are part of the trial protocol - so the reason for the missingness can be written down, which is what makes the case for treating it as close to MAR. That case comes from knowing how the data were generated, never from the data themselves. 312 is the enrolled count, and their stretch of the map is largely white; reading the red block as the trial arm reverses the mechanism entirely. Platelet count is the exception in the by-group table: it is missing for only a handful of registry rows, which is ordinary scattered missingness and far too little to fill the block. Whether missingness has a structure, and what that structure is, decides whether complete cases, imputation or a sensitivity analysis is the right next step.
Chapters that use this method
Watch next
Survival Analysis [Simply Explained]Sources and licences
This page is original writing