翻译自：http://www.sthda.com/english/wiki/survival-analysis-basics

生存分析基础

生存分析对应于一组统计方法，用于调查感兴趣事件发生所花费的时间。

生存分析可用于许多领域，例如：

• 用于患者生存时间分析的癌症研究，
• “事件历史分析”的社会学，
• 在工程中用于“故障时间分析”。

在癌症研究中，典型的研究问题如下：

• 某些临床特征对患者生存有何影响？
• 一个人生存3年的概率是多少？
• 两组患者的生存率是否存在差异？

内容

目标

本章的目的是描述生存分析的基本概念。在癌症研究中，大多数生存分析使用以下方法：

• Kaplan-Meier图 （Kaplan-Meier plots）可视化生存曲线
• 对数秩检验 （Log-rank test）以比较两组或更多组的生存曲线
• 用Cox比例风险回归（Cox proportional hazards regression）描述变量对生存的影响。下一章将讨论Cox模型：Cox比例风险模型。

在这里，我们将从解释生存分析的基本概念开始，包括：

• 如何生成和解释生存曲线，
• 以及如何量化和测试两组或更多组患者之间的生存差异。

然后，我们将继续使用Cox比例风险模型描述多元分析。

基本概念

在这里，我们首先定义生存分析的基本术语，包括：

• 生存时间和事件（Survival time and event）
• 删失（Censoring）
• 生存函数和风险函数（Survival function and hazard function）

癌症研究中的生存时间和事件类型（Survival time and event）

有不同类型的事件，包括：

• 复发（Relapse）
• 进展（Progression）
• 死亡（Death）

从“响应到治疗”（完全缓解）到发生关注事件的时间通常称为生存时间（或事件发生时间）。

癌症研究中两个最重要的标度包括：i）死亡时间；ii）无复发生存时间，对应于对治疗的反应与疾病复发之间的时间。也称为无病生存时间（disease-free survival time）和无事件生存时间（event-free survival time）。

删失（Censoring）

如上所述，生存分析着眼于直到发生感兴趣事件（复发或死亡）之前的预期持续时间。但是，在研究期间内某些人可能未观察到该事件，从而产生了所谓的删失（Censoring）现象。

审查可能以下列方式出现：

1. 患者在研究时间段内尚未（尚未）经历感兴趣的事件，例如复发或死亡；
2. 研究期间患者失去随访；
3. 患者经历了另一种事件，因此无法进行进一步的随访。

在生存分析中会遇到这种现象，称为右侧删失。

这里补充一下右侧删失和左侧删失的意思：
以右侧为例，当患者发生上述情况时，在时间轴这个时间点的右侧（即该时间点之后）数据点是缺失的，因此称为右侧删失。临床床上经常遇到的是右侧删失。
这样左侧删失也容易理解了，既被随访者由于某些原因在时间轴内某一点之前没能参与随访，因此在改时间点之前（既时间轴左侧）数据是缺失的，因此称为左侧删失。

比如研究者想跟踪调查青少年12岁至18岁之前的视力变化，如果某个被调查者在14岁才开始进行随访就会产生左侧缺失，如果某个被调查者在14岁由于玩游戏过度而住院无法继续参与随访就会产生右侧缺失。

生存和风险函数

使用两个相关的概率来描述生存数据：生存概率和危险概率。

生存函数，也被称为幸存者函数S(t) 是从时间起源（例如癌症诊断）到指定的未来时间t生存的概率。

风险函数，记为h(t)，是在时间t被观察的个人发生某项事件的概率。

注意，生存函数侧重于没有事件，危害函数着重于事件发生。

Kaplan-Meier生存估计

Kaplan-Meier（KM）方法是一种非参数方法，用于根据观察到的生存时间估算生存概率（Kaplan和Meier，1958年）。

时间点 t(i) 的生存概率 计算如下：

估计的概率

估计概率（S（t））是仅在每个事件发生时改变值的阶跃函数。也可以计算生存概率的置信区间。
KM生存曲线是KM生存概率与时间的关系曲线，它提供了一个有用的数据摘要，可用于估计中位生存时间等指标。

R中的生存分析

安装并加载所需的R包

我们将使用两个R包：

• survival 用于计算存活分析

• survminer 用于总结和可视化生存分析结果

• 安装软件包

install.packages(c("survival", "survminer"))

• 加载包

library("survival")
library("survminer")

示例数据集

我们将使用生存包中提供的肺癌数据。

data("lung")
head(lung)

  inst time status age sex ph.ecog ph.karno pat.karno meal.cal wt.loss
