Sensitivity analysis

What a request for a sensitivity analysis is actually asking for

A sensitivity analysis swaps out one analysis decision that could reasonably have gone the other way, reruns, and asks whether the conclusion turns over.

What it varies is an analysis decision. Every word there rules something out, and the three things it rules out get mixed in constantly:

The last row is the one most often confused, so to be explicit: a sensitivity analysis does not need a multiplicity adjustment. It does not ask whether one of these is real; it asks whether the same conclusion survives a change of setting. Adjusting it answers a question about significance when the question was about stability.

The dataset, and what it cannot be used to say

This page uses survival::rotterdam, 2,982 records from a breast cancer registry, of whom 339 received hormonal therapy. The dataset has no missing values at all (0 cells), so the very common sensitivity axis of swapping one missing-data method for another does not exist here; the four axes below are something else.

Sixteen defensible analyses

Four axes, two options each, every combination:

The grid is a full factorial rather than a hand-picked string of settings, because only a factorial design lets the spread be attributed to an axis. A string of settings yields a scatter of points, and which decision is driving the result is the question worth asking.

The table in one sentence: this dataset gives no single answer to whether there is an association. The analysis setting gives the answer.

The axis that drives the curve is the one no methods section mentions

Averaging each axis over the other three and taking the difference gives how far that axis moves the estimate on the log scale:

The dominant axis is whether node count is capped. The lower panel of the figure shows it directly: the dark dots in that row sit almost entirely in the left half.

This is worth pausing on. A methods section states the endpoint and lists the covariates — the two axes a reader can see, ranking second and third here. First place goes to a data-preparation decision that never appears at all.

The main analysis sits mid-curve, and that is the point

The pre-specified main analysis takes the least-processed option on each of the four axes: recurrence-free survival, core covariates, node count as recorded, follow-up complete.

Its hazard ratio is 0.884 (95% CI 0.760–1.028, p = 0.111), not statistically significant, at rank 12 of sixteen. On the same curve, 7 settings are significant.

The main analysis is neither the most conservative nor the most significant setting. Its standing does not come from where it sits on the curve; it comes from having been fixed before the curve existed. Were standing to come from position, then:

• picking the most significant setting could be written up as the specification with the clearest signal
• picking the most conservative one could be written up as the cautious approach

Both sentences are available, and both would be written after seeing the results. Picking a conservative setting after the fact and calling it caution is the same act as picking a significant one.

A robust result is an empty phrase

The sentence appears in papers constantly and carries no information on its own: it does not say which perturbations were run, what each one estimated, or whether any of them changed the conclusion.

Three things belong in the report:

1. Which perturbations were run, listed individually, including the ones that changed nothing
2. The effect estimate and interval under each, not just a p value, not a balance table, and certainly not a single sentence
3. Whether any of them turned the conclusion over, and if so which, and why the main analysis is still the main analysis

Robust does not mean correct

The curve above is what a result that is not robust looks like. The opposite case deserves more caution: sixteen settings all reaching significance would still not establish that the conclusion is correct.

The reason is visible in this page’s own data. Confounding by indication is present in all sixteen settings, and swapping the outcome definition, the covariates, the coding, or the time window does not touch it. A tidy, uniform curve can be tilted as a whole in one direction.

Influence diagnostics makes the same point for meta-analysis: a pooled analysis can be perfectly robust to dropping any single study while every included study shares one systematic bias — because that bias does not disappear when one study leaves.

A sensitivity analysis rules out only what it actually perturbed. About the assumptions it left alone, it has said nothing.

What picking the most significant setting costs, measured

The claim above is that picking afterwards and fixing in advance are different things. How different can be measured.

The method: shuffle the exposure column, which makes the true association exactly zero, and rerun the same 16 settings. Repeat 4,000 times (seed 20260826), recording two numbers each time — the smallest p value, and the p value of the pre-specified setting.

• Reporting only the most significant setting: 17.9% of runs claim significance
• Reporting only the pre-specified setting: 6.5%

The gap is 11.4 percentage points. It needs no inference to read, because the exposure was shuffled, so every significant result above is a false positive.

Same settings, same dataset, same true association of zero. The only thing that changes is which result gets written up. That is why pre-specification is not an administrative formality.

What to look for when reading a paper

