生存分析
2018-09-03
2018-09-03
1 基本概念
生存分析可以用来估计总体的生存率及其有关指标,如癌症患者治疗后的生存时间资料;还可以对不同处理组生存率进行比较,如比较哪种癌症治疗方法更有效果;同时还可探索生存时间资料的影响因素;预测不同的影响因素水平的个体生存时间和一定时间的生存率。近年来也越来越多人用在互联网数据挖掘中,例如去预测信息在社交网络的传播程度,或者去预测用户流失的概率。
1. 生存分析
是研究生存现象和响应时间数据及其统计规律的一门学科。将事件的结果(终点事件)和出现这一结果所经历的时间结合起来分析的一种统计分析方法。生存分析与其他多因素分析最大的区别,就是生存分析考虑了观测结果出现的时间长短。
2. 事件
指研究中规定的生存研究的终点,在研究开始之前就已经制定好。根据研究性质的不同,事件可以是患者的死亡、疾病的复发、仪器的故障,也可以是下岗工人的再就业等等。
3. 生存时间
指从某一起点到事件发生所经过的时间。生存时间不仅仅指医学中的存活,也可以是机器出故障前的正常运行时间等等。
终点事件与起始事件之间的时间间隔。
终点事件:研究者所关系的特定结局。
起始事件:反应研究对象生存过程的起始特征的事件。
终点事件与起始事件是相对而言的,都是由特定的研究目的所决定的,是整个研究过程的标尺,需要在设计时明确规定,并在研究期间严格遵守,不能随意改变。比如研究药品是否起作用和起作用的时间问题,可以将服药作为起始事件,终点事件设为痊愈。
4. 生存时间的类型
广义上指某个起点事件开始到某个终点事件发生所经历的时间,度量单位可以是年、月、日、小时等,常用符号t所示。这个时间也未必是通常意义上的时间,也可以是和时间相关的变量。比如距离等,具体要根据研究目的而定义。根据研究对象的结局,生存时间数据可分为两种类型:
完全数据:从观察起点到发生死亡事件所经历的时间。
不完全数据:生存时间观察过程的截止不是由于死亡事件,而是由其他原因引起的。不完全数据分为:删失数据,截断数据。
不完全主要原因:
- 失访:指失去联系;
- 退出:死于非研究因素或非处理因素而退出研究;
- 终止:设计时规定的时间已到而终止观察,但研究对象仍然存活。
删失的表现形式:
- 右删失:只知道实际寿命大于某数;
- 左删失:只知道实际寿命小于某数;
- 区间删失:只知道实际寿命在一个时间区间内。
截尾数据的产生:
往往是因为实验设计的要求使得数据天然具有上界或者下界。
- 左截尾:数据都大于某个值。
- 右截尾:数据都小于某个值。
如一个实验研究退休职工的生存情况,那么显然这些数据都是左截尾的,因为所有个体的年龄都大于退休年龄(如t≥60)
5. 条件生存概率
生存时间的分布一般不呈现正态分布。条件生存概率表示某时段开始存活的个体,到该时段结束时仍存活的可能性。
6. 生存率
生存概率又称为生存率或生存函数,表示观察对象经历t个单位时间段后仍存活的可能性。
生存率与条件生存概率不同。条件生存概率是单个时段的结果,而生存率是条件生存概率的累积,是多个时段的累积结果。例如,3年生存率是第1年存活,第2年也存活,第3年还存活的可能性。
7. 生存曲线
以观察时间为横轴,以生存率为纵轴,将各个时间点所对应的生存率连接在一起的曲线图。
生存曲线是一条下降的曲线,分析时应注意曲线的高度和斜率。平缓的生存曲线表示高生存率或较长生存期,陡峭的生存曲线表示低生存率或较短生存期。
8. 中位生存期
又称半数生存期,表示恰好有50%的个体尚存活的时间。中位生存期越长,表示疾病的预后越好。
9.生存函数S(t)
观察对象的生存时间T大于某时刻t的概率。
\[S(t)=\frac{时刻t的存活例数}{期初观察例数}\]
10. 死亡函数F(t)
观察对象的生存时间T不大于某时刻t的概率。
\[F(t)=1-S(t)\]
11. 死亡密度函数f(t)
观察对象在某时刻t的瞬时死亡率。
\[f(t)=\frac{\nabla{F(t)}}{\nabla{t}}\]
12. 风险函数h(t)
生存到时刻t的观察对象在时刻t的瞬时死亡率。 \[h(t)=\frac{f(t)}{S(t)}\]
2 生存分析的方法
1.描述法
根据样本观测值直接计算每一个时间点或时间区间上的生存函数、死亡函数、风险函数等,并采用列表或绘图的形式显示生存时间的分布规律。方法对数据分布无要求。
2.非参数法
生存时间的分布未知,估计生存函数;或检验危险因素对生存时间的影响。常用乘积极限法和寿命表法,可以用来估计生存函数,比较两组或多组生存分布函数。分析危险因素对生存时间的影响。
3.参数法
根据样本观测值来估计假设分布模型中的参数,获得生存时间的概率分布模型。生存时间经常服从的分布有:指数分布、Weibull分布、对数正态分布、对数Logistic分布、Gamma分布。不仅能处理非参数法解决的问题,同时还可建立生存时间与危险因素之间的关系模型。但是我们必须提前知道生存时间的分布。
4.半参数法
不需要对生存时间的分布做出假定,但是却可以通过一个模型来分析生存时间的分布规律,以及危险因素对生存时间的影响,最著名的就是COX回归。不仅能处理参数法解决的问题,还不需要事先知道生存时间的分布。
3 R实现
R实现生存分析涉及survival包和OIsurv包,同时数据集来自于KMsurv包。
# 在第一次导入包之前记得安装相关的包
library(KMsurv)
library(survival)
library(OIsurv)
# 查看数据集
data(aids)
head(aids)
## infect induct adult
## 1 0.00 5.00 1
## 2 0.25 6.75 1
## 3 0.75 5.00 1
## 4 0.75 5.00 1
## 5 0.75 7.25 1
## 6 1.00 4.25 1
3.1 生存类型数据
survival包中的许多函数将方法应用到Surv对象,这些对象是使用Surv()函数创建的生存类型对象。在这里,我们讨论了右删失Surv对象的构造和左截尾右删失Surv对象的构造。
具体形式:
Surv(time, event), Surv(time, time2, event, type)
1. 右截尾数据
在Surv()函数中仅需要两个参数:一个时间向量和一个表示哪些时间被观察和删失的向量。
# 查看数据集
data(tongue)
head(tongue)
## type time delta
## 1 1 1 1
## 2 1 3 1
## 3 1 3 1
## 4 1 4 1
## 5 1 10 1
## 6 1 13 1
attach(tongue)
my.surv.object <- Surv(time[type==1], delta[type==1])
my.surv.object
## [1] 1 3 3 4 10 13 13 16 16 24 26 27 28 30
## [15] 30 32 41 51 65 67 70 72 73 77 91 93 96 100
## [29] 104 157 167 61+ 74+ 79+ 80+ 81+ 87+ 87+ 88+ 89+ 93+ 97+
## [43] 101+ 104+ 108+ 109+ 120+ 131+ 150+ 231+ 240+ 400+
detach(tongue)
标识加号是指被右截的观察结果。Surv()中的第一个参数应该输入观察时间和右截尾时间的向量。第二个参数使用指标向量来表示事件是否被观察。
2. 左截尾右删失数据
第一个参数输入左截时间,第二个参数输入生存事件向量和删失时间,第三个参数输入右删失指示向量。
# 查看数据集
data(psych)
head(psych)
## sex age time death
## 1 2 51 1 1
## 2 2 58 1 1
## 3 2 55 2 1
## 4 2 28 22 1
## 5 1 21 30 0
## 6 1 19 28 1
attach(psych)
my.surv.object <- Surv(age, age+time, death)
my.surv.object
## [1] (51,52] (58,59] (55,57] (28,50] (21,51+] (19,47] (25,57]
## [8] (48,59] (47,61] (25,61+] (31,62+] (24,57+] (25,58+] (30,67+]
## [15] (33,68+] (36,61] (30,61+] (41,63] (43,69] (45,69] (35,70+]
## [22] (29,63+] (35,65+] (32,67] (36,76] (32,71+]
detach(psych)
3.2 Kaplan-Meier估计
生存分析的主要目的是估计生存函数,常用的方法有Kaplan-Meier法。Kaplan-Meier估计是生存函数S(t)的非参数极大似然估计(MLE)。该估计是一个在观测事件时间\(t_i\)上跳跃的阶跃函数。
具体形式:
survfit(formula, conf.int = 0.95, conf.type = “log”)
# 查看数据集
data(tongue)
attach(tongue)
my.surv <- Surv(time[type==1], delta[type==1])
survfit(my.surv~1)
## Call: survfit(formula = my.surv ~ 1)
##
## n events median 0.95LCL 0.95UCL
## 52 31 93 67 NA
Survfit()也有一些可选参数。例如,置信水平使用第二个参数conf.int更改(例如,90%的置信范围,conf.int=0.90)。conf.type参数描述了置信区间的类型。默认的是“log”,它等于对转换函数g(t) =log(t)。“log-log”选项使用g(t)=log(−log(t))。一个线性置信区间是使用参数conf.type=“plain”创建的。与R中的许多函数一样,survfit()函数返回可以通过正确的命令访问隐藏信息。下面我们考虑这个隐藏信息的几个元素,它们存储在一个列表中。若需要查看整个信息列表的结构,可对my.fit和summary(my.fit)应用str函数。
# 查看数据集
my.fit <- survfit(my.surv~1)
summary(my.fit)$surv# 生存函数值
## [1] 0.9807692 0.9423077 0.9230769 0.9038462 0.8653846 0.8269231 0.8076923
## [8] 0.7884615 0.7692308 0.7500000 0.7115385 0.6923077 0.6730769 0.6538462
## [15] 0.6340326 0.6142191 0.5944056 0.5745921 0.5547786 0.5342312 0.5061138
## [22] 0.4779963 0.4481216 0.4161129 0.3814368 0.3051494 0.2288621
summary(my.fit)$time
## [1] 1 3 4 10 13 16 24 26 27 28 30 32 41 51 65 67 70
## [18] 72 73 77 91 93 96 100 104 157 167
summary(my.fit)$n.risk
## [1] 52 51 49 48 47 45 43 42 41 40 39 37 36 35 33 32 31 30 29 27 19 18 16
## [24] 14 12 5 4
summary(my.fit)$n.event
## [1] 1 2 1 1 2 2 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
summary(my.fit)$sstd.err
## NULL
summary(my.fit)$lower
## [1] 0.9441432 0.8810191 0.8534195 0.8271685 0.7774159 0.7302341 0.7073724
## [8] 0.6849189 0.6628317 0.6410776 0.5984672 0.5775712 0.5569274 0.5365233
## [15] 0.5156073 0.4949392 0.4745101 0.4543127 0.4343412 0.4137057 0.3837502
## [22] 0.3546085 0.3240114 0.2916674 0.2571797 0.1692469 0.1010963
str(my.fit)
## List of 13
## $ n : int 52
## $ time : num [1:45] 1 3 4 10 13 16 24 26 27 28 ...
## $ n.risk : num [1:45] 52 51 49 48 47 45 43 42 41 40 ...
## $ n.event : num [1:45] 1 2 1 1 2 2 1 1 1 1 ...
## $ n.censor : num [1:45] 0 0 0 0 0 0 0 0 0 0 ...
## $ surv : num [1:45] 0.981 0.942 0.923 0.904 0.865 ...
## $ type : chr "right"
## $ std.err : num [1:45] 0.0194 0.0343 0.04 0.0452 0.0547 ...
## $ lower : num [1:45] 0.944 0.881 0.853 0.827 0.777 ...
## $ upper : num [1:45] 1 1 0.998 0.988 0.963 ...
## $ conf.type: chr "log"
## $ conf.int : num 0.95
## $ call : language survfit(formula = my.surv ~ 1)
## - attr(*, "class")= chr "survfit"
str(summary(my.fit))
## List of 15
## $ n : int 52
## $ time : num [1:27] 1 3 4 10 13 16 24 26 27 28 ...
## $ n.risk : num [1:27] 52 51 49 48 47 45 43 42 41 40 ...
## $ n.event : num [1:27] 1 2 1 1 2 2 1 1 1 1 ...
## $ n.censor : num [1:27] 0 0 0 0 0 0 0 0 0 0 ...
## $ surv : num [1:27] 0.981 0.942 0.923 0.904 0.865 ...
## $ type : chr "right"
## $ std.err : num [1:27] 0.019 0.0323 0.037 0.0409 0.0473 ...
## $ lower : num [1:27] 0.944 0.881 0.853 0.827 0.777 ...
## $ upper : num [1:27] 1 1 0.998 0.988 0.963 ...
## $ conf.type : chr "log"
## $ conf.int : num 0.95
## $ call : language survfit(formula = my.surv ~ 1)
## $ table : Named num [1:9] 52 52 52 31 147 ...
## ..- attr(*, "names")= chr [1:9] "records" "n.max" "n.start" "events" ...
## $ rmean.endtime: num 400
## - attr(*, "class")= chr "summary.survfit"
Kaplan-Meier估计可以使用plot(my.fit)绘制。
plot(my.fit,main="Kaplan-Meier estimate with 95% confidence bounds", xlab = "time",ylab="survival function")
上面分析不同的组被包含在一个单一的Surv对象中。例如,tongue数据集中的type变量描述了患者的DNA谱。我们可以通过对type变量上的Surv对象进行回归得到这些组的Kaplan-Meier估计:
my.fit1 <- survfit( Surv(time, delta) ~ type )
summary(my.fit1)
## Call: survfit(formula = Surv(time, delta) ~ type)
##
## type=1
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 1 52 1 0.981 0.0190 0.944 1.000
## 3 51 2 0.942 0.0323 0.881 1.000
## 4 49 1 0.923 0.0370 0.853 0.998
## 10 48 1 0.904 0.0409 0.827 0.988
## 13 47 2 0.865 0.0473 0.777 0.963
## 16 45 2 0.827 0.0525 0.730 0.936
## 24 43 1 0.808 0.0547 0.707 0.922
## 26 42 1 0.788 0.0566 0.685 0.908
## 27 41 1 0.769 0.0584 0.663 0.893
## 28 40 1 0.750 0.0600 0.641 0.877
## 30 39 2 0.712 0.0628 0.598 0.846
## 32 37 1 0.692 0.0640 0.578 0.830
## 41 36 1 0.673 0.0651 0.557 0.813
## 51 35 1 0.654 0.0660 0.537 0.797
## 65 33 1 0.634 0.0669 0.516 0.780
## 67 32 1 0.614 0.0677 0.495 0.762
## 70 31 1 0.594 0.0683 0.475 0.745
## 72 30 1 0.575 0.0689 0.454 0.727
## 73 29 1 0.555 0.0693 0.434 0.709
## 77 27 1 0.534 0.0697 0.414 0.690
## 91 19 1 0.506 0.0715 0.384 0.667
## 93 18 1 0.478 0.0728 0.355 0.644
## 96 16 1 0.448 0.0741 0.324 0.620
## 100 14 1 0.416 0.0754 0.292 0.594
## 104 12 1 0.381 0.0767 0.257 0.566
## 157 5 1 0.305 0.0918 0.169 0.550
## 167 4 1 0.229 0.0954 0.101 0.518
##
## type=2
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 1 28 1 0.9643 0.0351 0.8979 1.000
## 3 27 1 0.9286 0.0487 0.8379 1.000
## 4 26 1 0.8929 0.0585 0.7853 1.000
## 5 25 2 0.8214 0.0724 0.6911 0.976
## 8 23 1 0.7857 0.0775 0.6475 0.953
## 12 21 1 0.7483 0.0824 0.6031 0.929
## 13 20 1 0.7109 0.0863 0.5603 0.902
## 18 19 1 0.6735 0.0895 0.5190 0.874
## 23 18 1 0.6361 0.0921 0.4790 0.845
## 26 17 1 0.5986 0.0939 0.4402 0.814
## 27 16 1 0.5612 0.0952 0.4025 0.783
## 30 15 1 0.5238 0.0959 0.3658 0.750
## 42 14 1 0.4864 0.0961 0.3302 0.716
## 56 13 1 0.4490 0.0957 0.2956 0.682
## 62 12 1 0.4116 0.0948 0.2621 0.646
## 69 10 1 0.3704 0.0938 0.2255 0.608
## 104 8 2 0.2778 0.0904 0.1468 0.526
## 112 5 1 0.2222 0.0877 0.1025 0.482
## 129 4 1 0.1667 0.0815 0.0639 0.435
## 181 2 1 0.0833 0.0717 0.0155 0.449
detach(tongue)
3.3 累计风险函数
累积风险函数在连续数据上和生存函数的关系如下:
\[S(t)=exp{(-H(t))}\] 风险函数的最大似然估计可以通过kaplan-Meier估计的逆变换得到:
\[\hat{H(t)}=-log\hat{S(t)}\]
尽管survival包中没有函数自动计算任何一个表单,但是通过summary(survfit())返回的对象可以用来计算估计值:
data(tongue)
attach(tongue)
my.surv <- Surv(time[type==1], delta[type==1])
my.fit <- summary(survfit(my.surv~1))
hat.H <- -log(my.fit$surv)
hat.H
## [1] 0.01941809 0.05942342 0.08004271 0.10109612 0.14458123 0.19004360
## [7] 0.21357410 0.23767165 0.26236426 0.28768207 0.34032581 0.36772478
## [13] 0.39589566 0.42488319 0.45565485 0.48740355 0.52019337 0.55409493
## [19] 0.58918625 0.62692657 0.68099379 0.73815221 0.80269073 0.87679870
## [25] 0.96381008 1.18695363 1.47463570
上面得到的结果即是累积风险函数值。下图是累积风险函数值在时间上的走势。
plot(my.fit$time, hat.H, xlab="time", ylab="cumulative hazard value",main="cumulative hazards", type="s")
3.4 多组生存时间比较
survdiff(formula, rho=0)
给定两个或两个以上的样本,生存时间是否存在显著差异? 我们可以结合风险函数假设为:
\[H_0:对所有的t有 h_1(t)=h_2(t)=...=h_n(t)\\H_1:h_i(t_0) \neq h_j(t_0)\] 函数survdiff()用于这个假设检验。第一个参数是生存对象对一个分类的协变量,通常是一个变量指定哪个组对应哪个生存时间。survdiff()返回的对象提供了有用的信息。
data(btrial)
attach(btrial)
survdiff(Surv(time, death)~im)
## Call:
## survdiff(formula = Surv(time, death) ~ im)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## im=1 36 16 20.19 0.869 5.49
## im=2 9 8 3.81 4.599 5.49
##
## Chisq= 5.5 on 1 degrees of freedom, p= 0.02
上面假设检验结果P值等于0.02小于0.05,所以拒绝原假设,两组的样本的生存时间存在显著差异。
第二个survdiff()参数rho,默认等于0,参数rho用来设置权重,如果给生存曲线的第一部分赋予更大的权重,rho大于0。如果更大权重在生存曲线的后面部分,用小于0的。
3.5 Cox比例风险模型
如果考虑其他影响生存时间分布的因素,可以使用Cox回归模型(也叫比例风险模型)它是一种半参数回归模型,利用数学模型拟合生存分布与影响因子之间的关系,评价影响因子对生存函数分布的影响程度。该模型以生存结局和生存时间为应变量,可同时分析众多因素对生存期的影响,能分析带有截尾生存时间的资料,且不要求估计资料的生存分布类型。这里的前提是影响因素的作用不随时间改变,如果不满足这个条件,则应使用含有时间依存协变量的Cox回归模型。
现在我们考虑时间依存的协变量,如干预措施或可能的环境因素在研究过程中产生变化。为了在R中使用时间相关的协变量,我们大量地应用左截断。例如,如果有一个受干预的病人\(i\),然后我们将病人\(i\)分为两个独立的观察:干预前和干预后。假设患者干预发生在\(t_i = 45\)时,观察病人\(i\)在\(t_j = 58\)时发生的事件。我们将这个病人的病历分成两片段:0到45和45到58。干预的协变量可以在每个时间间隔被赋予不同的值。
我们考虑以下示例(模拟数据)。病人记录来自一个诊所。在患者清醒后的150天内获得。数据的事件是酗酒,有四个变量。
* event变量描述被观察的或被删失的复发时间;
* delta变量描述该事件是被观察(真实)还是被删失(假);
* gender变量是一个独立于时间的协变量;
* int是一个依赖时间的协变量表示患者是否进行了干预,有些患者在计划的干预时间前复发,所以里面有NA值。所有干预均发生在患者清醒后。也就是说,只要没有复发,干预协变量随着时间的推移而变化。数据relapse保存在OIsurv包中。
data(relapse)
head(relapse)
## event delta gender inter
## 1 150 FALSE 0 84
## 2 53 TRUE 1 50
## 3 12 TRUE 1 NA
## 4 150 FALSE 0 89
## 5 150 FALSE 1 77
## 6 135 TRUE 1 7
这些数据可以通过以下两个步骤进行建模:
1. 构造生存记录,包括左和右删失的观察。每一个接受干预的病人的生存记录分成两个:干预前和干预后。
2. 新的生存记录可以通过coxph()函数运行。
我们首先初始化新的向量来表示生存记录的变量。变量\(t_1\)和 \(t_2\)分别表示开始和结束时间,\(d\)分别表示观察到复发(TRUE)还是右删失(FALSE),\(g\)代表性别,\(i\)代表病人是否正在经历干预治疗。
N <- dim(relapse)[1]
t1 <- rep(0, N+sum(!is.na(relapse$int)))# 初始化开始时间
t2 <- rep(-1, length(t1))# 构建结束时间向量
d <- rep(-1, length(t1))# 时间是否删失
g <- rep(-1, length(t1))# 性别
i <- rep(FALSE, length(t1))# 初始化干预为FALSE
接下来,对每个病人进行检查,看他们是否接受过干预。如果他们没有,那么他们的记录被简单地复制到所使用的新变量中。如果有干预,他们的观察分为两部分:干预前和干预后。
j <- 1
for(ii in 1:dim(relapse)[1]){
if(is.na(relapse$int[ii])){
t2[j] <- relapse$event[ii]
d[j] <- relapse$delta[ii]
g[j] <- relapse$gender[ii]
j <- j+1
} else { # 若发生干预事件,分裂数据
g[j+0:1] <- relapse$gender[ii] # 每个时间的性别是一样的
d[j] <- 0 # 干预前没有复发
d[j+1] <- relapse$delta[ii] # 干预后是否复发
i[j+1] <- TRUE # 干预事件发生
t2[j] <- relapse$int[ii]-1 # 干预前的结束时间
t1[j+1] <- relapse$int[ii]-1 # 干预后的开始时间
t2[j+1] <- relapse$event[ii] # 干预后的结束时间
j <- j+2 # 两条记录增加
}
}
虽然许多患者有随时间变化的协变量(即干预),但每一个新的间隔存在不改变的协变量,生存对象就可以建立Cox PH模型。
mySurv <- Surv(t1, t2, d)# 构建生存对象
myCPH <- coxph(mySurv ~ g + i)# cox回归
test.ph <- cox.zph(myCPH)
test.ph
## rho chisq p
## g -0.0325 0.216 0.642
## iTRUE -0.0225 0.109 0.741
## GLOBAL NA 0.314 0.855
使用R语言survival包中的cox.zph函数进行模型诊断,从上面的结果可以看出,两个变量的P值都大于0.05,说明每个变量均满足PH检验,而模型的整体检验P值0.855,模型整体满足PH检验。
这个例子是一个简单的例子;只有一个与时间相关的协变量,它每次最多改变一次。在某些情况下,可能有许多随时间变化的协变量。在大多数情况下,可以使用相同的方法。每当一个协变量从一个时间单位变化到下一个时间单位,将间隔分割为两个,根据需要使用左截断和右删失。
4 本章汇总
名称 | 类别 | 功能 |
---|---|---|
KMsurv | 包 | 生存分析的数据集 |
survival | 包 | 生存分析包 |
OIsurv | 包 | 生存分析包 |
data() | 函数 | 载入R内置数据集 |
head() | 函数 | 查看数据 |
attach() | 函数 | 减少重复键入表名 |
detach() | 函数 | 解除attach函数 |
Surv() | 函数 | 生存对象产生 |
survfit() | 函数 | 生存函数拟合 |
survdiff() | 函数 | 生存时间比较 |
dim(relapse) | 函数 | 得到对象的大小 |
rep() | 函数 | 将对象重复多少次产生新对象 |
is.na() | 函数 | 判断是否为空值 |
coxph() | 函数 | cox回归拟合函数 |
cox.zph() | 函数 | cox回归模型检验 |
5 参考文献
[1] 彭非. 生存分析[M]. 中国人民大学出版社, 2004.
[2] http://www.openintro.org/stat/down/Survival-Analysis-in-R.pdf