Ⅱ型区间删失数据的生存分析

山西医科大学卫生统计教研室（030001） 梁 洁 王 彤崔 燕



Ⅱ型区间删失数据的生存分析

山西医科大学卫生统计教研室（030001） 梁 洁 王 彤△崔 燕

问题的提出

例：用两种疗法治疗94名乳腺癌患者随访数月后的结果，其中针对46名乳腺癌患者用放射性疗法（RT）治疗，另48名乳腺癌患者用放射性疗法加辅助性化学疗法（RCT）治疗，每4～6个月随访一次，观察乳腺收缩情况，比较两种疗法的优劣［1］。

表1　两种疗法治疗94名乳腺癌患者随访数月后结果

表1中数据包含了不同患者乳腺收缩的时间信息，但没有观察到精确时间，其中38名患者在随访结束时仍没有观察到乳腺收缩，此为右删失数据；另56名患者的乳腺收缩时间落在不同时间区间内，例如观察时间（6，10］，表示在第6个月随访时患者未出现乳腺收缩，但是在下一次随访，即第10个月时，患者出现了乳腺收缩，乳腺收缩情况出现在第6个月至第10个月两次随访之间，但精确的时间未知，此为Ⅱ型区间删失数据。当数据中包含精确数据，左/右删失数据及Ⅱ型区间删失数据时，称其为混合区间删失数据。

Ⅱ型区间删失主要发生在需要定期随访观察的医疗研究数据中，许多临床试验和纵向研究都存在这种现象，如肿瘤发病率的研究，AIDS的临床医学研究等。由于在实际随访中，有些患者不是按照预先确定的观测时间进行观察，而是选择在较方便的时间进行观察；有些患者可能会错过一个或多个观察后再继续进行观察。因此不同患者的随访观察时间是不同的，研究者只能得到患者在出现某个结局之前最后一次临床观察的时间和出现此结局之后第一次临床观察的时间，即Ⅱ型区间删失数据。如何正确处理此类数据是得出正确结论的关键。Ⅱ型区间删失数据既不同于可以精确测得的数据，又不同于缺失数据，我们要根据其提供的不完整的数据信息，来估计出相对稳定的回归模型参数，从而解决临床实际问题。

处理此类数据的传统方法有：

（1）忽略Ⅱ型区间删失数据，用剩余数据进行生存分析。显然这种方法会丢失数据信息，尤其在Ⅱ型区间删失比例较大的情况下（本例Ⅱ型区间删失数据的比例为59.6%），得到的结果通常是有偏的，不可信的。

（2）Ad hoc：用区间左端点、右端点或中点代替区间删失数据。Flygare（1985）［2］用两参数Weibull回归模型，分别求得参数的最大似然估计值（MLE）以及用区间中点替代删失数据后的参数估计值（MDE），并对两者进行比较，发现在删失比例较大时，MDE可能会造成对参数系统性的错估，MLE优于MDE。Wei （1989）［3］将中点填补法应用于删失比例较小且区间宽度较窄的Ⅱ型区间删失数据，但当删失比例较大、区间宽度较宽且有变异时，此法对经验生存函数估计会产生偏倚，导致估计参数有偏，并低估回归系数的标准误。

解决方法

为了更好地利用Ⅱ型区间删失数据信息，获得更精确的参数估计值，并克服复杂的计算过程，解决含有Ⅱ型区间删失数据的回归分析方法有：参数回归分析，半参数回归分析，非参数最大似然估计以及填补法。

1.参数回归分析

目前对于Ⅱ型区间删失数据的处理方法多为半参数回归分析和非参数回归分析方法，参数回归分析方法研究相对较少。Lindsey（1998）［5］总结了含有混合区间删失数据的参数回归分析，分别在9种不同的分布假设下，对3个存在严重混合区间删失数据的实例进行模拟。发现对严重混合区间删失数据进行参数估计时，改变分布假设后所得估计结果很稳健，并且与相应的非参数模型相比，参数模型包含了更多的数据信息，而且可以直接利用最大似然的理论来实现。在分析Ⅱ型区间删失数据或混合区间删失数据时，常用的参数回归模型包括线性回归模型和Weibull回归模型。

