生存分析

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