Skip to content

生存分析与纵向数据

许多生物学研究关心从一个明确起点到事件发生经历了多长时间。种子萌发、个体死亡、疾病复发、移植物排斥和动物首次繁殖都属于时间到事件资料。研究结束时,一部分观察对象尚未发生事件;另一些对象在研究中途失访,或在进入研究之前已经度过了一段风险期。它们提供了带有观测边界和时间次序的信息,生存分析据此描述事件随时间积累的过程。

纵向研究则在同一个体、样地或实验单位上反复测量响应。它能够把“不同个体在某一时点的差异”与“同一个体随时间的变化”区分开来,同时需要处理个体内相关、不规则访视、时间变化协变量和失访。时间到事件与纵向资料常在同一研究中相遇:例如反复测量的生物标志物既预测死亡,死亡又终止后续测量。分析应先确定时间轴、观测机制和希望估计的生物学量,再选择相应模型。

时间轴、事件与观测机制

时间原点与时间尺度

令随机变量 \(T\geq 0\) 表示从时间原点到事件的时间。时间原点必须对所有对象具有可比的生物学意义,例如随机分组、确诊、孵化或首次暴露;“进入数据库的日期”只有在它与风险过程一致时才是合适的原点。年龄、日历时间和入组后的随访时间可能同时影响风险。Cox 模型等方法只能沿一条主要时间尺度建立风险集,其余时间尺度可作为协变量、分层因素或更细致模型的一部分。

事件定义也须在看数据之前固定。全因死亡只有一个终点,疾病特异性死亡则依赖死因判定;“首次复发”“任一复发”和“复发次数”对应不同的观察过程。若把相互排斥的不同结局合成复合终点,应说明各组成事件及其临床或生物学重要性,并分别呈现频繁事件和重点终点的贡献。

删失、截断与风险集

删失(censoring)表示对象已经进入观察,但事件时间只能定位在一个范围内。最常见的右删失发生在研究结束或失访时,此时只知道 \(T>C\);左删失只知道 \(T\leq C\),例如首次检测时已感染而感染时刻未知;区间删失则只知道 \(L<T\leq R\),常见于按固定间隔检查萌发或发病的研究。把区间终点任意记在首次阳性或区间中点,会人为改变时间分布;有足够区间删失时应使用相应的区间删失方法。1

截断(truncation)改变的是对象能否进入样本。左截断又称延迟入组:只有在入组时仍未发生事件的对象才会被观察,入组前已经发生事件者不在样本中。对象应从自己的入组时刻开始进入风险集。左删失记录已观察对象的早期事件时刻边界,左截断则使早期失败对象未能进入样本,二者对应不同的观测机制。

时刻 \(t\) 的风险集由刚到 \(t\) 之前已经进入观察、尚未发生目标事件且仍在随访的对象组成。删失对象在删失时刻之前仍贡献信息,之后才退出风险集。标准生存方法通常要求在已纳入的协变量和既往信息条件下,删失与未来事件过程独立。行政性研究截止常较接近这一假设;因病情恶化而退出则可能是信息性删失,需要补充协变量、加权、敏感性分析或联合建模来表达其观测机制。

生存函数、风险函数与累积风险

生存函数给出超过时刻 \(t\) 仍未发生事件的概率:

\[ S(t)=P(T>t)=1-F(t). \]

连续时间下的风险函数(hazard)是条件于刚到 \(t\) 仍处于风险中时,单位时间内的瞬时事件率:

\[ h(t)=\lim_{\Delta t\to 0} \frac{P(t\leq T<t+\Delta t\mid T\geq t)}{\Delta t}. \]

风险率可以大于 1,因为它是率而非某一时间区间内的概率。累积风险为

\[ H(t)=\int_0^t h(u)\,\mathrm du, \]

因此生存函数也可由累积风险写成

\[ S(t)=\exp\{-H(t)\}. \]

这三个函数观察同一事件过程的不同侧面。\(S(t)\) 适合表达“存活到某时点的比例”,\(h(t)\) 适合描述仍处于风险中的对象此刻经历事件的强度,\(H(t)\) 则把沿途的风险率累积起来。风险比(hazard ratio,HR)比较两个条件风险率,与某时点风险比和中位生存时间之比分属不同效应尺度。风险随时间的形状和研究问题共同决定应报告哪一种尺度。

Kaplan–Meier 估计与组间比较

乘积极限估计