（1）线性回归模型：G′omez（2003）［6］假设自变量为离散型分布，在线性回归模型基础上，用普通最小二乘法与两步条件算法联合估计回归系数和协变量的边际分布。但对于自变量为连续性分布的Ⅱ型区间删失数据，其联合似然方程很复杂，此法不适用。G′omez （2005）［7］用优化技术对其进行改进，但是这种优化技术仍然不适用于自变量连续性分布的Ⅱ型区间删失数据。为解决这个问题，Calle（2005）［8］提出了半参数贝叶斯方法，用Dirichlet混合过程指定模型中除Ⅱ型区间删失协变量以外的所有成分并将其参数化，但是其迭代算法和抽样过程十分复杂耗时。丁邦俊（2012）［9］假设自变量为离散分布的Ⅱ型区间删失数据，应变量为可观测数据，用EM算法将区间删失的问题转化为缺失值问题：E步根据似然方程建立伪似然方程，用缺失值所对应的期望值进行填补；M步最大化伪似然方程得到下一步迭代的参数估计值。对应每一个E步，有相应的M步，经过迭代得到线性回归的参数极大似然估计值（MLE）并证明在一定条件下估计量分布渐近正态分布。Yang（2014）［10］提出在无信息区间删失条件下，结合对数样条的线性回归分析来解决自变量为连续性分布的Ⅱ型区间删失数据。将（li，ui）划分成（K +1）个子区间，用对数样条法构造时间（T）的分布函数，得到对数样条模型的对数似然函数，从而计算参数的极大似然估计值。需要注意的是，存在Ⅱ型区间删失数据时，对数似然函数可能是非凸的，如果用修正的Newton-Raphson法，得到的结果可能是局部最大值而不是全局最大值，且当删失区间长度固定且较宽时，会出现找不到极值的情况或遇到算法不收敛问题（由于样条光滑法适用于密集数据，区间较宽，会导致数据变得非常稀疏）。此时求出T的条件期望，用条件期望替代未知T，算出线性回归模型的最小二乘估计值。将此法与中点替代法相比，当删失比例较大时，由中点替代法得到的估计有偏，且偏倚随着删失区间的增宽而增加，而对数样条法可得到更加精确的估计值和较小的误差均方；与半参数贝叶斯方法相比，两者结果相近，但运算速度比其快100倍。当删失区间宽度固定且较小时，对数样条法表现最好，但当区间宽度较大时，三种方法的表现都不好，特别是对数样条法，其误差均方相当大，在程序运行过程中也会遇到算法不收敛的问题。

（2）Weibull回归模型：Pradhan（2014）［11］假设区间删失机制是无信息发生的，即删失信息与生存事件的发生相互独立，生存时间为Weibull分布时，计算在伽马先验分布下未知参数的贝叶斯估计。对于Weibull分布，当形状参数α已知，此时尺度参数β共轭先验分布为Gamma分布。但当形状参数也未知时，即使存在连续-离散的联合先验分布（连续指尺度参数，离散指形状参数），Weibull分布也没有连续共轭先验；且由于Gamma共轭先验不好解释，所以此法应用有很大的局限性。故当两参数都未知时，假设两参数先验分布为相互独立的Gamma分布，用重要抽样技术来实现贝叶斯估计。但贝叶斯方法计算复杂，且Pradhan没有对贝叶斯估计与EM计算的最大似然估计进行比较，也没有考虑协变量对贝叶斯估计的影响。

2.半参数回归分析

在处理Ⅱ型区间删失数据时，半参数回归方法比较常用，由于其同时考虑了发生结局的时间和影响结局的协变量，而且在计算上又优于复杂的非参数最大似然估计，故比较常用。

