R语言生存分析数据分析可视化案例

目标

本文的目的是对如何在R中进行生存分析进行简短而全面的评估。关于该主题的文献很广泛，仅涉及有限数量的（常见）问题。最近我们被客户要求撰写关于生存分析的研究报告，包括一些图形和统计输出。
可用的R包数量反映了对该主题的研究范围。 

---

 

 视频：R语言生存分析原理与晚期肺癌患者分析案例

R包

可以使用各种R包来解决特定问题。以下是本次用于读取，管理，分析和显示数据的软件包。
运行以下行以安装和加载所需的包。

if (!require(pacman)) install.packages("pacman")
pacman::p_load(tidyverse, survival  )

  

数据

该评价将基于orca数据集，数据集包含1985年1月1日至2005年12月31日期间芬兰最北部省份诊断为口腔鳞状细胞癌（OSCC）的338名患者的一部分。患者的随访始于癌症诊断之日，并于2008年12月31日死亡，迁移或随访截止日期结束。死亡原因分为两类：（1） ）OSCC死亡; （2）其他原因造成的死亡。
数据集包含以下变量：
id=序号，
sex=性别，类别1 =“女性”，2 =“男性”，
age=诊断癌症日期的年龄（年），
stage=肿瘤的TNM分期（因子）：1 =“I”，...， 4 =“IV”，5 =“unkn” 
time=自诊断至死亡或审查的随访时间（以年为单位），
event=结束随访的事件（因子）：1 =正常，2 =口腔癌死亡， 3 =其他原因造成的死亡。

 将数据从URL加载到R中。

head(orca)

  id    sex      age stage  time          event
1  1   Male 65.42274  unkn 5.081          Alive
2  2 Female 83.08783   III 0.419 Oral ca. death
3  3   Male 52.59008    II 7.915    Other death
4  4   Male 77.08630     I 2.480    Other death
5  5   Male 80.33622    IV 2.500 Oral ca. death
6  6 Female 82.58132    IV 0.167    Other death

summary(orca)

       id             sex           age         stage         time                   event    
 Min.   :  1.00   Female:152   Min.   :15.15   I   :50   Min.   : 0.085   Alive         :109  
 1st Qu.: 85.25   Male  :186   1st Qu.:53.24   II  :77   1st Qu.: 1.333   Oral ca. death:122  
 Median :169.50                Median :64.86   III :72   Median : 3.869   Other death   :107  
 Mean   :169.50                Mean   :63.51   IV  :68   Mean   : 5.662                       
 3rd Qu.:253.75                3rd Qu.:74.29   unkn:71   3rd Qu.: 8.417                       
 Max.   :338.00                Max.   :92.24             Max.   :23.258                       

生存数据分析

生存分析侧重于事件数据的时间。在我们的例子中，是诊断后的死亡时间。

为了定义失效时间随机变量，我们需要：
1。时间起源（诊断OSCC），
2。时间尺度（诊断后的年数，年龄），
3。事件的定义。我们将首先考虑总死亡率 。

图1：转换的框图。

         Alive Oral ca. death    Other death 
           109            122            107 

 
FALSE  TRUE 
  109   229 

以图形方式显示观察到的随访时间对于生存数据的分析非常有帮助。 

 

OSCC死亡更有可能在诊断后早期发生，而不是其他原因引起的死亡。类型怎么样？

 'Surv' num [1:338, 1:2]  5.081+  0.419   7.915   2.480   2.500   0.167   5.925+  1.503  13.333   7.666+ ...
 - attr(*, "dimnames")=List of 2
  ..$ : NULL
  ..$ : chr [1:2] "time" "status"
 - attr(*, "type")= chr "right"

然后将创建的生存对象用作生存分析的其他特定函数中的因变量。

---

估计生存函数

非参数估计

我们将首先介绍一类非参数估计 。

Kaplan–Meier 

生存曲线基于每个死亡时间的风险数量和事件数量。包的survfit()创建（估计）生存曲线 。 

Call: survfit(formula = Surv(time, all) ~ 1, data = orca)
 
         n     events     *rmean *se(rmean)     median    0.95LCL    0.95UCL 
   338.000    229.000      8.060      0.465      5.418      4.331      6.916 
    * restricted mean with upper limit =  23.3 

函数返回估计的生存曲线的摘要。 

   time n.risk n.event n.censor      surv     std.err     upper     lower
1 0.085    338       2        0 0.9940828 0.004196498 1.0000000 0.9859401
2 0.162    336       2        0 0.9881657 0.005952486 0.9997618 0.9767041
3 0.167    334       4        0 0.9763314 0.008468952 0.9926726 0.9602592
4 0.170    330       2        0 0.9704142 0.009497400 0.9886472 0.9525175
5 0.246    328       1        0 0.9674556 0.009976176 0.9865584 0.9487228
6 0.249    327       1        0 0.9644970 0.010435745 0.9844277 0.9449699

ggsurvplot()的survminer提供了估计的生存曲线的信息性说明。

 

默认的KM图表显示了生存函数。

生存曲线估算

生存曲线在精算师和人口统计学中非常普遍。它特别适用于分组数据。

为了在实际示例中显示此方法，我们首先需要创建聚合数据，即将后续分组并在每个层中计算风险。

基于分组的数据，我们估计会用生存曲线。 

      nsubs nlost nrisk nevent   surv    pdf hazard se.surv se.pdf se.hazard
0-1     338     0 338.0     64 1.0000 0.1893 0.2092  0.0000 0.0213    0.0260
1-2     274     4 272.0     41 0.8107 0.1222 0.1630  0.0213 0.0179    0.0254
2-3     229     9 224.5     21 0.6885 0.0644 0.0981  0.0252 0.0136    0.0214
3-4     199    12 193.0     20 0.6241 0.0647 0.1093  0.0265 0.0140    0.0244
4-5     167     9 162.5     13 0.5594 0.0448 0.0833  0.0274 0.0121    0.0231
5-6     145    14 138.0     13 0.5146 0.0485 0.0989  0.0279 0.0131    0.0274
6-7     118     5 115.5      8 0.4662 0.0323 0.0717  0.0283 0.0112    0.0254
7-8     105     8 101.0      9 0.4339 0.0387 0.0933  0.0286 0.0126    0.0311
8-9      88     7  84.5      1 0.3952 0.0047 0.0119  0.0288 0.0047    0.0119
9-10     80     4  78.0      8 0.3905 0.0401 0.1081  0.0288 0.0137    0.0382
10-11    68     4  66.0      5 0.3505 0.0266 0.0787  0.0291 0.0116    0.0352

Nelson-Aalen估计

  

图形比较

可以绘制不同的生存函数估计值来评估潜在的差异。