设不同事件时刻为 \(t_1<t_2<\cdots\),在 \(t_j\) 刚发生前有 \(n_j\) 个对象处于风险集,其中 \(d_j\) 个发生事件。Kaplan–Meier 乘积极限估计为

\[ \widehat S(t)= \prod_{t_j\leq t}\left(1-\frac{d_j}{n_j}\right). \]

每个事件时刻使曲线下降;删失使该对象在此后退出风险集,而曲线高度保持到下一个事件时刻。阶梯曲线由一连串条件存活比例相乘得到,这条经典推导利用了删失前的全部风险集信息。Kaplan 与 Meier 的原始论文建立了这一非参数估计及其大样本性质。2

Greenwood 公式给出常用方差近似:

\[ \widehat{\operatorname{Var}}\{\widehat S(t)\} =\widehat S(t)^2 \sum_{t_j\leq t}\frac{d_j}{n_j(n_j-d_j)}. \]

实际置信区间常在 log–log 等变换尺度上构造,使区间保持在 \([0,1]\)。中位生存时间是 \(\widehat S(t)\) 首次降至 0.5 或以下的时刻;随访期间曲线未越过 0.5 时应报告“尚未达到”,因为观测尾部之外缺少估计支持。图中应同时给出删失标记、若干关键时点的风险集人数和置信带。曲线尾端若只剩很少对象,阶梯仍可计算,但区间通常很宽。

Nelson–Aalen 累积风险估计为

\[ \widehat H(t)=\sum_{t_j\leq t}\frac{d_j}{n_j}. \]

\(\exp\{-\widehat H(t)\}\) 与 Kaplan–Meier 曲线通常接近,但有限样本中并不完全相同。前者在累积风险、复发事件和多状态模型中尤其自然。

对数秩检验与受限平均生存时间

对数秩(log-rank)检验在每个事件时刻按当时风险集计算各组期望事件数,再把“观察数减期望数”沿时间合并。它检验整段随访中各组生存曲线是否相同,效应大小需由生存概率、RMST 或回归模型另行给出;其权重使它对近似比例风险的组间差异较敏感。检验依赖组内非信息性删失,各组可以具有不同的删失时间分布。3

曲线明显交叉、处理效应延迟出现或只在早期存在时,单一 log-rank \(P\) 值和恒定风险比都可能掩盖时间结构。若替代权重、里程碑时点或分析窗口由观察曲线后才决定,会增加选择性解释。受限平均生存时间(restricted mean survival time,RMST)把预先规定的 \(\tau\) 之前的生存曲线面积作为效应量:

\[ \operatorname{RMST}(\tau)=\int_0^\tau S(t)\,\mathrm dt. \]

两组 RMST 之差可解释为在 \(0\)\(\tau\) 内平均无事件时间的差异。\(\tau\) 必须处于各组都有充分支持的随访范围,并在分析前确定;RMST 概括的是这一特定时间窗口。4

Cox 比例风险模型

部分似然与条件风险比

Cox 模型把协变量作用写成风险率上的乘法结构:

\[ h(t\mid\mathbf x)=h_0(t)\exp(\mathbf x^{\mathsf T}\beta), \]

其中 \(h_0(t)\) 是未指定具体形状的基准风险。若连续协变量 \(x_j\) 增加一个单位,其他协变量相同条件下的风险率乘以 \(e^{\beta_j}\)。模型在每个事件时刻比较事件对象与当时风险集中的对象,借助部分似然估计参数向量 \(\beta\),无需先指定 \(h_0(t)\) 的参数分布。Cox 的原始论文把这一条件似然思想用于带协变量的寿命资料。5

Cox 模型的半参数性指基准风险形状未作参数化,同时仍要求风险比在时间上保持恒定、协变量函数形式正确、删失在条件上非信息性,并恰当处理同一对象、家系或中心造成的相关性。事件时刻有大量并列值时,还应说明使用 Breslow、Efron 或精确方法。连续变量通常先用样条检查非线性,分组阈值则需要预设依据;交互项与分层因素也应由研究设计和生物机制决定。

比例风险诊断与时间变化效应

比例风险可通过分组的 log-minus-log 曲线、Schoenfeld 残差与时间的关系,以及预先设定的协变量—时间交互来检查。Grambsch 与 Therneau 将加权 Schoenfeld 残差与时间变化系数联系起来,使比例风险假设可以按变量诊断。6 一个显著的全局检验提示模型结构需要重新考察,但图形中的偏离形状、效应大小和时间范围同样重要;样本很大时微小偏离也可显著,事件很少时检验又可能缺乏能力。

