文章信息
- 翟映红, 陈琪, 韩贺东, 赵欣欣, 高永晴, 周详, 贺佳.
- Zhai Yinghong, Chen Qi, Han Hedong, Zhao Xinxin, Gao Yongqing, Zhou Xiang, He Jia.
- 联合模型介绍及在医学研究中的应用
- Introduction of joint model and its applications in medical research
- 中华流行病学杂志, 2019, 40(11): 1456-1460
- Chinese Journal of Epidemiology, 2019, 40(11): 1456-1460
- http://dx.doi.org/10.3760/cma.j.issn.0254-6450.2019.11.021
-
文章历史
收稿日期: 2019-05-24
2. 第二军医大学卫生统计学教研室, 上海 200433;
3. 安徽医科大学卫生管理学院, 合肥 230032
2. Department of Health Statistics, The Second Military Medical University, Shanghai 200433, China;
3. Department of Health Management, Anhui Medical Universite, Hefei 230032, China
1.背景:在医学随访研究中,针对每位受试者通常需要收集不同类型的测量数据,包括重复测量的纵向数据以及感兴趣事件的发生时间。典型的例子是涉及艾滋病的临床研究,其中生物标志物CD4+T淋巴细胞(CD4)计数往往有规则地以预定时间间隔测量。以Abrams等[1]的随机临床试验为例,主要研究目的是比较两种替代抗反转录病毒药物去羟肌苷和扎西他滨的有效性和安全性,研究中收集的数据包括纵向(研究开始以及之后2、6、12和18个月的CD4计数)和生存(患者出现死亡事件的时间)2类。这些数据的特点:①纵向测量数据是间歇性收集,并受测量误差影响;②感兴趣事件的发生将终止纵向测量过程,可能导致有用信息缺失;③纵向测量过程会影响生存过程。这些属性意味着分开进行重复测量和生存结果的分析是低效的,因为未充分利用重复测量过程和生存过程之间的潜在关系,并且忽略了测量误差,会导致两者之间联系的有偏估计。联合模型(Joint model)的出现,为同时分析这2类数据提供可能。1997年Wulfsohn和Tsiatis[2]首先提出了用于处理重复测量和生存结局数据的联合模型,联合模型最早是应用于艾滋病研究中[3]。近年来,联合模型得到进一步的发展,形式有更多扩展[4],已被广泛用于其他临床研究领域,包括癌症[5]、心血管疾病[6]和肾移植研究[7]等。
2.联合模型的基本结构:联合模型主要有2类模型框架:随机效应联合模型(random e ects joint model)[8]和基于G估计的边际结构联合模型(marginal structural joint models)[9]。随机效应联合模型是最常见的,本文主要对其展开介绍。
随机效应联合模型由纵向子模型和生存子模型组成。Ti表示第i个受试者被观察到的事件发生时间,是真实事件时间Ti*和潜在删失时间Ci的最小值,δi是事件发生指标。yi(t)表示第i个受试者在时间点t被观测到的纵向结果值,真实的未观察到的纵向结果值用mi(t)表示(受到测量误差的影响)。
(1)纵向子模型:对于纵向子模型,常用的是线性混合效应模型:

其中β是固定效应的系数,bi是随机效应的系数,xi(t)和zi(t)分别表示固定和随机效应设计矩阵的行向量。εi(t)是测量误差项,假设与bi无关,服从均值为0,方差为σ2的正态分布。线性混合效应模型通过假设观察到的纵向结果的水平yi(t)等于真实水平mi(t)加上随机误差项来解释测量误差问题。
(2)生存子模型:对于生存子模型,常采用Cox比例风险模型:

Mi(t)={mi(s),0≤s≤t}表示直至时间点t的真实未观察到的纵向过程轨迹,h0(·)表示基线风险函数,ωi是基线协变量(如治疗指标,疾病史等)的向量,γ是回归系数。纵向子模型和生存子模型通过参数α连接,α量化了潜在纵向结果对事件风险的影响。在标准Cox模型中基线风险函数h0(·)通常未指定,联合模型中需要明确地定义h0(·)。可以使用已知参数分布的风险函数(如Weibull、对数正态、Gompertz和Gamma分布),或更灵活的方法比如分段常数(piecewise- constant approaches)和回归样条(regression splines approaches)[10]。
3.联合模型的参数估计:联合模型参数估计的方法目前主要是极大似然方法和使用马尔科夫链-蒙特卡罗的贝叶斯方法,本文仅介绍联合模型中使用最多的极大似然估计[11]。联合模型的极大似然估计是最大化事件发生时间和纵向结果联合分布{Ti,δi,γi}的对数似然函数。为了定义这种联合分布,将假设时间无关随机效应的向量bi是纵向和生存过程的基础。这意味着这些随机效应既解释了纵向和事件结果之间的关联,也解释了纵向过程中重复测量之间的相关性。


