如何用R的ggsurvplot绘制带风险表、置信区间的累积发生率图

问题描述 #

需在R中绘制包含指定时间点风险表、置信区间、删失标记的累积发生率图。现有数据包含天数（time）、事件状态（status，T=事件发生，F=未发生，含NA值）。使用ggsurvfit可生成带所需元素的生存曲线（无事件比例），但需将Y轴转换为累积发生率；尝试base R绘图或cmprsk::cuminc时，添加风险表参数会报错。

解决方案 #

方法1：基于ggsurvfit转换Y轴（推荐，单一事件场景） #

利用ggsurvfit的绘图体系，通过Y轴转换将生存曲线转为累积发生率，同时保留所有所需元素：

# 加载依赖包
library(survival)
library(tidyverse)
library(ggsurvfit)

# 模拟含NA的数据集
time = c(1:10, NA)
status = c(T, F, T, F, T, F, T, F, T, F, NA)
df <- data.frame(time, status)

# 拟合生存模型
ci_surv <- Surv(time, status)
ci_survfit <- survfit2(ci_surv ~ 1)

# 绘制累积发生率图
ggsurvfit(ci_survfit) +
  labs(x = "Days", y = "Cumulative Incidence") +
  scale_x_continuous(breaks = c(2, 5, 10)) +
  # 转换Y轴：累积发生率 = 1 - 无事件比例
  scale_y_continuous(
    trans = "reverse",  # 反转原生存曲线的Y轴方向
    labels = ~ 1 - .x,  # 将刻度标签转换为累积发生率
    breaks = seq(0, 1, 0.2)
  ) +
  add_confidence_interval() +
  add_risktable(
    risktable_title = "Number at risk",
    times = c(2, 5, 10)  # 指定风险表显示的时间点
  ) +
  add_censor_mark()

方法2：基于cuminc+tidycmprsk（竞争风险场景） #

若后续涉及多种事件类型的竞争风险分析，可使用cmprsk::cuminc结合tidycmprsk绘图：

library(tidycmprsk)
library(cmprsk)

# 将状态转换为数值型（1=事件发生，0=删失）
df$status_num <- ifelse(df$status, 1, 0)

# 拟合竞争风险模型
ci_cuminc <- cuminc(Surv(time, status_num) ~ 1, data = df)

# 绘图并添加所需元素
ggcuminc(ci_cuminc) +
  labs(x = "Days", y = "Cumulative Incidence") +
  scale_x_continuous(breaks = c(2, 5, 10)) +
  add_confidence_interval() +
  add_risktable(times = c(2, 5, 10)) +
  add_censor_mark()

错误原因说明 #

• base R的plot.survfit函数不支持atrisk参数，该参数属于survminer包的ggsurvplot，因此直接用base plot会报错。
• 单一事件场景下，cuminc虽能生成累积发生率，但没必要，直接转换生存曲线更高效；若强行使用，需搭配tidycmprsk的ggcuminc而非base plot来添加风险表。

内容的提问来源于stack exchange，提问作者Tim