1. Where the sensitivity analyses live. A table in the appendix listing an effect estimate per setting and a sentence in the discussion saying the result was robust are entirely different things.
2. Whether there are effect estimates, or only p values and balance tables. With only p values, you cannot see whether the conclusion would turn over.
3. Whether the main and sensitivity roles were fixed in writing. Look for pre-specified, or according to the protocol, and for a registration number.
4. Whether the listed perturbations are only the ones that changed nothing. On paper, a complete set and a filtered set are indistinguishable.
5. What was not perturbed. This usually matters more than what was — see the section above.

Run it yourself

library(survival)

d0 <- rotterdam
d0$rfs_time  <- pmin(d0$rtime, d0$dtime)
d0$rfs       <- as.integer(d0$recur == 1 | d0$death == 1)
d0$nodes_cap <- pmin(d0$nodes, 10)
d0$lpgr      <- log(d0$pgr + 1)

grid <- expand.grid(
  outcome  = c("rfs", "os"),
  covars   = c("core", "plus"),
  nodesf   = c("raw", "cap"),
  followup = c("full", "trunc5"),
  stringsAsFactors = FALSE
)

fit_one <- function(d, g) {
  if (g$followup == "trunc5") {
    tt <- 5 * 365.25
    d$rfs <- ifelse(d$rfs_time > tt, 0L, d$rfs); d$rfs_time <- pmin(d$rfs_time, tt)
    d$death <- ifelse(d$dtime > tt, 0L, d$death); d$dtime <- pmin(d$dtime, tt)
  }
  nod <- if (g$nodesf == "raw") "nodes" else "nodes_cap"
  cov <- if (g$covars == "core") c("age", "size", "grade", nod)
         else c("age", "size", "grade", nod, "meno", "lpgr")
  lhs <- if (g$outcome == "rfs") "Surv(rfs_time, rfs)" else "Surv(dtime, death)"
  f <- as.formula(paste(lhs, "~ hormon +", paste(cov, collapse = " + ")))
  summary(coxph(f, data = d))$coefficients["hormon", c("coef", "Pr(>|z|)")]
}

res <- t(sapply(seq_len(nrow(grid)), function(i) fit_one(d0, as.list(grid[i, ]))))
grid$hr <- exp(res[, 1]); grid$p <- res[, 2]
print(grid[order(grid$hr), ], row.names = FALSE)

import itertools
import numpy as np
import pandas as pd
from lifelines import CoxPHFitter

url = "https://vincentarelbundock.github.io/Rdatasets/csv/survival/rotterdam.csv"
d0 = pd.read_csv(url)
d0["rfs_time"] = np.minimum(d0["rtime"], d0["dtime"])
d0["rfs"] = ((d0["recur"] == 1) | (d0["death"] == 1)).astype(int)
d0["nodes_cap"] = d0["nodes"].clip(upper=10)
d0["lpgr"] = np.log(d0["pgr"] + 1)
size = pd.get_dummies(d0["size"], prefix="size", drop_first=True).astype(float)
d0 = pd.concat([d0, size], axis=1)
size_cols = list(size.columns)


def fit_one(outcome, covars, nodesf, followup):
    d = d0.copy()
    if followup == "trunc5":
        tt = 5 * 365.25
        d["rfs"] = np.where(d["rfs_time"] > tt, 0, d["rfs"])
        d["rfs_time"] = d["rfs_time"].clip(upper=tt)
        d["death"] = np.where(d["dtime"] > tt, 0, d["death"])
        d["dtime"] = d["dtime"].clip(upper=tt)
    nod = "nodes" if nodesf == "raw" else "nodes_cap"
    cov = ["age", "grade", nod] + size_cols
    if covars == "plus":
        cov += ["meno", "lpgr"]
    t, e = ("rfs_time", "rfs") if outcome == "rfs" else ("dtime", "death")
    cph = CoxPHFitter().fit(d[[t, e, "hormon"] + cov], duration_col=t, event_col=e)
    return float(np.exp(cph.params_["hormon"])), float(cph.summary.loc["hormon", "p"])


rows = []
for g in itertools.product(["rfs", "os"], ["core", "plus"], ["raw", "cap"],
                           ["full", "trunc5"]):
    hr, p = fit_one(*g)
    rows.append((*g, hr, p))

out = pd.DataFrame(rows, columns=["outcome", "covars", "nodes", "followup", "hr", "p"])
print(out.sort_values("hr").to_string(index=False))

Common misuses

Rerun every number on this page

/opt/homebrew/bin/Rscript figures/scripts/B8-07-sensitivity-analysis.R

Report a problem on this page