其中θ表示参数矢量。由于生存子模型和纵向子模型具有相同的随机效应,因此这种类型的联合模型称为共享参数模型(shared parameter models)。据前提出的建模假设,以及这些条件独立假设,第i个患者的联合对数似然贡献可以表示为:

生存模型的似然函数可以表示为:


f{γi(tij)|bi;θ}是纵向响应的单变量正态密度函数,g(bi;θ)是随机效应的多元正态密度函数。高斯积分[12]和蒙特卡罗[13]可用于参数的求解。
4.实例应用:通过实例来演示联合模型在R中的实现过程,并将联合模型的结果与扩展的Cox模型结果进行对比。
(1)数据来源:试验数据来自Abrams等[1]随机临床试验(数据集AIDS可从R软件JM程序包中获得),研究人群为467例在抗反转录病毒治疗期间接受齐多夫定治疗失败或对齐多夫定治疗不耐受的晚期AIDS患者,患者被随机分配接受去羟肌苷或者扎西他滨治疗。记录研究开始以及之后2、6、12和18个月的CD4计数和患者出现死亡事件的时间。CD4计数随时间的减少表明患者免疫系统状况恶化,往往标示较差的预后。在下面的例子中,我们关注t时刻的平方根CD4计数与同一时间点的死亡风险之间的关联强度(P<0.05为差异有统计学意义)。
(2)建模分析:
① 联合模型:纵向子模型,建立纵向子模型:

ddIi是ddI组的虚拟变量,固定效应部分包括时间(obstime)的主效应、治疗与时间的交互作用、性别、prevOI和AZT;随机效应部分包括截距和时间项。
生存子模型,建立生存子模型:

生存子模型中,包括治疗方式、性别、AZT、prevOI和从纵向子模型中估计的时间-依赖性平方根CD4的真实变化轨迹mi(t)。这里假设基线风险函数是分段常数。最后,使用JM包中的函数jointModel()拟合联合模型,由函数summary()可获得详细的结果输出。
② 扩展的Cox模型:传统的Cox模型不能处理时依性协变量,扩展的Cox模型(extended Cox model,也称Andersen-Gill模型)是标准Cox模型的扩展,可以处理外源性的时依性协变量,但不能考虑到纵向过程的测量误差[14-15]。我们将平方根CD4计数作为外源性协变量放入扩展的Cox模型。

