如何用R转换生存分析累积风险数据以适配Excel阶梯图绘制

我正在为客户开展项目，要求所有图表在Excel中制作。其中需绘制累积风险图，数据通过Survival::survfit函数获取。问题在于Excel无法直接绘制阶梯图，需转换数据：让每个时间点出现两次，一次对应前一个累积风险值，一次对应当前时间点的累积风险值。手动在Excel操作虽简单但繁琐耗时，作为R语言新手，尚未找到解决方法。

现有代码与原始数据 #

# Load survival package
library(survival)

# Create survfit object
Survival_Function <- survfit(Surv(lung$time,
                                 lung$status == 2)~1)

# Extract cumulative hazards
cumhaz <- data.frame(Survival_Function$time, Survival_Function$cumhaz)

head(cumhaz)

运行后得到的原始数据：

Survival_Function.time Survival_Function.cumhaz
1                      5              0.004385965
2                     11              0.017601824
3                     12              0.022066110
4                     13              0.031034720
5                     15              0.035559606
6                     26              0.040105061

目标数据格式（适配Excel阶梯图） #

Survival_Function.time Survival_Function.cumhaz
1                      5              0.004385965
2                     11              0.004385965
3                     11              0.017601824
4                     12              0.017601824
5                     12              0.022066110
6                     13              0.022066110
7                     13              0.031034720
8                     15              0.031034720
9                     15              0.035559606
10                    26              0.035559606
11                    26              0.040105061

解决方案 #

方法1：基础R实现（适合新手理解）

无需额外安装包，用基础循环即可完成转换：

# 初始化结果数据框，先加入第一行原始数据
step_data <- cumhaz[1, ]

# 从第二行开始遍历原始数据
for (i in 2:nrow(cumhaz)) {
  # 添加「当前时间 + 上一行累积风险」的行
  step_data <- rbind(step_data, data.frame(
    Survival_Function.time = cumhaz$Survival_Function.time[i],
    Survival_Function.cumhaz = cumhaz$Survival_Function.cumhaz[i-1]
  ))
  # 添加当前行的原始数据
  step_data <- rbind(step_data, cumhaz[i, ])
}

# 查看转换后的数据
head(step_data, 11)

方法2：dplyr+tidyr实现（更简洁高效）

如果习惯用tidyverse系列工具，代码更简洁：

library(dplyr)
library(tidyr)

step_data <- cumhaz %>%
  # 生成前一个时间点的累积风险值
  mutate(prev_cumhaz = lag(Survival_Function.cumhaz)) %>%
  # 移除第一行（其前风险值为NA），拆分每行成两行
  slice(-1) %>%
  pivot_longer(
    cols = c(prev_cumhaz, Survival_Function.cumhaz),
    values_to = "Survival_Function.cumhaz"
  ) %>%
  # 保留需要的列，合并第一行原始数据并排序
  select(Survival_Function.time, Survival_Function.cumhaz) %>%
  bind_rows(cumhaz[1, ], .) %>%
  arrange(Survival_Function.time, Survival_Function.cumhaz)

# 查看结果
head(step_data, 11)

内容的提问来源于stack exchange，提问作者Asger Paludan-Müller