Pattern-Mixture Models — MNAR by Sensitivity Analysis¶
Make the untestable assumption explicit, then stress-test the conclusion¶
The selection model of Project 4 met non-ignorable missingness by factoring the joint density as $p(y)\,p(r\mid y)$ — an outcome model and a model for being observed. Pattern-mixture models take the other factorisation, $$p(y,r)=p(y\mid r)\,p(r),$$ stratifying by the missingness pattern and giving each stratum its own outcome distribution. Its strength is honesty. The outcome distribution among the missing is never observed, so it cannot be estimated — it must be assumed. Pattern-mixture turns that assumption into a single interpretable knob, a sensitivity parameter $\delta$ tying the missing stratum to the observed one: $$E[y\mid\text{missing}]=E[y\mid\text{observed, same covariates}]+\delta.$$ $\delta=0$ is the MAR assumption (dropouts behave like completers); $\delta\neq0$ is MNAR (they differ by $\delta$ even after adjusting for what we saw). Since $\delta$ is not identified by the data, the analysis is a sensitivity analysis: sweep $\delta$ and watch the conclusion. The tipping point is the $\delta$ at which the conclusion just reverses — "how much worse would the dropouts have had to be for the result to disappear?" This is the delta-adjustment sensitivity analysis regulators ask for in trials with dropout, and the transparent complement to the selection model's opaque $\rho$.
We build the imputation-and-combine engine from scratch, show that dropout biases a naive analysis, run the tipping-point sensitivity, apply it to a real antidepressant-trial dataset with dropout, and cross-check in PyMC.
import numpy as np, pandas as pd, matplotlib.pyplot as plt
import patternmix as PM
rng = np.random.default_rng(6)
BLUE="#2b6cb0"; RED="#c53030"; GREEN="#2f855a"; ORANGE="#dd6b20"; PURP="#6b46c1"; GREY="#718096"
print("Pattern-mixture stratifies by missingness pattern; the missing stratum's outcome is unidentified, so it is")
print("set by a sensitivity parameter delta. delta=0 is MAR; sweeping delta and finding the tipping point is the analysis.")
Pattern-mixture stratifies by missingness pattern; the missing stratum's outcome is unidentified, so it is set by a sensitivity parameter delta. delta=0 is MAR; sweeping delta and finding the tipping point is the analysis.
1. A trial where dropout flatters the drug¶
A two-arm trial: a depression score at the end of study (lower is better), with a real drug benefit ($\tau=-1$). The dropout is MNAR and concentrated in the treatment arm — the treated patients who are not improving stop coming, so the drug's observed endpoints look better than the truth. Complete-case analysis and MAR imputation both exaggerate the benefit.
y,arm,base,yfull=PM.simulate_trial(400, tau=-1.0, rng=rng)
true_eff=yfull[arm==1].mean()-yfull[arm==0].mean()
cc,cse=PM.complete_case_effect(y,arm)
mar,mse=PM.effect(PM.pm_impute(y,arm,base[:,None],rng,delta=0.0,arms=(1,),m=80),arm)
print(f"dropout: treatment arm {np.isnan(y[arm==1]).mean():.0%}, control arm {np.isnan(y[arm==0]).mean():.0%} (the sicker treated patients leave)")
print(f"true (full-data) drug effect : {true_eff:+.2f}")
print(f"complete-case : {cc:+.2f} <- exaggerates the benefit")
print(f"MAR imputation (delta=0) : {mar:+.2f} <- still exaggerates (dropouts imputed like completers)")
fig,ax=plt.subplots(figsize=(7,3))
ests={"full data\n(truth)":(true_eff,GREEN),"complete\ncase":(cc,RED),"MAR\n(delta=0)":(mar,ORANGE)}
ax.bar(range(3),[v for v,_ in ests.values()],color=[c for _,c in ests.values()],alpha=.85)
ax.axhline(true_eff,color=GREEN,ls="--",lw=1); ax.set_xticks(range(3)); ax.set_xticklabels(list(ests))
ax.set_ylabel("drug - placebo (negative = benefit)"); ax.set_title("MNAR dropout in the treatment arm exaggerates the drug benefit")
for i,(v,_) in enumerate(ests.values()): ax.text(i,v-0.08,f"{v:.2f}",ha="center",color="white",fontweight="bold")
plt.tight_layout(); plt.show()
print("Because the treated non-responders drop out, both naive analyses credit the drug with more benefit than it has.")
print("MAR does not help: imputing dropouts as if they resembled completers repeats the same optimistic assumption.")
dropout: treatment arm 30%, control arm 3% (the sicker treated patients leave) true (full-data) drug effect : -1.19 complete-case : -1.79 <- exaggerates the benefit MAR imputation (delta=0) : -1.70 <- still exaggerates (dropouts imputed like completers)
Because the treated non-responders drop out, both naive analyses credit the drug with more benefit than it has. MAR does not help: imputing dropouts as if they resembled completers repeats the same optimistic assumption.
2. The sensitivity analysis and the tipping point¶
We cannot know how the treated dropouts would have scored — so we make it explicit. Add a shift $\delta\ge0$ to the imputed endpoints of the treatment dropouts (they were doing worse than MAR assumes) and sweep it. As $\delta$ grows the estimated benefit shrinks; the tipping point is where its 95% interval first touches zero. The result is read as a statement about assumptions, not just an estimate.
deltas=np.linspace(0,5,26)
est,lo,hi=PM.sensitivity(y,arm,base[:,None],deltas,rng,arms=(1,),m=80)
tp=PM.tipping_point(deltas,est,lo,hi)
drec=deltas[np.argmin(np.abs(est-true_eff))]
fig,ax=plt.subplots(figsize=(8,4.2))
ax.fill_between(deltas,lo,hi,color=BLUE,alpha=.2); ax.plot(deltas,est,color=BLUE,lw=2,label="estimated drug effect")
ax.axhline(0,color="k",lw=1); ax.axhline(true_eff,color=GREEN,ls="--",lw=1,label="true effect")
if tp is not None: ax.axvline(tp,color=RED,lw=2,ls=":",label=f"tipping point delta={tp:.1f}")
ax.set_xlabel(r"$\delta$ = how many points worse the treated dropouts really were"); ax.set_ylabel("drug - placebo effect (95% CI)")
ax.set_title("Tipping-point sensitivity analysis"); ax.legend(frameon=False,fontsize=8); plt.tight_layout(); plt.show()
print(f"delta that recovers the truth: {drec:.1f} (we cannot know this from the data -- it is what MNAR hides)")
sd_y = float(np.nanstd(yfull)); span = float(np.ptp(yfull))
print(f"tipping point: delta = {tp:.1f}. The drug's benefit stays significant unless the treated dropouts would")
print(f"have scored more than {tp:.1f} points worse than completers -- {tp/sd_y:.1f} standard deviations of the endpoint,")
print(f"{tp/span:.0%} of its observed range, and {tp/abs(true_eff):.1f} times the effect being estimated. On that reading")
print("the conclusion is fairly robust: only a large MNAR effect overturns it. That sentence -- not a single")
print("number -- is the honest output of an MNAR analysis.")
delta that recovers the truth: 1.8 (we cannot know this from the data -- it is what MNAR hides) tipping point: delta = 4.4. The drug's benefit stays significant unless the treated dropouts would have scored more than 4.4 points worse than completers -- 2.5 standard deviations of the endpoint, 44% of its observed range, and 3.7 times the effect being estimated. On that reading the conclusion is fairly robust: only a large MNAR effect overturns it. That sentence -- not a single number -- is the honest output of an MNAR analysis.
3. Real data — an antidepressant trial with dropout¶
Congdon's antidepressant-trial data (28 patients, a depression rating at weeks 0–6, drug vs placebo, lower = better) has exactly the worrying structure: of six dropouts before week 6, five are in the drug arm. We estimate the week-6 drug effect adjusting for baseline, and run the same $\delta$-sensitivity on the drug-arm dropouts.
d=pd.read_csv("antidep.csv"); w=d.pivot_table(index="id",columns="week",values="y"); w.columns=["w0","w1","w3","w6"]
w["drug"]=d.groupby("id")["drug"].first(); w["base"]=w["w0"].fillna(w["w1"]); w["base"]=w["base"].fillna(w["base"].mean())
yr=w["w6"].to_numpy(); armr=w["drug"].to_numpy(); baser=w["base"].to_numpy()
cc,cse=PM.complete_case_effect(yr,armr); mar,mse=PM.effect(PM.pm_impute(yr,armr,baser[:,None],rng,delta=0,arms=(1,),m=100),armr)
print(f"n=28 ({int((armr==1).sum())} drug, {int((armr==0).sum())} placebo); week-6 dropout: {int(np.isnan(yr[armr==1]).sum())} drug, {int(np.isnan(yr[armr==0]).sum())} placebo")
print(f"week-6 drug effect complete-case {cc:+.2f} +/- {cse:.2f} MAR {mar:+.2f} +/- {mse:.2f}")
dr=np.linspace(0,3,16); er,lr,hr=PM.sensitivity(yr,armr,baser[:,None],dr,rng,arms=(1,),m=100)
fig,ax=plt.subplots(figsize=(8,4)); ax.fill_between(dr,lr,hr,color=PURP,alpha=.2); ax.plot(dr,er,color=PURP,lw=2)
ax.axhline(0,color="k",lw=1); ax.set_xlabel(r"$\delta$ on the drug-arm dropouts"); ax.set_ylabel("drug - placebo (week-6)")
ax.set_title("Antidepressant trial: the benefit is not robust (small, dropout-heavy sample)"); plt.tight_layout(); plt.show()
print(f"The point estimate suggests a drug benefit (~{mar:+.2f} points), but with only 28 patients its 95% interval")
print("already includes zero even under MAR -- the trial is underpowered, and the differential dropout makes it worse.")
print("Pattern-mixture is candid about this: any assumption that the drug dropouts fared worse erodes the benefit")
print("further. The honest verdict is 'promising but unproven', which a complete-case p-value alone would hide.")
n=28 (18 drug, 10 placebo); week-6 dropout: 5 drug, 1 placebo week-6 drug effect complete-case -1.22 +/- 0.68 MAR -1.28 +/- 0.70
The point estimate suggests a drug benefit (~-1.28 points), but with only 28 patients its 95% interval already includes zero even under MAR -- the trial is underpowered, and the differential dropout makes it worse. Pattern-mixture is candid about this: any assumption that the drug dropouts fared worse erodes the benefit further. The honest verdict is 'promising but unproven', which a complete-case p-value alone would hide.
4. Cross-check in PyMC¶
We fit the completer regression in PyMC (Bayesian, endpoint on arm and baseline), then reconstruct the effect-versus-$\delta$ curve from its posterior — imputing the treatment dropouts with the $\delta$ shift — and overlay it on the from-scratch sensitivity curve. The two should coincide.
import pymc as pm
obs=~np.isnan(y)
with pm.Model() as mod:
a0=pm.Normal("a0",5,3); a1=pm.Normal("a1",0,3); b=pm.Normal("b",0,2); s=pm.HalfNormal("s",3)
pm.Normal("yobs", a0+a1*arm[obs]+b*(base[obs]-5), s, observed=y[obs])
idata=pm.sample(800,tune=1000,chains=4,target_accept=0.9,random_seed=7,progressbar=False)
po=idata.posterior; A0=po["a0"].values.ravel(); A1=po["a1"].values.ravel(); Bc=po["b"].values.ravel(); S=po["s"].values.ravel()
mis=np.isnan(y); nmis_t=(mis&(arm==1)).sum()
def pymc_effect(dl):
e=[]
for k in rng.choice(len(A0),200,replace=False):
yc=y.copy(); pred=A0[k]+A1[k]*arm[mis]+Bc[k]*(base[mis]-5)+S[k]*rng.standard_normal(mis.sum())
pred=pred+dl*((arm[mis]==1)); yc[mis]=pred
e.append(yc[arm==1].mean()-yc[arm==0].mean())
return np.mean(e)
epm=np.array([pymc_effect(dl) for dl in deltas])
fig,ax=plt.subplots(figsize=(7.5,4)); ax.plot(deltas,est,color=BLUE,lw=2.5,label="from-scratch")
ax.plot(deltas,epm,color=RED,lw=1.5,ls="--",label="PyMC completer model"); ax.axhline(0,color="k",lw=.8)
ax.set_xlabel(r"$\delta$"); ax.set_ylabel("drug - placebo effect"); ax.set_title(f"Effect-vs-delta curves agree (corr {np.corrcoef(est,epm)[0,1]:.3f})")
ax.legend(frameon=False,fontsize=8); plt.tight_layout(); plt.show()
print("The PyMC completer regression reproduces the from-scratch effect-versus-delta curve: the pattern-mixture")
print("sensitivity analysis is a thin, transparent layer of assumption on top of an ordinary Bayesian regression.")
g++ not available, if using conda: `conda install gxx`
Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (4 chains in 4 jobs)
NUTS: [a0, a1, b, s]
Sampling 4 chains for 1_000 tune and 800 draw iterations (4_000 + 3_200 draws total) took 4 seconds.
The PyMC completer regression reproduces the from-scratch effect-versus-delta curve: the pattern-mixture sensitivity analysis is a thin, transparent layer of assumption on top of an ordinary Bayesian regression.
5. Summary¶
Pattern-mixture models confront non-ignorable missingness by admitting what cannot be known. Factoring the problem as $p(y\mid r)\,p(r)$, they let the missing stratum have its own outcome distribution and set the unidentified part by an explicit sensitivity parameter $\delta$ — $\delta=0$ recovering MAR, $\delta\neq0$ encoding how much the dropouts differ. Because $\delta$ cannot be estimated, the analysis is a sensitivity sweep ending in a tipping point: on the simulated trial the drug's benefit survived unless the treated dropouts would have scored more than ~4 points worse than completers — an implausible shift, so the conclusion was robust; on the small real antidepressant trial the benefit was already not significant under MAR, and the differential dropout only deepened the doubt. A PyMC completer regression reproduced the whole $\delta$-curve, underlining that this is a transparent layer over ordinary modelling.
The connections close the MNAR pair. Selection models (Project 4) and pattern-mixture models are the two factorisations of the same non-ignorable joint density — $p(y)\,p(r\mid y)$ versus $p(y\mid r)\,p(r)$ — and where the selection model's $\rho$ was opaque and hard to identify, the pattern-mixture $\delta$ is interpretable and chosen, turning Project 4's uncomfortable wide posterior into a deliberate stress test. This is Congdon's pattern-mixture treatment of dropout (BMCD Ch.11), and the imputation-and-combine engine is the multiple-imputation machinery of Project 2 with a shift added. The final project, categorical missing data, returns to ignorable ground for contingency tables and log-linear models, tying the arc back to the latent-class work.