定量纵向数据缺失值处理方法的模拟比较研究

陈丽嫦 衡明莉 王 骏 陈平雁△

【提 要】 目的 比较末次观测结转法(LOCF)、重复测量的混合效应模型法(MMRM)、多重填补法(MI)在处理纵向缺失数据中的统计性能。方法 以双臂设计、4次访视、3种访视间相关程度为应用背景，采用Monte Carlo模拟技术，产生模拟完整纵向数据后考虑两种缺失比例和三种缺失机制，即完全随机缺失(MCAR)、随机缺失(MAR)和非随机缺失(MNAR)的缺失数据集。以完整纵向数据的分析结果为基准，评价不同处理方法的统计性能，包括Ⅰ类错误、检验效能、组间疗效差的估计误差及其95%置信区间(95%CI)宽度。结果 所有情况下，MMRM和MI均可控制Ⅰ类错误，检验效能略低于完整数据；LOCF大多难以控制Ⅰ类错误，检验效能变异较大。多数情况下MMRM和MI的点估计误差较低，LOCF则表现不稳定。所有情况下，MI的95%CI最宽，MMRM次之，LOCF最窄。结论 MCAR和MAR缺失机制下，MMRM与MI的统计性能相当，受各种因素影响较有规律，可根据实际情况选择其中一个作为主要分析。LOCF因填补方法的特殊性使得变异较小，精度较高，但其最大的缺陷是不够稳健且不能有效控制I类错误，需谨慎使用。基于MNAR缺失机制对缺失数据进行敏感性分析以考察试验结果的稳健性是必要的。

纵向数据在临床试验中颇为常见，如对受试者做多个访视点的重复观察记录。纵向数据如果存在缺失值会导致分析结果产生潜在偏倚[1-3]。目前，有关定量纵向数据缺失值的处理方法有多种，但哪种方法更具优良特性尚无定论。本研究针对定量纵向数据缺失资料，采用模拟研究方法，考虑完全随机缺失(missing completely at random，MCAR)、随机缺失(missing at random，MAR)和非随机缺失(missing not at random，MNAR)三种不同缺失机制，比较常用的四种缺失数据处理方法，即末次观测结转法(last observation carried forward，LOCF)、重复测量的混合效应模型法(mixed model for repeated measurements，MMRM)、基于马尔可夫链蒙特卡罗(Markov chain Monte Carlo，MCMC)的多重填补法和基于线性回归的多重填补法在纵向缺失数据方面的统计性能。

原理与方法

1.末次观测结转法(LOCF)

LOCF属于单一填补法，是指将最近一次的观察数据填补缺失值。在所有缺失值填补完成后，形成一个完整的数据集，再按照既定的分析方法进行分析。

2.重复测量的混合效应模型(MMRM)

该方法为一般混合效应回归模型(general mixed-effects regression model,MRM)的一种特殊形式，由Mallinckrodt 等人在2001年定义[4]。MRM模型如下：

Yi=Xiβ+Ziνi+εi

(1)

其中，Yi为第i个受试者ni×1维反应向量；Xi为ni×p维已知固定效应设计矩阵；β为p×1维未知固定效应参数；Zi为ni×r维已知随机效应设计矩阵；νi为r×1维随机效应参数，服从N(0,∑v)分布，εi为ni×1维随机误差，服从N(0,∑εi)分布。νi和εi相互独立。故式(1)中，Yi服从均数为Xiβ，方差协方差矩阵为Zi∑vZ′i+∑εi的多元正态分布。

MMRM模型将组别、访视时间以及二者的交互作为固定因素，受试者内误差作为随机效应。有研究表明，无论方差协方差矩阵的真实情况如何，受试者内误差采用非结构化(unstructured,UN)的协方差矩阵可以控制Ⅰ类错误[5-6]。

3.多重填补

由Rubin[7]提出的多重填补法旨在解决调查研究中无响应的情况，也适用于处理临床研究的缺失值。该方法通过对每个缺失值填补多次，形成多个完整的数据集，对每个完整数据集按既定的分析方法处理后再使用统一的方法[8]进行合并得出综合的结论。

(1)马尔可夫链蒙特卡罗(MCMC)[9]

(2)线性回归法

根据已有观测数据，建立缺失值与协变量的回归方程，基于此方程，从参数的后验预测分布模拟出新的方程用于缺失值的填补。假设Yj是一个含有缺失值的连续性变量，建立回归方程

(2)

模拟研究

1.完整数据模拟