（1）Cox比例风险（Cox PH）模型：Satten （1996）［12］提出Cox比例风险模型，用Monte-Carlo结合EM算法来估计回归方程，此法需要进行大量计算。Goggins（1998）［13］在Satten的基础上，将Monte-Carlo模拟放在EM算法的E步，从而避免大量的计算。Goetghebeur（2000）［14］将EM算法与Breslow估计相结合：M步用最大化标准Cox偏似然来估计回归系数，然后用Breslow估计基线风险；E步估计风险集大小，以及在每一个可能的时间集中事件发生的个数。E步同Turnbull的自相合算法，但计算更简单也不需要迭代。Betensky（2002）［15］发展了局部似然估计方法，平滑基线风险函数，其计算比较复杂。Cai （2003）［16］假设一个分段线性基线风险，用惩罚样条方法估计限制最大似然估计。Heller（2011）［17］在估计方程的基础上，用一个反概率权重来选择明确顺序的事件时间对，并将其应用于结肠癌数据分析中。此法不需要估计基线风险，用标准的结构方程技术产生估计值及其近似分布，如果恰当定义选择概率模型和比例风险模型，就能建立无偏估计方程，且计算简单，适用于删失比例较大的临床数据。Sun（2013）［18］基于修正的局部似然产生两种估计过程：第一种估计过程是以区间的左端点填补Ⅱ型区间删失数据，转化为右删失数据来处理，方法简单易行，但是只利用了区间左端点的信息；第二种估计过程利用了两个区间端点，构建出常见的局部似然方程。两种方法都在估计方程理论基础上进行估计，避免了求基线累积风险函数，经模拟研究两种方法都有效可行。比较而言，在删失区间较宽的情况下，第二种估计过程可能会更稳定，效果更好。

（2）半参数AFT，PO，AH模型：由于Cox PH模型并不能很好地拟合所有生存数据，对于某些生存数据要考虑线性转换模型。Odell（1992）［19］将基于Weibull的加速失效时间（AFT）回归模型，应用于Framingham的冠心病研究数据中，计算模型参数的最大似然估计值并与中点替代法进行比较，发现前者的估计结果更好。Murphy（1997）［20］提出应用比例优势（PO）模型，用伪似然方法估计参数。Shen（1998）［21］用sieve极大似然估计和单调样条方法对PO模型的基线风险和回归系数进行估计。Rabinowitz（2000）［22］基于条件logistic回归拟合PO模型，不需要估计无限维度的冗余参数，风险集越小，其估计效果越好。Zeng（2006）［23］提出了用可加风险（AH）模型。Zhu （2008）［24］提出将Ⅱ型区间删失数据通过某种转换方式转换为Ⅰ型区间删失数据，进而用处理Ⅰ型区间删失数据的方法构建AH模型。Wang（2010）［25］在AH模型条件下，用一些基于估计方程的方法来处理Ⅱ型区间删失数据，这种方法比较稳健，可以应用于有信息和无先验信息的区间删失数据，且不用估计基线方程，实现过程简单快速。

（3）其他模型：Lin（2010）［26］提出了半参数probit模型，作为一种替代PH，PO，AH，AFT的模型，从频率论和贝叶斯的观点，进行Ⅱ型区间删失数据回归系数及基线生存函数的平滑估计，简单易行。此法基于似然且不需要对观察过程进行假设，所以对于任何区间删失数据都适用，可以进行拓展。Shao（2013）［27］用二分模型来描述这种含有复杂数据的生存过程，用灵活的变系数模型拟合Ⅱ型区间删失数据。通过交叉验证法选择带宽，用局部多项式方法来拟合模型中未知方程，讨论估计值的稳定性和渐近性分布，并提出用bootstrap做推断，构建BIC类的模型选择方法，合理解释了模型中的参数和非参数部分，易于理解，计算和推断易于完成。

3.非参数最大似然估计

在多数情况下，含Ⅱ型区间删失数据的随机变量总体分布未知，所以非参数回归分析方法的研究和应用较多。Ⅱ型区间删失数据生存函数的非参数最大似然估计没有闭集解，常用的获得非参数最大似然估计值的方法有：

