第十章:把低疗效留成空白

本章造假目标

八名受试者的三个月随访已经完成,数值越高表示状态越好。治疗组第 3 月的四个结果依次为 10、9、4、3,对照组依次为 7、6、8、9;完整数据下两组均值为 6.56.57.57.5,治疗组落后 1 分。结果表只按第 3 月完整病例计算,验收要求治疗组至少领先 1.51.5 分。允许的操作只有把治疗组的第 3 月结果改记为“缺失”,不得改动任何可见数值、早期随访、组别和计划人数。任务是求出通过验收所需的最少空白数,构造最终报告,再用缺失指示、选择偏差、逆概率加权、多重插补和敏感性曲线恢复被空白遮住的信息。

缺失会同时改变信息量和进入统计量的观测集合。若缺失位置与结局有关,完整病例均值就带有选择条件。本章先完成两格空白的构造,再把这个有限样本操作写成概率模型。

两格为什么是最少空白

后台完整随访表为

编号组别基线第 1 月第 2 月第 3 月
T1治疗78910
T2治疗6789
T3治疗8864
T4治疗9753
C1对照7777
C2对照6666
C3对照8888
C4对照9999

本章数值任务的目标量是所有计划受试者在第 3 月的组间均值差。完整数据给出

θ^full=YˉTYˉC=10+9+4+347+6+8+94=6.57.5=1.\widehat\theta_{\mathrm{full}} =\bar Y_T-\bar Y_C =\frac{10+9+4+3}{4} -\frac{7+6+8+9}{4} =6.5-7.5=-1.

固定隐藏数量 kk 后,要让剩余治疗均值尽可能大,就应保留最大的 4k4-k 个结果:若保留值小于某个隐藏值,交换二者便会提高均值。各个候选值为

空白数 kk保留的治疗结果最大治疗均值最大报告效应
0010,9,4,310,9,4,36.56.51-1
1110,9,410,9,423/37.66723/3\approx7.6671/60.1671/6\approx0.167
2210,910,99.59.522

只留一个空白时,隐藏最小值 3 已经给出全部单次操作中的最大效应,仍未达到 1.51.5。再隐藏结果 4,完整病例报告变为

θ^cc=10+927+6+8+94=9.57.5=2.\widehat\theta_{\mathrm{cc}} =\frac{10+9}{2}-\frac{7+6+8+9}{4} =9.5-7.5=2.

两格通过验收;单格最优值仍不合格,所以两格也是最少数量。最终结果表若只展示均值,会写成治疗组 9.59.5、对照组 7.57.5、组间差 2.02.0。相应分母是 2 与 4,必须与均值同时报告。

造假动作留下的影子

可见疗效数值保持原样,变化发生在“谁进入第 3 月均值”。两个空白全部落在治疗组,并且对应第 2 月状态最低的两人;治疗组观测率为 50%50\%,对照组为 100%100\%。后台导出、随访日志、数据版本与个体时间线能够恢复 T3、T4 的身份及原始结果。

怎样把空白位置写成数据

一个结果看不见时,它的值未知;该位置是否可见仍然可以记录。对第 3 月结果 YiY_i 定义观测指示