比例风险出现偏离时,可按问题选择不同处理:若某因素只影响基准风险而其系数不是研究目标,可按它分层;若效应确实随时间改变,可使用分段效应或 \(\beta(t)\);也可改用灵活参数模型、里程碑效应或 RMST。报告时应直接给出不同时间段的效应或预测曲线,说明平均风险比所概括的时间范围。

时间依赖协变量与时间变化系数是两个不同概念。前者指 \(x_i(t)\) 的取值随时间更新,后者指同一协变量的回归作用 \(\beta(t)\) 随时间改变。时间依赖协变量通常把每个对象的数据拆成 \((\text{start},\text{stop}]\) 区间;在某个事件时刻进入风险集的只能是该时刻之前已经知道的协变量值。用事件之后才测得的值解释事件之前的风险,会产生前视偏倚;把“必须存活到接受处理”期间错误记入处理组,则形成不死时间偏倚。官方 survival 说明以风险集和计数过程形式详细展示了这一时间顺序原则。7

参数模型与灵活参数模型

参数生存模型为时间分布或风险函数指定形状。指数分布假定风险恒定;Weibull 分布允许风险单调增加或下降;Gompertz 模型常用于随年龄近似指数变化的死亡率;log-normal 与 log-logistic 分布可以产生先升后降的风险。指定分布后可以直接估计绝对生存概率、分位数和更远时间的外推,但外推结果强烈依赖尾部分布假设。

加速失效时间模型(accelerated failure time,AFT)常写为

\[ \log T=\mathbf x^{\mathsf T}\gamma+\sigma\varepsilon. \]

在相应参数化下,\(e^{\gamma_j}\) 是时间比:大于 1 表示事件时间被拉长。它与 Cox 模型的风险比回答不同尺度的问题,某些分布如 Weibull 可同时具有 PH 与 AFT 表示,另一些则不能。

灵活参数生存模型使用限制性立方样条描述 log 累积风险、log 累积优势或其他尺度,并可加入协变量与时间的交互,兼顾平滑的基准曲线和非比例风险。Royston–Parmar 模型是这一思路的重要实现。8 节点数量、边界行为和时间尺度需要检查。若模型用于超出数据范围的长期生存或成本效果推断,还应结合外部自然史资料并比较多个合理尾部假设。

竞争风险与多状态过程

当一种事件发生后会阻止另一种事件再发生时,存在竞争风险。例如个体死于其他原因后便不可能再观察目标疾病死亡。第 \(k\) 类原因特异风险描述尚未发生任何事件者立即经历该事件的速率;相应累积发生函数(cumulative incidence function,CIF)为

\[ F_k(t)=P(T\leq t,J=k) \]

它也可由仍无任何事件的生存概率与原因特异累积风险表示为

\[ F_k(t)=\int_0^t S(u-)\,\mathrm dH_k(u). \]

CIF 同时受第 \(k\) 类风险和其他竞争事件风险影响。把竞争事件当成普通右删失,再用 \(1-\widehat S_{KM}(t)\) 估计目标事件的现实概率,相当于假想竞争事件发生后仍可继续经历目标事件,通常会高估累积发生概率;各原因这样算出的曲线甚至可能相加超过 1。

原因特异 Cox 模型在竞争事件发生时将对象移出风险集,回答当前无事件者的瞬时原因特异风险怎样随协变量改变。Fine–Gray 模型直接回归亚分布风险,以保留已发生竞争事件的对象于修改后的风险集中来连接 CIF。9 两类模型的系数属于不同风险集和不同估计目标,应按病因机制或预后决策问题选择;P 值则用于量化给定模型下的证据。稳妥报告通常给出各结局的绝对累积发生概率,并明确所用模型。

疾病过程还可扩展为无病、患病、复发和死亡等多个状态。Aalen–Johansen 估计把各转移的风险增量组合成状态占有概率;只有两个状态时它退化为 Kaplan–Meier,在单次竞争风险结构中则给出各原因 CIF。10 多状态模型要求明确允许的转移、每条转移的时间尺度,以及 Markov 或半 Markov 假设。把多个状态压成单一“失败”虽然简单,却可能丢失疾病进程本身。

复发事件与终末事件

同一个体可多次感染、发作、产卵或住院。只分析首次事件具有清楚的随机化比较含义,也避开了后续事件受先前事件影响的复杂性,但会丢弃疾病负担和事件间隔的信息。描述阶段可画平均累积函数,表示截至时间 \(t\) 每个对象平均累积经历的事件数,并同时报告仍处于观察中的人数。