1    3  306      2  74   1       1       90       100     1175      NA
2    3  455      2  68   1       0       90        90     1225      15
3    3 1010      1  56   1       0       90        90       NA      15
4    5  210      2  57   1       1       90        60     1150      11
5    1  883      2  60   1       0      100        90       NA       0
6   12 1022      1  74   1       1       50        80      513       0

• 机构：机构代码
• 时间：以天为单位的生存时间
• 状态：审查状态1 =审查，2 =失效
• 年龄：岁
• 性别：男= 1女= 2
• ph.ecog：ECOG成绩得分（0 =好5 =死）
• ph.karno：医师对Karnofsky成绩评分（差= 0-良好= 100）
• pat.karno：Karnofsky表现评分，按患者评分
• meal.cal：用餐时消耗的卡路里
• wt.loss：最近六个月的体重减轻

计算生存曲线：survfit（）

我们想按性别计算生存率。

生存软件包中的功能survfit（）可用于计算kaplan-Meier生存估计。其主要论​​点包括：

• 使用Surv（）函数创建的生存对象
• 以及包含变量的数据集。

要计算生存曲线，请键入：

fit <- survfit(Surv(time, status) ~ sex, data = lung)
print(fit)

Call: survfit(formula = Surv(time, status) ~ sex, data = lung)
        n events median 0.95LCL 0.95UCL
sex=1 138    112    270     212     310
sex=2  90     53    426     348     550

默认情况下，函数print（）显示生存曲线的简短摘要。它打印观察值，事件数，中位生存率和中位值的置信度限制。

如果要显示生存曲线的更完整摘要，请键入以下内容：

# Summary of survival curves
summary(fit)
# Access to the sort summary table
summary(fit)$table

访问survfit（）返回的值

函数survfit（）返回变量列表，包括以下组件：

• n：每条曲线中的对象总数。
• 时间：曲线上的时间点。
• n。风险：时间t处有风险的受试者人数
• n.event：在时间t发生的事件数。
• n。审查者：在时间t退出事件而不发生风险的被审查者的数量。
• 下，上：曲线的置信度上限和下限。
• 分层：表示曲线估计的分层。如果地层不为NULL，则结果中有多条曲线。层次的水平（一个因子）是曲线的标签。

可以按以下方式访问组件：

d <- data.frame(time = fit$time,
                  n.risk = fit$n.risk,
                  n.event = fit$n.event,
                  n.censor = fit$n.censor,
                  surv = fit$surv,
                  upper = fit$upper,
                  lower = fit$lower
                  )
head(d)

  time n.risk n.event n.censor      surv     upper     lower
1   11    138       3        0 0.9782609 1.0000000 0.9542301
2   12    135       1        0 0.9710145 0.9994124 0.9434235
3   13    134       2        0 0.9565217 0.9911586 0.9230952
4   15    132       1        0 0.9492754 0.9866017 0.9133612
5   26    131       1        0 0.9420290 0.9818365 0.9038355
6   30    130       1        0 0.9347826 0.9768989 0.8944820

可视化生存曲线

我们将使用功能ggsurvplot（）[在Survminer R包中]生成两组受试者的生存曲线。

也可以显示：

• 使用参数conf.int = TRUE对生存函数的95％置信度限制。
• 使用option risk.table按时间划分处于风险中的个体的数量和/或百分比。risk.table的允许值包括：
  • TRUE或FALSE指定是否显示风险表。默认值为FALSE。
  • “绝对”或“百分比”：分别显示按时间划分处于风险中的受试者的绝对数和百分比。使用“ abs_pct”显示绝对数字和百分比。
• 对数秩检验的p值，使用pval = TRUE比较组。
• 使用参数surv.median.line在中值生存时的水平/垂直线。允许的值包括c（“ none”，“ hv”，“ h”，“ v”）之一。v：垂直，h：水平。

# Change color, linetype by strata, risk.table color by strata
ggsurvplot(fit,
          pval = TRUE, conf.int = TRUE,
          risk.table = TRUE, # Add risk table
          risk.table.col = "strata", # Change risk table color by groups
          linetype = "strata", # Change line type by groups
          surv.median.line = "hv", # Specify median survival
          ggtheme = theme_bw(), # Change ggplot2 theme
          palette = c("#E7B800", "#2E9FDF"))

可以使用以下参数进一步定制该图：

• conf.int.style =“步骤”以更改置信区间带的样式。
• xlab更改x轴标签。
• break.time.by = 200以时间间隔将x轴断开200。
• risk.table =“ abs_pct”，以显示处于风险中的个人的绝对数量和百分比。
• risk.table.y.text.col = TRUE，risk.table.y.text = FALSE在风险表图例的文本注释中提供条形而不是名称。
• ncensor.plot = TRUE，以绘制时间t处受检对象的数量。正如Marcin Kosinski所建议的那样，这是对生存曲线的一个很好的附加反馈，因此人们可以认识到：生存曲线看起来如何，风险集的数量是多少，以及风险集变小的原因是什么？是由事件还是受审查事件引起的？
• legend.labs更改图例标签。