Turnbull算法：Turnbull（1976）［28］提出自相合算法，结合EM算法得到单调收敛的经验分布函数最大似然估计值，与Newton-Raphson算法相比，此法运算简单，容易实现，便于理解，所以Turnbull的自相合算法在处理Ⅱ型区间删失数据时，常被采用。Gentlemen （1994）［29］应用凸优化技术验证自相合算法得到的最大似然估计值是唯一的。Shen（2014）［30］在应用Turnbull自相合算法的前提下，考虑了自变量的权重，获得条件分布函数的非参数估计值，在带宽固定的情况下，也容易实现。但是自相合算法在删失比例较大的情况下，估计值偏倚会增大。

其他求非参数最大似然估计值的算法：Wang （2008）［31］提出Ⅱ型区间删失数据降维的非参数最大似然算法，用CNM（constrained Newton method），SBN （subspace-based Newton），及ICM（iterative convex minorant algorithm）算法有效地剔除迭代过程中出现的冗余支撑区间，找到梯度最大的支撑区间，从而有效提高计算效率。用迭代的方法估计非参数最大似然估计值的缺点是难以得到估计值的渐近分布。Yuan （2014）［32］用保序回归的方法得到一个近似的非参数最大似然估计值，并对其渐近分布和收敛速度做出评价。

上述几种算法中自相合算法简单易懂容易完成，对于小/中等样本数据，或右删失比例较大的数据，此法是很好的选择，但目前还没有方法来证明其估计值即为非参数最大似然估计值；ICM算法迭代步骤较少，计算所用时间短，收敛速度快，而且它们的整体收敛性可以被证明；保序回归算法可以得到估计值的渐近分布，但计算较复杂。

4.填补法

在用非参数回归分析方法分析Ⅱ型区间删失数据时，除了用以上算法求含有Ⅱ型区间删失数据变量的最大似然估计值之外，还可以用填补法将Ⅱ型区间删失数据转换成精确数据或者Ⅰ型区间删失数据，进而用成熟的方法求解其非参数最大似然估计值。Pan（2000）［33］用两种多重填补法：PMDA（poor man｀s data augmentation）与ANDA（asymptotic normal data augmentation），将Ⅱ型区间删失数据问题转化为右删失数据问题，容易实现，与非参数最大似然估计进行比较，在小样本和中等样本情况下，多重填补法得到的结果更好。当删失比例较大时，推荐ANDA填补法。Xiao （2012）［34］利用多重填补的思想，在Weibull分布条件下用迭代单点填补法——分位数填补法处理Ⅱ型区间删失数据，即用删失区间的条件分位数填补删失数据，估计结果稳定，方法简单易实现，但在删失比例较大时估计可能会有偏。Deng（2014）［35］用无偏转换的思想，将Ⅱ型区间删失数据进行无偏转换，转换后的值与真值有相同的条件期望，得到回归函数的最邻近估计，其有很强的稳定性和渐近正态性。Han（2013）［36］建立删失数据的伪观察值POS（pseudo-observations）［37］，通过PH或PO等适当的连接函数，将问题转化为广义线性回归模型，这种方法计算比较简单，不需要建立许多假设，可考虑将此法应用于双边区间删失［38］及其他类型的区间删失数据分析。

结 论

生存分析的首要任务是估计生存函数。估计的生存函数可用来：评价生存变量的模型假设是否正确；估计生存概率；比较不同组之间的生存情况；预测患者未来的生存时间等等。参数模型有光滑数据的内在性质，又可从邻近点获取数据信息，与非参数模型相比，Ⅱ型区间删失对于参数模型的影响要小，因此如果有足够的先验信息提示数据可用参数模型，那么参数估计要比非参数估计更简单高效。然而在实际研究中，常常没有足够的先验信息提示哪一种参数模型假设是合理的，因而半参数回归分析和非参数回归分析被更多地应用于分析含Ⅱ型区间删失数据的实例中。而对于生存函数的半参数和非参数估计，与右删失数据的处理方法相比，无论在理论还是实践中，对Ⅱ型区间删失数据的处理方法要更复杂。