(3)结果分析:联合模型和扩展Cox模型的结果都显示,平方根CD4计数与死亡风险之间存在明显关联(P<0.000 1)。联合模型中平方根CD4计数的单位减少对应于死亡风险增加exp(-α)=1.29倍(95%CI:1.20~1.39),扩展的Cox模型中平方根CD4计数的单位减少对应于死亡风险exp(-α′)=1.18倍增加(95%CI:1.12~1.24)。点估计值和区间估计值二者都存在不可忽略的差异,这是由于扩展的Cox分析中未考虑到测量误差,导致的回归系数衰减。
(4)联合模型的诊断:当在实践中使用联合模型时,必须验证模型的假设,评估这些假设的标准工具是残差图(residual plots)[11]。
① 纵向子模型的残差:在标准线性混合效应模型中,经常使用个体特异性(条件)残差[subject-specific(conditional)residuals]和边际(人群平均)残差[the marginal(population averaged)residuals]。两种都可用于检查联合模型纵向子模型的假设。
图 1中的残差围绕在零附近波动,表明纵向模型的拟合比较好;图 2中当拟合值较小时,残差多为正值,可能是固定效应的设计矩阵的形式不合适所导致。
![]() |
图 1 个体特异性残差与拟合值的图像 |
![]() |
图 2 边际残差与拟合值的图像 |
② 生存子模型的残差:联合模型的相对风险子模型的标准残差类型是鞅残差(martingale residuals)。鞅残差可以被视为第i个受试者在时间t观察到的事件(这里指死亡)数量与基于拟合模型的同一事件的预期事件数之间的差异。
从图 3中可以观察到,对于较小的拟合值,鞅残差多表现为负值,表明需要进一步研究目前的纵向方程是否很好的拟合纵向变化过程,可以在调整其他基线协变量的同时,检查鞅残差的整体变化趋势。
![]() |
图 3 鞅残差与拟合值的图像 |
5.讨论:联合模型同时建模纵向数据和生存数据,相对于单独分析的主要优点是能够正确处理噪声和未能完全观察到的时依性协变量信息,这使得能够无偏估计纵向过程和生存过程之间的关系。而且,通过使用随机效应可以很容易地将个体差异纳入其中,这些特点使其在医学动态预测中也具有独特优势。随着个性化医疗的发展,在临床视角对患者进行个体动态预测的兴趣日益增加,基于联合模型的动态预测是目前比较流行的方法。例如Garre等[7]利用肾移植患者血清肌酐浓度倒数的变化轨迹来预测移殖失败的时间;Proust-Lima和Taylor[16]使用前列腺癌患者放治后的前列腺特异性抗原变化轨迹动态预测患者前列腺癌的复发;Andrinopoulou等[17]使用严重主动脉瓣狭窄患者纵向脑尿钠肽数据实现生存的动态预测等。针对多变量纵向数据和生存数据的联合模型、复发事件和纵向数据的联合模型、多变量生存数据和纵向数据的联合模型也被提出并广泛应用到医学研究中。同时,应用随机效应联合模型有2个注意事项:首先,随机效应联合模型依赖关于随机效应的假设。然而,最近的研究表明,关于随机效应的假设的稳健性随着每位受试者的观察数量而增加[18];其次,模型计算的复杂性是一个挑战。需要进一步研究开发参数估计方法,以便通过联合建模方法来分析更大的数据集[19]。
随着研究的逐步深入,联合模型形式出现多种扩展,其中包括处理多种失效类型,考虑分类纵向结果,通过潜在类别而不是随机效应将2种结果联系起来等。此外,有学者对联合模型中的2个子模型进行替代研究,用非线性混合效应模型[20]、广义线性混合效应模型[21]替代标准联合模型中的线性混合效应模型拟合纵向数据;用加速失效时间模型、转变模型,治愈模型替代Cox模型等[22-25]。R中的软件包JM、joineR和JMbayes;SAS中的PROC NLMIXED;Stata中的stjm命令,都可用于拟合基本的联合模型。国内目前关于联合模型的介绍还比较少,联合模型的应用范围也不够广泛,应将联合模型与我国的医疗卫生数据结合起来,在医学研究中发挥其作用。
利益冲突 所有作者均声明不存在利益冲突
[1] |
Abrams DI, Goldman AI, Launer C, et al. A comparative trial of didanosine or zalcitabine after treatment with zidovudine in patients with human immunodeficiency virus infection[J]. N Engl J Med, 1994, 330(10): 657-662. DOI:10.1056/nejm199403103301001 |
[2] |
Wulfsohn MS, Tsiatis AA. A joint model for survival and longitudinal data measured with error[J]. Biometrics, 1997, 53(1): 330-339. DOI:10.2307/2533118 |
[3] |
Tsiatis AA, Degruttola V, Wulfsohn MS. Modeling the relationship of survival to longitudinal data measured with error. Applications to survival and CD4 counts in patients with AIDS[J]. J Am Stat Assoc, 1995, 90(429): 27-37. DOI:10.2307/2291126 |
[4] |
Król A, Mauguen A, Mazroui Y, et al. Tutorial in joint modeling and prediction:a statistical software for correlated longitudinal outcomes, recurrent events and a terminal event[J]. J Stat Softw, 2017, 81(3). DOI:10.18637/jss.v081.i03 |
[5] |
Ibrahim JG, Chu H, Chen LM. Basic concepts and methods for joint models of longitudinal and survival data[J]. J Clin Oncol, 2010, 28(16): 2796-2801. DOI:10.1200/jco.2009.25.0654 |
[6] |
Andrinopoulou ER, Rizopoulos D, Jin RY, et al. An introduction to mixed models and joint modeling:analysis of valve function over time[J]. Ann Thorac Surg, 2012, 93(6): 1765-1772. DOI:10.1016/j.athoracsur.2012.02.049 |
[7] |
Garre FG, Zwinderman AH, Geskus RB, et al. A joint latent class changepoint model to improve the prediction of time to graft failure[J]. J Roy Stat Soc A, 2008, 171(1): 299-308. DOI:10.1111/j.1467-985x.2007.00514.x |
[8] |
Tsiatis AA, Davidian M. Joint modeling of longitudinal and time-to-event data:an overview[J]. Stat Sin, 2004, 14(3): 809-834. DOI:10.2307/24307417 |
[9] |
Robins JM, Hernán MA, Brumback B. Marginal structural models and causal inference in epidemiology[J]. Epidemiology, 2000, 11(5): 550-560. DOI:10.1097/00001648-200009000-00011 |
[10] |
Klein JP, van Houwelingen HC, Ibrahim JG, et al. Handbook of survival analysis[M]. Boca Raton: CRC Press, 2016.
|
[11] |
Rizopoulos D. Joint models for longitudinal and time-to-event data:with applications in R[M]. Boca Raton: CRC Press, 2012.
|
[12] |
Henderson R, Diggle P, Dobson A. Joint modelling of longitudinal measurements and event time data[J]. Biostatistics, 2000, 1(4): 465-480. DOI:10.1093/biostatistics/1.4.465 |
[13] |
Song X, Davidian M, Tsiatis AA. A semiparametric likelihood approach to joint modeling of longitudinal and time-to-event data[J]. Biometrics, 2002, 58(4): 742-753. DOI:10.1111/j.0006-341x.2002.00742.x |
[14] |
Andersen PK, Gill RD. Cox's regression model for counting processes:a large sample study[J]. Ann Stat, 1982, 10(4): 1100-1120. DOI:10.1214/aos/1176345976 |
[15] |
Andersen PK, Borgan O, Gill RD, et al. Statistical models based on counting processes[M]. New York: Springer Science & Business Media, 2012.
|
[16] |
Proust-Lima C, Taylor JMG. Development and validation of a dynamic prognostic tool for prostate cancer recurrence using repeated measures of posttreatment PSA:a joint modeling approach[J]. Biostatistics, 2009, 10(3): 535-549. DOI:10.1093/biostatistics/kxp009 |
[17] |
Andrinopoulou ER, Rizopoulos D, Geleijnse ML, et al. Dynamic prediction of outcome for patients with severe aortic stenosis:application of joint models for longitudinal and time-to-event data[J]. BMC Cardiovasc Disor, 2015, 15: 28. DOI:10.1186/s12872-015-0035-z |
[18] |
Gould AL, Boye ME, Crowther MJ, et al. Joint modeling of survival and longitudinal non-survival data:current methods and issues. Report of the DIA Bayesian joint modeling working group[J]. Stat Med, 2015, 34(14): 2181-2195. DOI:10.1002/sim.6141 |
[19] |
McCrink LM, Marshall AH, Cairns KJ. Advances in joint modelling:a review of recent developments with application to the survival of end stage renal disease patients[J]. Int Stat Rev, 2013, 81(2): 249-269. DOI:10.1111/insr.12018 |
[20] |
Wu L. A joint model for nonlinear mixed-effects models with censoring and covariates measured with error, with application to AIDS studies[J]. J Am Stat Assoc, 2002, 97(460): 955-964. DOI:10.1198/016214502388618744 |
[21] |
Pulkstenis EP, Ten Have TR, Landis JR. Model for the analysis of binary longitudinal pain data subject to informative dropout through remedication[J]. J Am Stat Assoc, 1998, 93(442): 438-450. DOI:10.1080/01621459.1998.10473693 |
[22] |
Tseng YK, Hsieh F, Wang JL. Joint modelling of accelerated failure time and longitudinal data[J]. Biometrika, 2005, 92(3): 587-603. DOI:10.1093/biomet/92.3.587 |
[23] |
Zeng D, Lin DY. Maximum likelihood estimation in semiparametric regression models with censored data[J]. J Roy Stat Soc B, 2007, 69(4): 507-564. DOI:10.1111/j.1369-7412.2007.00606.x |
[24] |
Yu MG, Taylor JMG, Sandler HM. Individual prediction in prostate cancer studies using a joint longitudinal survival-cure model[J]. J Am Stat Assoc, 2008, 103(481): 178-187. DOI:10.1198/016214507000000400 |
[25] |
Law NJ, Taylor JMG, Sandler H. The joint modeling of a longitudinal disease progression marker and the failure time process in the presence of cure[J]. Biostatistics, 2002, 3(4): 547-563. DOI:10.1093/biostatistics/3.4.547 |