ggsurvplot(
   fit,                     # survfit object with calculated statistics.
   pval = TRUE,             # show p-value of log-rank test.
   conf.int = TRUE,         # show confidence intervals for 
                            # point estimaes of survival curves.
   conf.int.style = "step",  # customize style of confidence intervals
   xlab = "Time in days",   # customize X axis label.
   break.time.by = 200,     # break X axis in time intervals by 200.
   ggtheme = theme_light(), # customize plot and risk table with a theme.
   risk.table = "abs_pct",  # absolute number and percentage at risk.
  risk.table.y.text.col = T,# colour risk table text annotations.
  risk.table.y.text = FALSE,# show bars instead of names in text annotations
                            # in legend of risk table.
  ncensor.plot = TRUE,      # plot the number of censored subjects at time t
  surv.median.line = "hv",  # add the median survival pointer.
  legend.labs = 
    c("Male", "Female"),    # change legend labels.
  palette = 
    c("#E7B800", "#2E9FDF") # custom color palettes.
)

生存分析

Kaplan-Meier图可以解释如下：

横轴（x轴）表示以天为单位的时间，纵轴（y轴）表示生存的可能性或生存的人口比例。线代表两组的存活曲线。曲线中的垂直下降表示事件。曲线上的垂直刻度线表示此时已对患者进行检查。

• 在零时，生存概率为1.0（或100％的参与者还活着）。
• 在时间250，性别= 1的存活概率约为0.55（或55％），性别= 2的存活概率约为0.75（或75％）。
• 性别= 2的中位生存期约为270天，性别= 2的中位生存期约为426天，这表明性别= 2的生存期高于性别= 1

可以使用以下代码获得每组的中位生存时间：

summary(fit)$table

      records n.max n.start events   *rmean *se(rmean) median 0.95LCL 0.95UCL
sex=1     138   138     138    112 325.0663   22.59845    270     212     310
sex=2      90    90      90     53 458.2757   33.78530    426     348     550

每组的中位生存时间代表生存概率S（t）为0.5的时间。

性别= 1（男性）的中位生存时间为270天，而性别= 2（女性）则为426天。与男性相比，女性肺癌似乎具有生存优势。但是，要评估这种差异是否具有统计学显着性，需要进行正式的统计检验，这将在下一节中讨论。

注意，置信极限在曲线的尾部很宽，因此很难进行有意义的解释。这可以通过以下事实来解释：在实践中，通常有一些患者在随访结束时迷失了随访或存活。因此，在随访结束之前在x轴上缩短图可能是明智的（Pocock等，2002）。

可以使用参数xlim缩短生存曲线，如下所示：

ggsurvplot(fit,
          conf.int = TRUE,
          risk.table.col = "strata", # Change risk table color by groups
          ggtheme = theme_bw(), # Change ggplot2 theme
          palette = c("#E7B800", "#2E9FDF"),
          xlim = c(0, 600))

生存分析

注意，可以使用参数fun指定三个经常使用的转换：

• “ log”：幸存者功能的对数转换，
• “事件”：绘制累积事件（f（y）= 1-y）。也称为累积发生率
• “ cumhaz”绘制了累积危害函数（f（y）= -log（y））

例如，要绘制累积事件，请键入：

ggsurvplot(fit,
          conf.int = TRUE,
          risk.table.col = "strata", # Change risk table color by groups
          ggtheme = theme_bw(), # Change ggplot2 theme
          palette = c("#E7B800", "#2E9FDF"),
          fun = "event")

生存分析

累积性危险是常用来估计危险概率。定义为：

H(t)=−log(survivalfunction)=−log(S(t))

累积危险（H（t））可以解释为死亡的累积力。换言之，如果事件是可重复的过程，则它对应于时间t时每个个体预期的事件数。

要绘制累积危害，请键入以下内容：

ggsurvplot(fit,
          conf.int = TRUE,
          risk.table.col = "strata", # Change risk table color by groups
          ggtheme = theme_bw(), # Change ggplot2 theme
          palette = c("#E7B800", "#2E9FDF"),
          fun = "cumhaz")

生存分析

Kaplan-Meier生命表：生存曲线摘要

如上所述，您可以使用函数summary（）来获得生存曲线的完整摘要：

summary(fit)