Andersen–Gill 模型把 Cox 计数过程扩展到重复事件,常在总时间尺度上让对象每次事件后重新进入风险集,并用个体聚类稳健方差处理同一对象内相关;它估计的是事件强度关系。11 Prentice–Williams–Peterson 模型按事件次序分层,只有已经经历第 \(k-1\) 次事件者才进入第 \(k\) 次事件的风险集,可采用从研究起点计时的总时间或从上一次事件起计时的间隔时间。12 Frailty 模型则用共享随机效应表示未观测的个体易感性,所得系数是条件于该易感性的效应。

死亡等终末事件会中止复发观察,需要在模型中另行表达。高复发风险者若更早死亡,观察到的复发次数反而可能较少。分析必须说明估计目标是首次事件、总事件负担、事件强度、某次序事件还是存活期间的复发过程,并考虑联合 frailty、复发—终末事件模型或多状态模型。事件定义、事件后是否继续处于风险、风险间隔、聚类单位和终末事件处理都属于结果的一部分。

纵向资料与变化轨迹

重复横断面与个体随访

重复横断面调查在多个时点抽取不同个体,能够描述总体均值怎样变化;纵向设计在同一单位上反复测量,可以进一步估计个体内变化、个体间轨迹差异以及早期状态与后续结局的关联。两者即使每个时点样本量相同,抽样单位和可回答的问题也不同。

纵向数据通常整理成长格式:每行对应一个“单位—时点”观察,至少保留单位标识、实际测量时间、响应、随时间变化的协变量和基线协变量。访视编号适合表示程序次序,实际年龄或距随机化时间则保留真实间隔;计划在第 3 月测量的记录可能因访视窗口而具有不同实际间隔。分析前应检查每个对象的测量次数、时间分布和失访位置,同时绘制个体轨迹图、分组均值与区间和各时点分布,以呈现个体变化方向。

时间函数、基线与变化量

时间可以作为分类因素,让每个访视有独立均值;也可用线性项、分段线性项或样条描述平滑轨迹。把时间中心化在基线、处理开始或其他有意义时刻,截距就对应那个时刻的平均响应。处理与时间的交互表示组间轨迹差异;若时间为分类因素,交互允许各访视差异独立变化,若时间为线性项,它只比较平均斜率。

同一个随时间变化的协变量同时包含两种信息。令个体 \(i\) 在时点 \(t\) 的协变量为 \(x_{it}\),可分解为

\[ x_{it}=\bar x_i+(x_{it}-\bar x_i). \]

\(\bar x_i\) 描述长期平均水平不同的个体之间怎样不同,\(x_{it}-\bar x_i\) 描述某个体偏离自身平均水平时响应怎样改变。只放入原始 \(x_{it}\) 往往把个体间与个体内关联混在一个系数里。时间变化暴露还可能受先前响应影响;其因果效应解释需要处理时间变化混杂和反馈。

基线值既可作为后续响应的协变量,也可作为纵向响应的一部分。随机试验中,比较随访值并调整基线通常比只比较“随访减基线”的变化量更精确;变化量模型回答的是变化差异,并更直接继承两次测量的误差。约束纵向模型可把基线也纳入响应,同时利用随机化所蕴含的组间共同基线均值。观察性研究的基线调整则需要结合时间顺序和混杂结构解释。无论采用哪种方式,都应预先说明估计的是某时点均值差、平均斜率差、曲线面积还是变化量差。

个体内相关与模型选择

同一个体的相邻测量通常比相隔很久的测量更相似。复合对称结构假定任意两个时点相关相同;AR(1) 假定等间隔访视的相关随间隔阶数指数衰减;连续时间指数相关适合实际间隔不规则的情形;无结构协方差为每对时点分别估计协方差,灵活但参数数目随访视数迅速增长。随机截距和随机斜率则用个体轨迹的层级差异生成相关性。随机效应协方差与给定随机效应后的残差协方差是两个层次,必要时可同时存在。

经典重复测量方差分析适合设计平衡、时点固定且缺失很少的连续响应;其球形性假设要求任意两时点差值具有相同方差。违反球形性时可修正自由度,多变量方差分析以更完整的数据和足够样本换取对球形性的放宽。它们提供传统“组别、时间、组别×时间”分解;不规则时间、不同测量次数和层级嵌套出现后,混合模型通常更自然。