根据一项治疗黄斑水肿，以最后一个访视点较基线的最佳矫正视力变化值为主要评价指标的双臂阳性对照临床试验结果，设置每组样本量为160例，包括基线在内有4个访视点。假设观测数据服从多元正态分布。设置四种疗效变化模式，所有疗效变化模式中，对照组4个访视点均数为(57,57,60,62)。E0：两组各时间点的总体均数相同；E00：两组总体均数在第1和第4次访视相同，在第2和第3次访视不同，试验组4个访视点均数为(57,58,61,62)；E1：总体均数两组基线相同，试验组后3次访视比对照组大，试验组4个访视点均数为(57,62,63,66)；E11：总体均数两组基线相同，试验组在第2和第3次访视稍大于对照组，在第4个访视点较对照组有更大的提高，试验组4个访视点均数为(57,58,61,66)。上述4种疗效变化模式中两组的方差均相同，4个访视点的方差为(100,110,130,130)；其相应的相关阵如下：

其中C1、C2和C3分别代表各访视点评价指标具有低、中和高相关性。运用SAS 9.4对以上12种组合各模拟生成2000个完整数据集。

2.缺失数据的构造

根据三种缺失机制(MCAR、MAR、MNAR)，通过随机删除完整数据中的部分数据构造缺失数据。设定每个受试者均有基线后第一次测量，即仅可能在第3次和/或第4次访视出现缺失值。缺失模式为单调缺失，即本次缺失后，往后访视的数据一并缺失。

MCAR缺失机制下，假设每个访视点的缺失独立服从概率为p的二项分布。

MAR缺失机制下，假设某一访视数据的缺失与缺失前一次的观测结果有关，即

logit(p(yi=missing|y0,y1,…,yi-1))=a+byi-1

(3)

MNAR缺失机制下，假设某一访视数据的缺失与本次测试结果有关，即

logit(p(yi=missing|y0,y1,…,yi-1))=a+byi

(4)

设定的缺失机制参数见表1，试验组约有10%～26%的缺失比例，对照组约有20%～26%的缺失比例，试验组缺失比例均比对照组低。若两组缺失机制参数设定相同，则缺失比例接近，相差约为0%～5%；若两组缺失机制参数设定不同，两组相差10%～15%。

表1 缺失机制参数设定*

*：p代表二项分布的概率，a和b代表函数(3)、(4)的截距和斜率参数。

综上，每个完整数据集均对应构造18种缺失数据集，即三种缺失机制(用M1、M2、M3表示)、三种相关阵(C1、C2、C3)和两种缺失比例(用D0和D1分别表示相近和相差较大)。

3.缺失值处理方法

分别采用LOCF、MMRM、MCMC多重填补和回归多重填补，对构造的缺失数据进行分析处理。LOCF填补后，将基线观测值作为协变量进行协方差分析，比较两组主要评价指标是否存在统计学差异。MMRM模型纳入基线观测值、组别、访视点、访视点和组别的交互作用作为固定效应，受试者内误差作为随机效应(采用非结构协方差矩阵)。MCMC和回归法均进行5次填补，对填补后的数据进行协方差分析。本次研究还将对完整数据进行协方差分析，以作为各方法比较的基准。

4.评价指标

评价指标包括Ⅰ类错误、检验效能，组间疗效的估计误差及其95%置信区间宽度。

结 果

1.Ⅰ类错误和检验效能

如图1所示，在所有设定情况中，MMRM以及多重填补法(MCMC、回归法)均可控制Ⅰ类错误且变化平稳。LOCF表现不稳定，多数情况下难以控制Ⅰ类错误。

在所有设定情况中，MMRM的检验效能略高于多重填补法(MCMC、回归法)，但MMRM以及多重填补法的检验效能均低于完整数据，其变化趋势根据不同的相关系数矩阵和缺失比例情况与完整数据的变化基本相同。随着相关系数的增加，检验效能增大，两组缺失比例差异的扩大将降低检验效能。LOCF在E1模式下的所有设定情况下，检验效能略高于完整数据，E11模式下表现不稳定，近一半情况下与完整数据相近，另一半情况下与MMRM和多重填补法相当。

图1 四种处理方法各种情况下的Ⅰ类错误及检验效能**：M1、M2、M3分别代表完全随机缺失、随机缺失、非随机缺失；C1、C2、C3分别代表三种相关系数矩阵；D0代表两组缺失机制参数设置相同，缺失比例接近，D1表示两组缺失机制参数设置不同，缺失比例相差较大(下同)。

2.组间疗效差的估计误差

