1. 使用 $\color{green}{survival}$ 包绘制生存曲线图

1.1 Kaplan-Meier分析和Log-rank test

# 安装并加载所需的R包
# install.packages("survival")
library(survival)

# 查看内置数据集
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

# Kaplan-Meier 分析，
KM <- survfit(Surv(time, status) ~ sex, data = lung, type = 'kaplan-meier', conf.type = 'log')
KM

summary(KM) # 展示每个时间点的细节
## Call: survfit(formula = Surv(time, status) ~ sex, data = lung, type = "kaplan-meier", 
##     conf.type = "log")
## 
##                 sex=1 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    11    138       3   0.9783  0.0124       0.9542        1.000
##    12    135       1   0.9710  0.0143       0.9434        0.999
##    13    134       2   0.9565  0.0174       0.9231        0.991
##    15    132       1   0.9493  0.0187       0.9134        0.987

summary(KM)$table # 结果提取，包含各组生存状态、中位生存时间、置信区间等统计
##       records n.max n.start events    rmean se(rmean) median 0.95LCL 0.95UCL
##  sex=1     138   138     138    112 326.0841  22.91156    270     212     310
##  sex=2      90    90      90     53 460.6473  34.68985    426     348     550

# 对数秩检验（Log-rank test）比较组间生存曲线，可以获得组间生存率的差异p值。
# Log-rank test是非参数检验，近似于卡方检验，零假设是组间没有差异
survdiff(Surv(time, status) ~ sex, data = lung)
## 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.001

★ p=0.001，组间差异非常显著，即两组患者之间的生存率具有显著差异。

1.2 绘制生存曲线

# 绘制生存曲线，反映了尚在世的患者数量比例和时间的关系
plot(KM, main = 'Kaplan-Meier ', xlab = 'Time (days)', ylab = 'Overall survival', 
    lwd = 2, col = c('blue', 'red'))
 
legend(x = 'topright', col = c('blue', 'red'), lwd = 2, legend = c('1: male', '2: female'))
 
# 绘制累积风险曲线，反映了疾病风险和时间的关系，与累积的去世患者数量有关
plot(KM, main = 'Cumulative hazard', xlab = 'Time (days)', ylab = 'Cumulative hazard', 
    lwd = 2, col = c('blue', 'red'), fun = 'cumhaz')
 
legend(x = 'topright', col = c('blue', 'red'), lwd = 2, legend = c('1: male', '2: female'))

survival-1

2. 使用 $\color{green}{survminer}$ 包绘制生存曲线图

2.1 基本用法

# 安装并加载所需的R包
# install.packages("survminer")
# if(!require(devtools)) install.packages("devtools")
# devtools::install_github("kassambara/survminer")
library(survminer)

# 查看内置数据集
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

# Fit survival curves
require("survival")

# 使用survfit()函数拟合KM生存曲线
fit <- survfit(Surv(time, status) ~ sex, data = lung)

# 使用ggsurvplot()函数绘制基础KM生存曲线
ggsurvplot(fit)

survminer-1

2.2 其他参数

ggsurvplot(fit, data = lung,  
           main = "Survival curve", # 添加标题
           surv.median.line = "hv", # 添加中位数生存时间线
           font.main = c(16, "bold", "darkblue"), # 设置标题字体大小、格式和颜色
           font.x = c(14, "bold.italic", "red"), # 设置x轴字体大小、格式和颜色
           font.y = c(14, "bold.italic", "darkred"), # 设置y轴字体大小、格式和颜色
           font.tickslab = c(12, "plain", "darkgreen"), # 设置坐标轴刻度字体大小、格式和颜色
           legend = c(0.9, 0.4), # 通过坐标指定图例位置
           legend.title = "Sex", legend.labs = c("Male", "Female"), # 更改图例标题和标签
           size = 1,  # 改变线条大小
           linetype = "strata", # 按组更改线型
           break.time.by = 250, # 将时间轴以250作为打断
           # palette = "Dark2" # 使用 brewer 调色板“Dark2”
           palette = c("#E7B800", "#2E9FDF"), # 自定义调色板
           ggtheme = theme_bw(), # 修改主题
           conf.int = TRUE, # 添加置信区间
           pval = TRUE, # 添加p值
           risk.table = TRUE, risk.table.y.text.col = TRUE, # 添加风险表,并按层更改风险表y 文本颜色  
           tables.height = 0.2, # 设置风险表的高度
           tables.theme = theme_cleantable(), # 设置风险表的主题
           xlim = c(0, 1050) # 更改 x 轴范围
           )

survminer-2.png

2.3 多组生存曲线

head(colon)
##   id study      rx sex age obstruct perfor adhere nodes status differ extent surg node4 time etype
## 1  1     1 Lev+5FU   1  43        0      0      0     5      1      2      3    0     1 1521     2
## 2  1     1 Lev+5FU   1  43        0      0      0     5      1      2      3    0     1  968     1
## 3  2     1 Lev+5FU   1  63        0      0      0     1      0      2      3    0     0 3087     2
## 4  2     1 Lev+5FU   1  63        0      0      0     1      0      2      3    0     0 3087     1
## 5  3     1     Obs   0  71        0      0      1     7      1      2      2    0     1  963     2
## 6  3     1     Obs   0  71        0      0      1     7      1      2      2    0     1  542     1

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

ggsurvplot(fit2, pval = TRUE, 
           break.time.by = 800,
           risk.table = TRUE,
           risk.table.col = "strata",
           risk.table.height = 0.5, 
           palette = "Dark2", 
           ggtheme = theme_bw(),
           legend.labs = c("A", "B", "C", "D", "E", "F") # 更改图例标签
           )

# 绘制分面生存曲线
ggsurvplot_facet(fit, colon, 
                 facet.by = c("rx", "adhere"),
                 palette = "npg", 
                 pval = TRUE,
                 surv.median.line = "hv",  # 增加中位生存时间
                 conf.int = TRUE) # 增加置信区间