线性或广义线性混合模型通过随机效应描述个体特异轨迹,估计条件于随机效应的关系;GEE 直接估计总体平均关系,并用工作相关结构提高效率和构造稳健方差。两者对应不同的目标和系数解释,在非线性链接下差异尤其明显。本章前一页已系统说明随机截距、随机斜率与个体条件效应,以及 GEE 的工作相关与总体平均解释。Laird 与 Ware 的经典随机效应框架说明了不平衡、非等间隔纵向资料怎样分离个体内与个体间变异;Liang 与 Zeger 则建立了相关响应的 GEE 路径。1314

模型选择应回到估计目标。研究个体生长轨迹、异质性和个体预测时,混合模型通常更合适;研究干预在总体平均上的影响时,GEE 可能更直接。聚类数很少时,GEE 的三明治(sandwich)标准误会有明显小样本偏差;随机效应分布、协方差结构和残差方差若严重错设,混合模型的推断也会受影响。比较模型时应综合拟合轨迹、残差、自相关、个体影响和 AIC,并让相关结构的复杂度与数据量相配。

缺失、脱落与访视过程

纵向缺失需要同时记录原因和时间模式。完全随机缺失(missing completely at random,MCAR)表示在所考虑的资料中,观察概率不依赖已观察或未观察的值;随机缺失(missing at random,MAR)表示给定已经观察到的历史和协变量后,观察概率不再依赖当前未观察值;非随机缺失(missing not at random,MNAR)表示即使给定这些信息,缺失仍依赖未观察值。MAR 是关于给定已观测信息后缺失概率的条件性假设。在似然或贝叶斯推断中,“可忽略”通常还要求缺失机制参数与结局模型参数可分离。Rubin 的框架奠定了这些定义及其可忽略条件。15

正确指定的似然型混合模型可在给定已观察资料的 MAR 下使用所有已观测结局;普通未加权 GEE 的一致性通常需要更强的 MCAR 条件,MAR 情形可考虑逆概率加权 GEE、多重插补或适当的联合模型。完整病例分析会丢弃部分轨迹,其偏倚取决于缺失机制;末次观测结转假定退出后的状态保持不变,通常低估不确定性。国家研究委员会关于临床试验缺失资料的报告强调,应在设计阶段减少缺失、保存所有随机化对象的信息、明确主要假设并进行敏感性分析。16

MAR 与 MNAR 对未观测值作出不同假设,通常无法仅凭已观察数据区分。对信息性失访可采用选择模型、模式混合模型、共享参数模型或 delta 调整,考察在一组科学上合理的不可识别参数下结论是否改变。死亡后的生活质量因死亡而失去定义,需要复合结局、存活者估计目标、while-alive 过程或其他明确策略。

访视过程本身也可能携带信息。病情恶化者更频繁就诊会使测量时点与潜在结局相关;只在已到访记录上拟合轨迹,可能把“更常被观察”误作“总体变化”。研究应保存计划访视和实际访视、未访原因、治疗停止、事件与死亡,并在必要时联合建模访视强度或进行加权分析。

联合纵向—事件模型

当一个带测量误差的纵向标志物既与事件风险相关,又因事件或健康恶化而信息性终止时,分成两个独立模型常不充分。联合模型用纵向子模型估计潜在真实轨迹 \(m_i(t)\),再把它与事件子模型连接。例如共享参数模型可写为

\[ h_i(t)=h_0(t)\exp\left\{ \gamma^{\mathsf T}\mathbf w_i+ \alpha m_i(t) \right\}. \]

纵向子模型和事件子模型通过共享随机效应或潜在轨迹相关联。\(\alpha\) 表示给定模型后当前潜在标志物水平与事件风险的关联;也可连接轨迹斜率、累积暴露或多个标志物。联合模型可用于修正内源性时间依赖协变量的测量误差、处理与结局相关的失访,以及在新测量到来后更新个体动态风险预测。Tsiatis 与 Davidian 系统整理了这类共享随机效应模型及其估计目标。17

联合模型中的 \(\alpha\) 描述指定模型下的关联,其因果含义还需要相应设计和识别假设。纵向轨迹函数、随机效应分布、删失、基准风险和连接形式都可能影响结果;死亡前测量稀少时,斜率关联尤其难以识别。应把联合模型与较简单分析对照,展示个体拟合轨迹与事件校准,并对连接结构和缺失假设做敏感性分析。

分析流程与完整报告