组间疗效差的估计误差如图2所示，在MCAR和MAR缺失机制下，多重填补法(MCMC、回归法)的估计误差较小，在0附近波动；MMRM与多重填补法相近，但在部分MAR缺失机制下估计误差稍大，会低估组间疗效。疗效变化模式、访视点间相关性和两组缺失比例对MMRM和多重填补法几乎没有影响。LOCF在MCAR和MAR缺失机制下的估计误差较大且不稳定，受疗效变化模式的影响较大，在E11疗效变化模式会高估疗效，其余模式下低估疗效；多数情况下，各访视点之间的相关性减小和/或两组缺失比例差距的增加会增大估计误差。

图2 四种处理方法各种情况下的组间疗效估计误差

在MNAR缺失机制下，MMRM和多重填补法表现一致，估计误差随着各访视点相关性的增加以及两组缺失比例的接近而减少。如果两组的缺失比例接近，多数情况下，MMRM和多重填补法的估计误差最小。疗效变化模式对MMRM和多重填补法几乎没有影响。LOCF受到疗效变化模式的影响较大，在E11疗效变化模式下出现高估组间差异以及两组缺失比例的增大反而降低组间疗效估计误差的情况。访视点间相关性对其影响不大。

3.组间疗效差估计的95%置信区间宽度

组间疗效差估计的95%置信区间宽度如图3所示，所有情况下，多重填补法的95%置信区间宽度最大，MMRM次之但与其相近，LOCF最窄。95%置信区间宽度均随着相关系数的增强而变窄，缺失机制对其影响不大。MMRM和多重填补法两种方法中，缺失比例差异的增加降低了两组疗效差值的95%置信区间宽度。

图3 四种处理方法各种情况下的组间疗效的95%置信区间宽度

讨 论

本文共设定四种疗效变化模式，每种变化模式下根据三种缺失机制、三种各访视点疗效相关系数和两种两组缺失比例情况设定18种缺失数据集，对每种情形的三种缺失机制分别采用四种缺失值处理方法(LOCF、MMRM、MCMC、多重填补的回归法)进行处理分析。

LOCF法简单、容易理解，但多数情况下，Ⅰ类错误难以控制，检验效能和估计误差表现不稳定，MAR和MNAR缺失机制增加其不稳定性。访视点之间相关系数的变化和两组缺失比例变化相较于MMRM和多重填补法没有固定规律的影响。这可能因为LOCF受到疗效变化模式影响更大，疗效变化越不稳定，估计误差越大(因LOCF假设缺失值的填补值为最后一次观测值的概率为100%，该假设同时降低了估计的变异)。多项研究也表明，LOCF方法不够稳健，降低了估计的变异，并不总是保守的[10-13]。在使用该方法时，需注意其前提假设的合理性，谨慎使用其作为主要分析。

MMRM的处理方法无需对缺失数据进行填补，纳入所有观测的数据建模进行分析，符合意向性原则，本方法在MCAR和MAR机制下各项统计性能优异、稳定，在Ⅰ类错误、检验效能及置信区间宽度上表现优于多重填补，亦有研究表明该方法统计性能优于多重填补[14]。该方法对缺失机制的假设为MAR，相关系数越大和两组缺失比例差异越小，估计误差越小；相关系数越大和两组缺失比例差异越大，置信区间宽度越窄；疗效变化模式对其影响较小。根据模拟的结果，仍需注意以下两点：(1) 由于MMRM中纳入了各种固定效应、组别和访视的交互效应，组别效应的统计检验并不与所关注的最后一个访视点的点估计及其置信区间完全一致，反而受疗效变化模式影响较大。(2) MNAR缺失机制下，估计误差将增加，应使用其他方法对偏离MAR缺失机制假设的情况进行敏感性分析。

多重填补法对一个缺失数据填补多次，相较于其他方法考虑了填补数据的变异，缺失机制假设为MAR。本研究考察的MCMC和回归法统计性能相近。考虑多重填补的方法置信区间最宽，变异程度的增加可能使得其在MAR缺失机制下，估计误差略小于MMRM。有研究指出多重填补法高估变异程度[15]，填补和分析之间存在冲突[16]。在本研究中，多重填补法和MMRM性能相当，可根据实际情况选择其中一个作为主要分析，另一个作为敏感性分析。需要注意缺失机制对各种方法的影响，建议采用基于MNAR缺失机制下的其他分析方法，探索试验结果的稳健性。

本研究虽然探索了四种方法在处理合计72种情形的效果，但疗效变化模式、相关系数矩阵、各组缺失比例等参数设定并不能涵盖所有的可能组合。本研究也未考虑具体的比较类型(如非劣和等效性)，故本研究结论具有一定的局限性。