Ri={1,Yi 被观测,0,Yi 缺失.R_i= \begin{cases} 1,&Y_i\text{ 被观测},\\ 0,&Y_i\text{ 缺失}. \end{cases}

按照 T1–T4、C1–C4 的顺序,修改后的指示向量为

R=(1,1,0,0,1,1,1,1).\mathbf R=(1,1,0,0,1,1,1,1).

治疗组与对照组的观测率分别为

p^obs,T=24=0.5,p^obs,C=44=1.\widehat p_{\mathrm{obs},T}=\frac24=0.5, \qquad \widehat p_{\mathrm{obs},C}=\frac44=1.

再比较治疗组第 2 月状态:

Y2,T,R=1=9+82=8.5,Y2,T,R=0=6+52=5.5.\overline Y_{2,T,R=1} =\frac{9+8}{2}=8.5, \qquad \overline Y_{2,T,R=0} =\frac{6+5}{2}=5.5.

空白位置与组别及已观测的早期状态存在明显对应。由 (X,R,Yobs)(\mathbf X,\mathbf R,\mathbf Y_{\mathrm{obs}}) 可以直接看见这种缺失图案;解释它怎样产生,还需要指定缺失机制

Pr(RY,X),\Pr(\mathbf R\mid \mathbf Y,\mathbf X),

其中 X\mathbf X 包含组别、早期状态及其他完整记录,Y\mathbf Y 表示本来应当取得的全部结局。一次实现的 R\mathbf R 可以与多种概率机制相容,因而图案与机制需要分开表述。

适用边界:纵向研究还要区分单调与非单调图案

若某人从时点 tt 起退出,随后结果都缺失,则 Rit=0R_{it}=0 会推出 Ri,t+1=0R_{i,t+1}=0,这种图案称为单调缺失。间断缺访后又恢复随访会产生非单调图案。本例第 2 月全部可见,第 3 月有两人缺失,属于单调图案;这项图案描述本身不决定 MCAR、MAR 或 MNAR。

三种缺失机制分别限制哪项依赖

给定一个缺失图案 r\mathbf r,把完整结局分成可见部分 Yobs\mathbf Y_{\mathrm{obs}} 与缺失部分 Ymis\mathbf Y_{\mathrm{mis}}。三类机制规定 R\mathbf R 的条件分布可以依赖哪些信息。

MCAR:缺失与数据状态无关

完全随机缺失(missing completely at random, MCAR)要求

Pr(R=rYobs,Ymis,X)=Pr(R=r).\Pr(\mathbf R=\mathbf r \mid \mathbf Y_{\mathrm{obs}},\mathbf Y_{\mathrm{mis}},\mathbf X) =\Pr(\mathbf R=\mathbf r).

例如,实验室在随机时刻断电,受损样本与组别、既往状态和检测值都无关。MCAR 下,完整病例可看作原样本的随机子样本,主要损失通常是精度。

本例空白集中在治疗组的低早期状态者,这种实现会削弱 MCAR 的可信度;单个有限样本仍可能偶然产生不均衡图案,不能只凭观测率差异在逻辑上排除 MCAR。

MAR:可见信息解释缺失概率

随机缺失(missing at random, MAR)要求在给定可见信息后,缺失图案的概率不再随缺失结局改变:

Pr(R=rYobs,Ymis,X)=Pr(R=rYobs,X).\Pr(\mathbf R=\mathbf r \mid \mathbf Y_{\mathrm{obs}},\mathbf Y_{\mathrm{mis}},\mathbf X) = \Pr(\mathbf R=\mathbf r \mid \mathbf Y_{\mathrm{obs}},\mathbf X).

假设退出规则只读取组别和第 2 月结果,例如“治疗组第 2 月不高于 6 时停止随访”。T3、T4 的第 2 月值为 6 与 5,这些可见信息已经决定缺失概率;这样的规则可以满足 MAR。MAR 允许缺失与已观测数据有很强关系。

MNAR:缺失仍依赖不可见结局

非随机缺失(missing not at random, MNAR)允许在给定全部可见信息后,R\mathbf R 仍随 Ymis\mathbf Y_{\mathrm{mis}} 改变。本章规定的操作先读取第 3 月结果,再把其中的 4 与 3 改为空白;选择规则直接使用了随后不可见的值,对应 MNAR 机制。

从最终可见表出发,分析者只知道两个结果缺失,无法观察选择规则是否读取过 4 与 3。联系记录、退出原因、辅助变量和敏感性分析因此承担关键作用。

方法来路:Rubin 为什么研究缺失机制

Rubin 在 1976 年研究缺失过程何时可以从统计推断中忽略。他给出的条件把数据模型与缺失模型连接起来:对似然或贝叶斯推断,MAR 还需配合两套参数的相异性(distinctness),使数据模型参数与缺失机制参数能够分别变化。MCAR、MAR 与不可忽略缺失随后成为描述研究流程和分析假设的共同语言。

适用边界:MAR 不能替任何分析方法自动担保

MAR 是关于缺失机制的条件独立假设。完整病例、加权、似然和插补使用这一假设的方式不同,还各自需要模型设定、正值性或参数可区分性等条件。观测数据通常也无法直接检验 MAR 对缺失值所施加的那部分限制。

完整病例均值偏向哪一部分人

令组别变量为 G{T,C}G\in\{T,C\}。以治疗组为例,目标均值是

μT=E(YG=T),\mu_T=\mathbb{E}(Y\mid G=T),

完整病例均值对应的总体量是

μT,obs=E(YR=1,G=T).\mu_{T,\mathrm{obs}} =\mathbb{E}(Y\mid R=1,G=T).

pT=Pr(R=1G=T)p_T=\Pr(R=1\mid G=T),并记

μT,mis=E(YR=0,G=T).\mu_{T,\mathrm{mis}} =\mathbb{E}(Y\mid R=0,G=T).

全期望公式给出

μT=pTμT,obs+(1pT)μT,mis.\mu_T =p_T\mu_{T,\mathrm{obs}} +(1-p_T)\mu_{T,\mathrm{mis}}.

移项后得到完整病例均值相对目标均值的选择偏差:

μT,obsμT=(1pT)(μT,obsμT,mis).\mu_{T,\mathrm{obs}}-\mu_T =(1-p_T) \left( \mu_{T,\mathrm{obs}}-\mu_{T,\mathrm{mis}} \right).

后台完整数据中,治疗组完整者与缺失者的末次均值分别为

YT,R=1=9.5,YT,R=0=4+32=3.5.\overline Y_{T,R=1}=9.5, \qquad \overline Y_{T,R=0}=\frac{4+3}{2}=3.5.

观测率为 pT=1/2p_T=1/2,所以本例的有限样本差值满足

9.56.5=12(9.53.5)=3.9.5-6.5 =\frac12(9.5-3.5) =3.

完整病例治疗均值比四人的完整均值高 3 分。组间效应也随之从 1-1 上升到 22,整整增加 3 分。

适用边界:完整病例是否有效取决于目标参数

若在每个组内都有 YRY\perp R,组内完整病例均值可以无偏估计组内目标均值。某些回归系数在更宽的条件下也可能由完整病例有效估计。判断依据是目标参数、分析模型与缺失依赖的组合,单独报告一个机制标签还不够。

为什么可见数据无法确定缺失者均值

最终表对治疗组只显示 10、9 和两个空白。以下三种补全拥有完全相同的可见部分:

对缺失者的假设缺失者均值完整治疗均值组间效应
与完整者均值相同9.59.59.59.522
缺失者均值为 5.55.55.55.57.57.500
后台真实值 4,34,33.53.56.56.51-1

可见数据只确定结果 10、9 和两个空白的位置,并未确定空白中结局的分布。关于缺失者的结论必须来自额外假设、辅助观测或外部证据。这种“多个完整数据分布对应同一观测数据分布”的现象称为不可识别性。后面的加权和插补都要明确说明它们补入了哪项假设。

逆概率加权怎样校正不同观测机会

XiX_i 包含组别、早期结果及其他完整记录。在 MAR 下,若这些变量足以解释观测概率,可写成

πi=Pr(Ri=1Xi),Pr(Ri=1Xi,Yi)=πi.\pi_i =\Pr(R_i=1\mid X_i), \qquad \Pr(R_i=1\mid X_i,Y_i)=\pi_i.

还需要正值性条件 0<πi10<\pi_i\leq1:目标总体中的每类人都要有正概率被观测。若 πi\pi_i 已知,则

E(RiYiπiXi,Yi)=YiπiE(RiXi,Yi)=Yi.\mathbb{E}\left( \left. \frac{R_iY_i}{\pi_i} \right|X_i,Y_i \right) = \frac{Y_i}{\pi_i}\mathbb{E}(R_i\mid X_i,Y_i) =Y_i.

因此,观测概率已知时,Horvitz–Thompson 形式的均值估计量

μ^HT=1ni=1nRiYiπi\widehat\mu_{\mathrm{HT}} =\frac1n\sum_{i=1}^n \frac{R_iY_i}{\pi_i}

用观测概率的倒数补偿选择机会。实际分析通常用模型估计值 π^i\widehat\pi_i 替代 πi\pi_i,也常把估计权重归一化,得到 Hájek 比率形式

μ^H=iRiYi/π^iiRi/π^i.\widehat\mu_{\mathrm{H}} = \frac{\sum_i R_iY_i/\widehat\pi_i} {\sum_i R_i/\widehat\pi_i}.

后者的分母会随样本变化,有限样本通常带有比率估计偏差,却常能改善稳定性。用于组间比较时,应在每个目标组内完成相应加权,或在统一模型中明确组别结构。

若某类人的 πi\pi_i 接近 0,权重会很大,方差随之上升;若 πi=0\pi_i=0,该类人的结局从不出现,任何有限权重都无法恢复其均值。本章的实际规则直接读取末次结局,已经违反上面的 MAR 条件。即使分析者只用组别和第 2 月状态拟合观测概率,治疗组第 2 月不高于 6 的层内也没有完整结局,实际重叠性为零。只依靠最终可见数据,IPW 无法还原 4 与 3。

方法来路:逆概率权重从哪里来

Horvitz 与 Thompson 在 1952 年研究不等概率抽样:入样机会较小的单位需要代表更多总体单位,由此得到按入样概率倒数加权的无偏总量估计。缺失数据把“是否入样”换成“结果是否被观测”,同一条件期望恒等式随之得到逆概率加权方法。

适用边界:极端权重需要事先处理规则

权重截断或稳定化可以降低方差,同时会改变偏差。阈值若在查看结果后才确定,又会增加分析选择。报告中应给出观测概率模型、权重分布、极端权重处理规则及相应敏感性结果。

单次均值填补怎样压缩方差

用治疗组完整病例均值 9.59.5 填补 T3、T4,治疗组数据会变成

10,9,9.5,9.5.10, 9, 9.5, 9.5.

先看一般情形。rr 个可见结果的均值为 Yˉobs\bar Y_{\mathrm{obs}}、样本方差为 sobs2s_{\mathrm{obs}}^2。再用同一均值填入 nrn-r 个空白。填补值相对均值的偏差全为 0,所以离均差平方和保持为

i:Ri=1(YiYˉobs)2=(r1)sobs2.\sum_{i:R_i=1} (Y_i-\bar Y_{\mathrm{obs}})^2 =(r-1)s_{\mathrm{obs}}^2.

若把 nn 个数都当成真实观测,填补后的样本方差成为

sfill2=r1n1sobs2.s_{\mathrm{fill}}^2 =\frac{r-1}{n-1}s_{\mathrm{obs}}^2.

只要 1<r<n1<r<n,乘数 (r1)/(n1)(r-1)/(n-1) 就小于 1。本例 r=2,n=4r=2,n=4,且

sobs2=(109.5)2+(99.5)221=0.5,s_{\mathrm{obs}}^2 =\frac{(10-9.5)^2+(9-9.5)^2}{2-1} =0.5,

因而

sfill2=2141×0.5=16.s_{\mathrm{fill}}^2 =\frac{2-1}{4-1}\times0.5 =\frac16.

若继续把四个值视为四次无误差观测,治疗组均值的朴素标准误会从

0.52=0.5\sqrt{\frac{0.5}{2}}=0.5

下降到

1/64=1240.204.\sqrt{\frac{1/6}{4}} =\frac1{\sqrt{24}} \approx0.204.

均值填补既把填补值放在样本中心,又把预测值计作已知事实。第一步压缩数据离散程度,第二步虚增有效分母;回归系数、相关系数和标准误也会受到类似影响。

多重插补怎样把补法差异放回方差

多重插补用多份随机补全表达缺失值的不确定性。一次完整流程包括:

  1. 用可见结局、组别、早期随访和辅助变量建立插补模型;

  2. 抽取模型参数,再从相应预测分布抽取缺失值;

  3. 对每份完整数据使用同一个预定分析,得到点估计与完整数据方差;

  4. 合并各份结果。

第二步若只重复填入同一组条件均值,插补间差异会消失,参数估计与个体预测的不确定性仍未进入结果。

mm 份完整数据给出

θ^1,,θ^m,U1,,Um.\widehat\theta_1,\ldots,\widehat\theta_m, \qquad U_1,\ldots,U_m.

合并点估计、平均组内方差和插补间方差分别为

θ=1mj=1mθ^j,\overline\theta =\frac1m\sum_{j=1}^m\widehat\theta_j, U=1mj=1mUj,B=1m1j=1m(θ^jθ)2.\overline U =\frac1m\sum_{j=1}^mU_j, \qquad B =\frac1{m-1}\sum_{j=1}^m (\widehat\theta_j-\overline\theta)^2.

全方差公式把总不确定性分成“给定一份完整数据后的分析方差”与“不同合理完整数据带来的估计差异”。Rubin 合并规则据此写成

T=U+(1+1m)B.T =\overline U+ \left(1+\frac1m\right)B.

U\overline U 保留每份数据内部的抽样不确定性,BB 衡量缺失信息造成的补全差异,B/mB/m 则反映只使用有限份插补的额外模拟误差。

用三份插补演示合并算术

假设三份插补把 T3、T4 的第 3 月结果分别补为

(6,5),(5,4),(4,3).(6,5),\qquad(5,4),\qquad(4,3).

对照组样本方差为 sC2=5/3s_C^2=5/3。每份数据都用两独立组均值差及其方差

Uj=sT,j24+sC24U_j=\frac{s_{T,j}^2}{4}+\frac{s_C^2}{4}

进行分析,得到

jj插补值sT,j2s_{T,j}^2θ^j\widehat\theta_jUjU_j
11(6,5)(6,5)17/317/30011/611/6
22(5,4)(5,4)26/326/30.5-0.531/1231/12
33(4,3)(4,3)37/337/31-17/27/2

于是

θ=0.5,U=95362.639,B=0.25.\overline\theta=-0.5, \qquad \overline U=\frac{95}{36}\approx2.639, \qquad B=0.25.

总方差和标准误为

T=9536+(1+13)0.25=107362.972,T1.724.T =\frac{95}{36} +\left(1+\frac13\right)0.25 =\frac{107}{36} \approx2.972, \qquad \sqrt T\approx1.724.

这组三份补全只用于展示合并算术。正式插补要从预先规定、与分析目标相容的预测模型中随机抽取;如此小的样本还需采用多重插补自由度与小样本校正,不能直接套用正态临界值。

方法来路:多重插补要解决单次填补的哪项缺陷

Rubin 在调查无回答研究中发展并系统化多重插补。它保留完整数据分析方法的便利,同时让多份合理补全之间的差异进入方差;1987 年的专著给出了系统理论与应用框架。这个方法传播的是指定插补模型内的不确定性,模型没有覆盖的 MNAR 机制仍需另做敏感性分析。

Δ\Delta 曲线寻找验收线与翻转点

把完整病例均值 9.59.5 作为缺失者均值的参照,再令两个缺失结果共同偏移 Δ\Delta

Ymis(Δ)=9.5+Δ.Y_{\mathrm{mis}}^{(\Delta)}=9.5+\Delta.

两名缺失者占治疗组的一半,对照组没有缺失,所以

θ^(Δ)=2+0.5Δ.\widehat\theta(\Delta) =2+0.5\Delta.

这条直线给出两个不同问题的临界值。报告刚好失去“至少领先 1.51.5 分”的验收资格时,

2+0.5Δ=1.5Δacc=1.2+0.5\Delta=1.5 \quad\Longrightarrow\quad \Delta_{\mathrm{acc}}=-1.

组间差方向翻转时,

2+0.5Δ=0Δ0=4.2+0.5\Delta=0 \quad\Longrightarrow\quad \Delta_0=-4.

后台缺失者真实均值为 3.53.5,相对参照的偏移是 Δtrue=6\Delta_{\mathrm{true}}=-6。代回可得

θ^(6)=2+0.5(6)=1,\widehat\theta(-6) =2+0.5(-6) =-1,

正好恢复完整数据结论。

Δ 敏感性曲线同时显示验收阈值、方向翻转点与后台完整数据。

Δ\Delta 敏感性曲线同时显示验收阈值、方向翻转点与后台完整数据。

更一般地,若治疗组与对照组缺失比例分别为 fT,fCf_T,f_C,两组都采用均值偏移,则简单均值差满足

θ^(ΔT,ΔC)=θ^(0,0)+fTΔTfCΔC.\widehat\theta(\Delta_T,\Delta_C) =\widehat\theta(0,0) +f_T\Delta_T-f_C\Delta_C.

Δ\Delta 是一项无法由观测数据直接估计的敏感性参数。其取值范围应由量表意义、既往研究、退出原因或专家知识约束;若分析方法含非线性链接、时间趋势或交互项,敏感性曲线通常也不再是直线。

把两格空白放回可复核分析链

两个空白把报告效应从 1-1 推到 22。治疗组一半末次结果消失,对照组全部可见;缺失者第 2 月均值为 5.55.5,完整者为 8.58.5;后台末次均值又相差 6 分。均值优势与缺失图案之间存在可检查的对应。

可复核报告至少需要:

  1. 按组别和时点报告计划人数、观测人数、退出时间与缺失原因;

  2. 保存每个结局的 RitR_{it},比较完整者与缺失者的已观测历史;

  3. 预先说明目标参数、主要分析、缺失机制假设和正值性要求;

  4. 同时给出完整病例、加权或插补分析,并说明各自使用的信息;

  5. 报告观测概率、权重分布、插补模型和随机种子;

  6. 给出关键 MNAR 偏移范围、验收阈值与方向翻转点;

  7. 对照原始随访、联系日志、数据版本和分析代码。

这些证据会把“治疗组均值 9.59.5”还原为带选择条件的量 E(YR=1,G=T)\mathbb{E}(Y\mid R=1,G=T),并让两格空白重新进入结论的不确定性。本章选择了进入末次均值的受试者;下一章选择进入回归方程的解释变量,研究比较条件怎样随模型公式改变。

本章知识链

  1. 隐藏一个治疗组结果时,最大报告效应只有 1/61/6;隐藏 4 与 3 后效应达到 22,所以两格是最少改动。

  2. 缺失指示记录可见图案,条件分布 Pr(RY,X)\Pr(\mathbf R\mid\mathbf Y,\mathbf X) 描述生成机制;同一图案可与多种机制相容。

  3. MCAR、MAR 与 MNAR 分别限制缺失概率对完整数据、可见数据和不可见结局的依赖;MAR 还不能单独保证某种分析有效。

  4. 完整病例均值的选择偏差等于缺失比例乘以完整者与缺失者的均值差;本例偏差为 3 分。

  5. 逆概率加权依赖 MAR、正确的观测概率模型与正值性;本例的规则读取缺失结局,低早期状态层内也没有完整结局。

  6. 单次均值填补把样本方差乘以 (r1)/(n1)(r-1)/(n-1);多重插补用 U\overline UBB 分别保留组内和补全之间的不确定性。

  7. 本例 Δ=1\Delta=-1 时失去验收资格,Δ=4\Delta=-4 时效应方向翻转,后台真实偏移 6-6 恢复效应 1-1

思考与练习

  1. 若验收线改为治疗组至少领先 0.50.5 分,求最少需要隐藏几个治疗组结果,并说明最有利的隐藏顺序。

  2. 写出第 2 月和第 3 月各自的观测指示向量,验证本例的缺失图案是否单调;再构造一个非单调图案。

  3. 为同一随访分别设计一个 MCAR、MAR 和 MNAR 机制,写出每个机制允许缺失概率依赖的变量。

  4. 从全期望公式推导完整病例选择偏差恒等式,并用治疗组数据复算 3 分偏差。

  5. 最终可见治疗结果为 10、9 和两个空白。另给出两种与可见数据相容的缺失者均值,计算相应完整治疗均值和组间效应,并说明观测数据为何不能选定其中一种。

  6. 在已知 πi=Pr(Ri=1Xi)\pi_i=\Pr(R_i=1\mid X_i) 且 MAR 成立时,证明 E(RiYi/πi)=E(Yi)\mathbb{E}(R_iY_i/\pi_i)=\mathbb{E}(Y_i);比较 Horvitz–Thompson 形式与归一化 Hájek 形式的分母。

  7. 证明用 rr 个完整病例的均值填补到总数 nn 后,样本方差变为 (r1)sobs2/(n1)(r-1)s_{\mathrm{obs}}^2/(n-1)。讨论这项变化为何没有包含预测误差。

  8. 复算三份插补中的三个治疗组样本方差、三个 UjU_jU\overline UBB 和 Rubin 总方差。

  9. 若敏感性参照效应为 1.21.2、治疗组缺失比例为 0.250.25、对照组无缺失,写出 θ^(Δ)\widehat\theta(\Delta) 并求方向翻转点。

  10. 为一项有四次随访的工程可靠性试验设计缺失审计表,使计划观测、实际观测、退出原因、权重和敏感性参数都能被复核。

专题导航