还可以使用功能surv_summary（）[在survminer程序包中]获取生存曲线的摘要。与默认的summary（）函数相比，surv_summary（）创建一个数据框，其中包含来自survfit结果的不错的摘要。

res.sum <- surv_summary(fit)
head(res.sum)

  time n.risk n.event n.censor      surv    std.err     upper     lower strata sex
1   11    138       3        0 0.9782609 0.01268978 1.0000000 0.9542301  sex=1   1
2   12    135       1        0 0.9710145 0.01470747 0.9994124 0.9434235  sex=1   1
3   13    134       2        0 0.9565217 0.01814885 0.9911586 0.9230952  sex=1   1
4   15    132       1        0 0.9492754 0.01967768 0.9866017 0.9133612  sex=1   1
5   26    131       1        0 0.9420290 0.02111708 0.9818365 0.9038355  sex=1   1
6   30    130       1        0 0.9347826 0.02248469 0.9768989 0.8944820  sex=1   1

函数surv_summary（）返回包含以下列的数据帧：

• 时间：曲线有阶跃的时间点。
• n。风险：处于t风险的受试者人数。
• n.event：在时间t发生的事件数。
• n.censor：审查事件的数量。
• surv：估计生存概率。
• std.err：生存标准误差。
• upper：置信区间的上限
• 较低：置信区间的下限
• 分层：表示曲线估计的分层。层次的水平（一个因子）是曲线的标签。

在生存曲线已拟合一个或多个变量的情况下，surv_summary对象包含表示变量的额外列。这使得可以按层次或某些因素组合来考虑ggsurvplot的输出。

surv_summary对象还有一个名为“表”的属性，其中包含有关生存曲线的信息，包括具有置信区间的生存中位数，以及受试者总数和每条曲线中的事件数。要访问属性“表”，请输入以下命令：

attr(res.sum, "table")

比较生存曲线的Log-Rank检验：survdiff（）

log-rank test 是比较两条或更多条生存曲线的最广泛使用的方法。零假设是两组之间的生存率没有差异。对数秩检验是一种非参数检验，它不对生存分布做出任何假设。本质上，对数秩检验将每个组中观察到的事件数与如果原假设为真（即，生存曲线相同）时的预期事件数进行比较。对数秩统计量近似作为卡方检验统计量分布。

函数survdiff（）[在生存包中]可用于比较两个或更多生存曲线的对数秩检验。

survdiff（）可以如下使用：

surv_diff <- survdiff(Surv(time, status) ~ sex, data = lung)
surv_diff

Call:
survdiff(formula = Surv(time, status) ~ sex, data = lung)
        N Observed Expected (O-E)^2/E (O-E)^2/V
sex=1 138      112     91.6      4.55      10.3
sex=2  90       53     73.4      5.68      10.3
 Chisq= 10.3  on 1 degrees of freedom, p= 0.00131 

该函数返回组件列表，包括：

• n：每组的科目数。
• obs：每组中事件的加权观测数量。
• exp：每个组中加权的预期事件数。
• chisq：用于检验相等性的卡方统计量。
• 阶层：（可选）每个阶层中包含的主题数。

生存差异的对数秩检验给出p值为p = 0.0013，表明性别群体的生存差异显着。

拟合复杂的生存曲线

在本节中，我们将使用多个因素的组合来计算生存曲线。接下来，我们将结合多种因素来研究ggsurvplot（）的输出

1. 使用结肠数据集拟合（复杂）生存曲线

require("survival")
fit2 <- survfit( Surv(time, status) ~ sex + rx + adhere,
                data = colon )

2. 使用survminer可视化输出。下图显示了根据rx＆粘附力值按性别分面的生存曲线。

# Plot survival curves by sex and facet by rx and adhere
ggsurv <- ggsurvplot(fit2, fun = "event", conf.int = TRUE,
                     ggtheme = theme_bw())

ggsurv$plot +theme_bw() + 
  theme (legend.position = "right")+
  facet_grid(rx ~ adhere)

概要

生存分析是用于数据分析的一组统计方法，其中感兴趣的结果变量是事件发生之前的时间。

生存数据通常根据两个相关功能进行描述和建模：

• 幸存者函数代表一个人从起源时间到超过时间t的某个时间生存的概率。通常用Kaplan-Meier方法估算。对数秩检验可用于测试组（例如治疗组）的生存曲线之间的差异。

• 危害函数给出了一次事件的瞬时可能性，并给出了到那个时间的生存率。它主要用作诊断工具或用于指定生存分析的数学模型。

在本文中，我们演示了如何使用两个R包（survival（用于分析）和 survminer（用于可视化））的组合来执行和可视化生存分析。