R 语言做竞争风险回归分析

喜欢我就快点关注我吧

1引言

此篇内容来自粉丝的投稿,感谢这位粉丝的分享!此外 俊粉投稿分享群 已成立,欢迎愿意分享的有志之士加入。

2内容

说到生存分析,相必大家想到的第一个都是 COX 回归和 KM 曲线,不过相关的内容曾老师已经分享过很多了,如果有需求可见曾老师这个专栏:生存分析 。

在生存分析的时候,我们一般只在意一个结局,例如患者死亡或者疾病复发。传统的生存模型认为删失事件和结局事件独立,两者不存在竞争风险 但是临床生存数据常常有多个结局,如果存在某种事件可能会影响结局事件发生的概率或者是完全阻碍其发生,那么两者就存在竞争风险,而此时如果我们将竞争事件作为删失事件进行传统的单因素回归就可能高估了累计生存率

因此:

  • 只有一个结局事件:用 COX — 肿瘤生存时间
  • 有多个竞争结局事件:用竞争风险模型 —慢性肾病患者死亡与透析,心肌梗死患者导致的死亡与其他死因

R 语言如何完成竞争风险回归?

进行竞争风险回归的核心思想就是关注一种结局,而排除另外一种的影响,目前有两种方法能满足这种需求

  • cause-specific hazards :即直接对感兴趣的结局事件单独计算生存率
  • Subdistribution hazard : 在考虑竞争事件的情况下计算结局事件的生存率

从统计学上考虑,cause-specific hazards 方法更倾向于计算单个因素对预后或治疗效果的影响;而 Subdistribution hazard 方法能更好构建预后列线图,评估对健康的影响(其实好像论文中 Subdistribution hazard 用的多一点) 。

下面直接上例子!我们使用 MASS 包的 Melanoma 数据:

  • time survival time in days, possibly censored.
  • status 1 died from melanoma, 2 alive, 3 dead from other causes.
  • sex 1 = male, 0 = female.
  • age age in years.
  • year of operation.
  • thickness tumor thickness in mm.
  • ulcer 1 = presence, 0 = absence.

能看到它有一个结局事件:患黑色素瘤死亡,一个竞争时间:死于其他原因,一个删失事件:仍存活

我们先稍微调整一下它 status 的编码方式,一般 0 是删失,1 是结局事件,2 是竞争事件:


1.Subdistribution hazard 方法:

(1)单因素计算累计风险:

在这个例子中,如果我们要估计患者在排除竞争事件后的患黑色素死亡率,应该怎么做呢?

我们使用 tidycmprsk 包中的 cuminc  函数计算竞争风险背景下的累计死亡率:


library(tidycmprsk)
library(gtsummary)
library(ggsurvfit)
cuminc(Surv(time, status) ~ 1, data = Melanoma)

然后我们可以用 ggsurvfit 包中的 ggcuminc() 函数可视化,默认以 status = 1 为结局事件:


cuminc(Surv(time, status) ~ 1, data = Melanoma) %>%
  ggcuminc() +
  labs(x = "Days") +
  add_confidence_interval() +
  add_risktable()

如果我们想同时查看黑色素瘤的累计死亡和其他原因的累计死亡,我们可以使用 outcome 参数 ggcuminc(outcome=):


cuminc(Surv(time, status) ~ 1, data = Melanoma) %>%
  ggcuminc(outcome = c("1""2")) +
  ylim(c(01)) +
  labs(x = "Days")

如果我们想看变量对累计生存的影响:tbl_cuminc 来自于 tidycmprsk 包:


cuminc(Surv(time, status) ~ ulcer, data = Melanoma) %>%
  tbl_cuminc(
    times = 1826.25,
    label_header = "**{time/365.25}-year cuminc**") %>%
  add_p()

也可以用 ggcuminc 可视化:


cuminc(Surv(time, status) ~ ulcer, data = Melanoma) %>%
  ggcuminc() +
  labs( x = "Days") +
  add_confidence_interval() +
  add_risktable()

(2)多因素回归计算 HR

使用 tidycmprsk 包中的 crr  函数进行对性别和年龄进行累计风险回归,评估这俩变量在排除竞争事件影响后与结局事件关系:


crr(Surv(time, status) ~ sex + age, data = Melanoma)

能看到,Sex 是显著的危险因素,而 age 无显著影响。

使用 gtsummary 包让结果更加简洁:


crr(Surv(time, status) ~ sex + age, data = Melanoma) %>%
  tbl_regression(exp = TRUE)

结语: 最常见的生存分析里也有着大学问!基于数据的类型选择最适合的统计方法,就好像从纷乱的算法世界中找到打开满载生物学意义的宝箱的钥匙,当你拿对了钥匙,才能一窥自然的奥秘。

参考文献:

  • 1.https://www.emilyzabor.com/tutorials/survival_analysis_in_r_tutorial.html
  • 2.聂志强, 欧艳秋, 曲艳吉, 袁海云, 刘小清. 临床生存数据新视角:竞争风险模型[J]. 中华流行病学杂志, 2017, 38(8): 1127-1131

3结尾

欢迎来投稿!


欢迎加入生信交流群。加我微信我也拉你进 微信群聊 老俊俊生信交流群 (微信交流群需收取 20 元入群费用,一旦交费,拒不退还!(防止骗子和便于管理)) 。QQ 群可免费加入, 记得进群按格式修改备注哦。

老俊俊微信:

知识星球:


往期回顾目录

听说你还想添加 KEGG 注释? BioSeqUtils 延伸序列及提取启动子序列 ClusterGVis 对接单细胞啦 About BioSeqUtils 基因特征序列提取 R 包 BioSeqUtils 关于 chatGPT 无法连接的问题 问问 chatGPT 来提取最长转录本及基因序列 GseaVis 一键对接 GSEA 软件结果并可视化 同学你又在画 GSEA? R 语言中的 S3/S4 类

声明:文中观点不代表本站立场。本文传送门:https://eyangzhen.com/190552.html

(0)
联系我们
联系我们
分享本页
返回顶部