在半参数回归分析中，有时既要估计回归系数，也要估计基线生存函数。但对于Ⅱ型区间删失数据而言，估计基线生存函数或风险函数的方法和计算是比较难的，由于其复杂的数据结构以及与样本量近似的参数个数特点，寻找一个合适的理论依据也很难。有些方法不估计基线生存函数，直接估计回归参数，比如在PH模型下的边际似然方法、Monte Carlo EM算法，但这些方法都是在假设事件发生时间的发生次序与观察到的数据一致的前提下估计的，而且计算都比较复杂。

非参数回归分析常用的Turnbull自相合算法，只考察了自变量相互独立，删失区间宽度固定的情况，没有考虑删失区间较宽或有重叠或区间宽窄不一的情况，而且删失比例的增大对结果的精度有影响。

由于上述方法计算比较复杂，或需要许多假设，而这些假设在实际应用中很难满足，所以目前在进行含Ⅱ型区间删失数据的生存分析时，一方面可以考虑用其他模型来估计生存函数，比如半参数probit模型，变系数模型，两者估计过程中假设较少，而且应用比较广泛；另一方面考虑填补法，将Ⅱ型区间删失数据的问题转换成精确数据或Ⅰ型区间删失数据。但是简单地用区间中点或端点进行填补会造成结果的偏倚，为了在转换过程中更好地利用区间信息，得到精确有效的估计结果，可以使用Deng的无偏转换填补思想和POS方法，这两种方法更好地利用了区间信息，容易实现且估计结果更加稳定。

在实际应用中，对于删失比例较小，删失区间宽度适中、变异较小的Ⅱ型区间删失数据，上述许多方法可得到较为满意的估计结果，但是对于删失比例较大，删失区间较宽，变异较大，协变量个数较多的Ⅱ型区间删失数据，如何从众多方法中做出较好的选择并进行成功的应用还需进一步研究讨论。

参考文献

［1］Sun J.The Statistical Analysis of Interval-censored Failure Time Data.American：springer，2006.

［2］FlygareME.Maximum Likelihood Estimation For The 2-Parameter Weibull Distribution Based On Interval-Data.Ieee Transations on Reliability，1985，R-34，1.

［3］Wei LJ，Lin DY，Weissfeld L.Regression analysis of multivariate incomplete failure time data by modeling marginal distributions.Journal of the American Statistical Association，1989，84：1065-1073.

［4］Goggins WB.A Proportional Hazards Model for MultivariateInterval-Censored Failure Time Data.Biometrics，2000，5（6）：940-943.

［5］Lindesey JK.A study of interval censoring in parametric regressionmodels.Lifetime Data Analysis，1998，4：329-354.

［6］G′omez G.Inference for a linear regression model with anintervalcensored covariate.Statist Med，2003，22：409-425.

［7］G′omez G，Calle ML，Oller R.Frequentist and Bayesian approaches for interval-censored data.Statistics Papers，2005，2：139-173.

［8］Calle ML，G′omez G.A semiparametric hierarchical method for a regression model with an interval-censored covariate.Australian &New Zealand Journal of Statistics，2005，47：351-364.

［9］Ding BJ.Regression model with an interval-censored data covariant.Chinese Journal of Applied Probability and Statistics，2012，28（2）：.

［10］Yang Y.A Logspline Estimation for a Linear Regression Model withan Interval-Censored Continuous Covariate.Taylor & Francis，2014，43（10）：2521-2539.

［11］Pradhan B.Analysis of Interval-Censored Data with Weibull LifetimeDistribution.Sankhya：The Indian Journal of Statistics，2014，76-B（1）：120-139.

［12］Satten GA.Rank-based inference in the proportional hazards model for interval-censored data.Biometrika，1996，83：355-370.

［13］Goggins WB.A Markov chain Monte Carlo EM algorithm for ana-lyzinginterval-censored data under the Cox proportional hazards model.Biometrics，1998，54（4）：1498-1507.