一项时间到事件分析应从研究设计重新构造时间轴:核对时间原点、事件日期、入组日期、最后确认无事件日期以及竞争事件,识别延迟入组和同日并列事件;随后用事件数、删失数、随访分布、风险表和非参数曲线了解资料支持的时间范围。模型结果至少给出效应估计、置信区间、绝对风险和 P 值。

Cox 或参数模型还应报告时间尺度、协变量编码与非线性、交互或分层、并列事件处理、聚类处理、比例风险检查及其后续策略。竞争风险分析要说明估计的是原因特异风险、亚分布风险还是 CIF;复发事件分析要写明事件次序、总时间或间隔时间、事件后风险状态和终末事件。用于预测时,训练与验证必须按对象而非记录拆分,并在预定时间点检查校准、时间依赖判别指标和 Brier 分数。

纵向分析应报告计划与实际测量时间、每个对象的测量次数、时间函数、基线处理、固定效应与随机效应、残差或工作相关结构,以及目标是个体条件效应还是总体平均效应。缺失比例应按访视和组别呈现,说明退出原因、主要缺失假设、插补或加权模型及 MNAR 敏感性分析。模型图应展示估计轨迹及区间,同时保留原始数据的时间分布;一条平滑均值曲线无法替代对个体差异、尾部稀疏和失访过程的说明。

时间分析的共同原则是让风险集和信息时间保持一致。谁在某时刻真正处于风险中、哪些协变量在那一刻已经可知、谁仍可能被测量,决定了估计量所比较的人群。明确这些条件后,经典生存曲线、回归模型、纵向轨迹和现代联合模型才会成为同一观察过程的连续描述。

参考资料与延伸阅读


  1. Penn State Eberly College of Science. 5.2 - Special Types of Event Times

  2. Kaplan, E. L. & Meier, P. (1958). Nonparametric Estimation from Incomplete Observations. Journal of the American Statistical Association, 53, 457–481. 

  3. Penn State Eberly College of Science. Lesson 5: Sample Size Calculation and Power(其中生存函数估计与 log-rank 检验部分). 

  4. Royston, P. & Parmar, M. K. B. (2013). Restricted mean survival time: an alternative to the hazard ratio for the design and analysis of randomized trials with a time-to-event outcome. BMC Medical Research Methodology, 13, 152. 

  5. Cox, D. R. (1972). Regression Models and Life-Tables. Journal of the Royal Statistical Society: Series B, 34, 187–220. 

  6. Grambsch, P. M. & Therneau, T. M. (1994). Proportional hazards tests and diagnostics based on weighted residuals. Biometrika, 81, 515–526. 

  7. Therneau, T. M., Crowson, C. S. & Atkinson, E. J. Using Time Dependent Covariates and Time Dependent Coefficients in the Cox Model. survival package vignette. 

  8. Royston, P. & Parmar, M. K. B. (2002). Flexible parametric proportional-hazards and proportional-odds models for censored survival data. Statistics in Medicine, 21, 2175–2197. 

  9. Fine, J. P. & Gray, R. J. (1999). A Proportional Hazards Model for the Subdistribution of a Competing Risk. Journal of the American Statistical Association, 94, 496–509. 

  10. Therneau, T. M. et al. Multi-state models and competing risks. survivalVignettes tutorial. 

  11. Andersen, P. K. & Gill, R. D. (1982). Cox's Regression Model for Counting Processes: A Large Sample Study. The Annals of Statistics, 10, 1100–1120. 

  12. Prentice, R. L., Williams, B. J. & Peterson, A. V. (1981). On the regression analysis of multivariate failure time data. Biometrika, 68, 373–379. 

  13. Laird, N. M. & Ware, J. H. (1982). Random-Effects Models for Longitudinal Data. Biometrics, 38, 963–974. 

  14. Liang, K.-Y. & Zeger, S. L. (1986). Longitudinal data analysis using generalized linear models. Biometrika, 73, 13–22. 

  15. Rubin, D. B. (1976). Inference and missing data. Biometrika, 63, 581–592. 

  16. National Research Council. (2010). The Prevention and Treatment of Missing Data in Clinical Trials. National Academies Press. 

  17. Tsiatis, A. A. & Davidian, M. (2004). Joint Modeling of Longitudinal and Time-to-Event Data: An Overview. Statistica Sinica, 14, 809–834. 

页面讨论

使用 GitHub 登录后可参与整页讨论;评论独立保存在 GitHub Discussions 中。