［14］Goetghebeur E.Semiparametric Regression Analysis of Interval-CensoredData.Biometrics，2000，5（6）：1139-1144.

［15］Betensky R，JCLL，Wand AMP.A local likelihood proportional hazards model for intervalcensored data.Statistics in medicine，2002，（21）：263-275.

［16］Cai T，Betensky RA.Hazard regression for interval-censored data with penalized spline.Biometrics，2003，59：570-579.

［17］Glenn Heller.Proportional hazards regression with interval censored data using an inverse probability weight.Lifetime Data Anal，2011，17：373-385.

［18］Sun J.Simple estimation procedures for regression analysis of interval-censored failure time data under the proportional hazards model.Lifetime Data Anal，2015，21（1）：138-55.

［19］Odell PM.Maximum Likelihood Estimation for Interval-Censored D.Biometrics，1992，48（3）：951-959.

［20］Murphy SA，Rossini AJ，Van der Vaart AW.Maximum likelihood estimation in the proportional odds model.Journal of the American Statistical Association，1997，92：968-976.

［21］Shen X.Proportional odds regression and sieve maximum likelihoodestimation.Biometrika Trust，1998，85（1）：165-177.

［22］Daniel Rabinowitz RABA.Using Conditional logistic Regression to Fit ProportionalOdds Models to Interval Censored Data.Biometrics，2000，56：511-518.

［23］Zeng D.Semiparametric additive risks model for interval-censoreddata.Statistica Sinica，2006，16：287-302.

［24］Zhu L.A transformation approach for the analysis of interval-censored failure time data1.Lifetime Data Anal，2008，14：167-178.

［25］Wang L.Regression analysis of case II interval-censored failuretime data with the additive hazards model.Statistica Sinica，2010，20：1709-1723.

［26］Lin X.A semiparametric probit model for case 2 interval-censoredfailure time data.Statistics in Medicine，2010，DOI：10.1002/ sim.3832.

［27］Shao F.Semiparametric varying-coefficient model for intervalcensored data with a cured proportion.Statistics in Medicine，2013，DOI：10.1002/ sim.6054.

［28］Turnbull BW.The Empirical Distribution Function with Arbitrarily-Grouped，Censored and Truncated Data.Journal of the Royal Statistical Society，1976，38（3）：290-295.

［29］Gentleman R，Geyer CJ.Maximum likelihood for interval censored data：Consistency and computation.Biometrika，1994，81：618-623.

［30］Wang Y.Dimension-reduced nonparametric maximum likelihood computation for interval-censored data.Computational Statistics & Data Analysis，2008，52：2388-2402.

［31］Shen P.A Generalization of Turnbull′s Estimator for Interval-Censored and Doubly Truncated Data.Taylor & Francis，2014，43：14，2958-2972.

［32］Yuan A.Approximate Nonparametric Maximum Likelihood Estimationfor Interval Censoring Model Case II.International Journal of Statistics and Probability，2014，3（3）.

［33］Pan W.A Multiple Imputation Approach to Cox Regression with-Interval-Censored Data.Biometrics，2000（56）：199-203.

［34］Xiao X.Study of an imputation algorithm for the analysis of intervalcensored data.Taylor &Francis，2012，84（3）：477-490.

［35］Deng WL.Nonparametric regression with interval-censored data.Acta Mathematica Sinica，English Series，2014，30（8）：1422-1434.

［36］Seungbong Han.A Semiparametric Regression Method for Interval.Communications in Statistics-Simulation and Computation，2013，1 （43）：18-30.

［37］Andersen PK.Pseudo-observations in survival analysis.Statistical Methods in Medical Research，2009，0（0）：1-29.

［38］Seungbong Han.A flexible semiparametric modeling approach for doublycensored data with an application to prostate cancer.Statistical Methods in Medical Research，2013，0（0）：1-18.

（责任编辑：郭海强）

通信作者：△王彤，E-mail：wtstat1@ sina.com