导语
GUIDE ╲
生存分析是指将终点事件和出现此事件所经历的时间结合起来分析的一种统计方法,研究生存现象和现象的响应时间数据及其规律,在肿瘤等疾病研究中运用广泛。在R中进行生存分析常用的包有survival包以及survminer包。[A Robust 6-lncRNA Prognostic Signature for Predicting the Prognosis of Patients With Colorectal Cancer Metastasis]中R包survminer用于确定高风险和低风险人群的最佳截点。[Change in Neutrophil to Lymphocyte Ratio During Immunotherapy Treatment Is a Non-Linear Predictor of Patient Outcomes in Advanced Cancers]采用survminer包进行单因素和多因素分析。今天我们来介绍survminer包。
survminer包应用
01
ggsurvplot()
install.packages(“survminer”);library(survminer);require("survival")
fit <- survfit(Surv(time, status) ~ sex, data = lung)#也可以用surv_fit
fit #Fig 1
# Basic survival curves
ggsurvplot(fit, data = lung) #Fig 2
ggsurvplot(fit, data = lung,
surv.median.line = "hv", #用于在中位生存值绘制水平/垂直线,c(“none”,“hv”,“h”,“v”)。
legend.title = "Sex",legend.labs = c("Male", "Female"), #改变图例
pval = "The hot p-value is: 0.003",conf.int = TRUE, #加p-value和置信区间
risk.table = TRUE, #加入风险表
tables.height = 0.2, #生存图下的表的高度,[0~1]。
palette = c("#E7B800", "#2E9FDF"),#曲线颜色
ggtheme = theme_bw(), #改变ggplot2背景theme_linedraw,theme_light,theme_gray等。
font.main = c(16, "bold", "darkblue"),font.x = c(14, "bold.italic", "red"),font.y = c(14, "bold.italic", "darkred"),font.tickslab = c(12, "plain", "darkgreen")# 改变字体大小类型颜色等
) #Fig 3
Fig 1
Fig 2
Fig 3
fit <- survfit( Surv(time, status) ~ sex + rx + adhere,data = colon )
ggsurvplot(fit, data = colon,
fun = "cumhaz", #定义生存曲线变换的任意函数,"event"累积事件f(y) = 1-y, "cumhaz"累积风险函数f(y) = -log(y),"pct"生存概率(百分比)。
conf.int = TRUE, risk.table = TRUE, risk.table.col="strata",
ggtheme = theme_bw())
#fun = "event"--Fig 4.1
#fun = "cumhaz"--Fig 4.2
#fun = "pct"--Fig 4.3
Fig 4.1
Fig 4.2
Fig 4.3
02
arrange_ggsurvplots():在同一页面展示多个曲线图
splots <- list()#创建一个列表并设置成不同的背景色
splots[[1]] <- ggsurvplot(fit, data = lung, risk.table = TRUE, ggtheme = theme_bw())
splots[[2]] <- ggsurvplot(fit, data = lung, risk.table = TRUE, ggtheme = theme_grey())
arrange_ggsurvplots(splots,
print = TRUE,#为TRUE时展示排列好的图形
ncol = 2, nrow = 1,#行和列的数量
risk.table.height = 0.4#这里设置risk.table为TRUE,风险表的高度默认0.25,有多个表时适当增加高度。
) #Fig 5
Fig 5
03
surv_median() & surv_pvalue() & surv_summary()
fit1 <- surv_fit(Surv(time, status) ~ sex, data = colon)
fit2 <- surv_fit(Surv(time, status) ~ adhere, data = colon)
fit.list <- list(sex = fit1, adhere = fit2)
surv_median(fit1)
surv_median(fit.list)
surv_median(fit.list, combine = TRUE) #Fig 6
surv_pvalue(fit1, data = colon, method = "survdiff",
get_coord = TRUE) #Fig 7
#pval:p值
#method: 计算p值的方法
#pval.txt: 用于注释图形的文本
#pval.x, pval.y: 注释图形的x&y的坐标
#method.x, method.y: 计算p值的方法的坐标x & y coordinates of pvalue method
res.sum <- surv_summary(fit1, data = colon);head(res.sum) #Fig 8
#time: 在曲线上有一阶跃的时间点。
#n.risk: 当前时间的受试者人数。
#n.event: 当前时间的事件数目。
#n.censor: 删失事件数目。
#surv: 生存评估。
#std.err: 标准误差。
#upper/lower: 置信区间的上下限。
#strata: 生存曲线的类别
Fig 6
Fig 7
Fig 8
04
ggsurvplot_df():从任何包含生存数据的数据框绘制生存曲线
head(surv_summary(fit2, colon)) #Fig 9
ggsurvplot_df(surv_summary(fit2, colon), conf.int = TRUE,
legend.title = "Adhere", legend.labs = c("0", "1"))#Fig 10
Fig 9
Fig 10
05
ggsurvevents():事件的时间分布
surv <- Surv(lung$time, lung$status)
ggsurvevents(surv,
normalized = TRUE,#展示相对事件数目
type = "radius") #Fig 11.1
#"cumulative"累积事件数,"radius"给定半径内的事件数
ggsurvevents(surv, normalized = TRUE, type = "cumulative")#Fig 11.2
ggsurvevents(surv, normalized = TRUE, type = "fraction")#Fig 11.3
Fig 11.1
Fig 11.2
Fig 11.3
06
ggcoxfunctional():展示解释连续变量的图,有助于在cox模型中正确选择函数变量形式。
data(mgus)
res.cox <- coxph(Surv(futime, death) ~ mspike + log(mspike) + I(mspike^2) +age + I(log(age)^2) + I(sqrt(age)), data = mgus)
res.cox #Fig 12
ggcoxfunctional(res.cox, data = mgus,
point.col = "blue", point.size = 0.5, point.shape = 8, point.alpha = 0.5,#点的颜色大小形状和可见度
title = "Title", caption = "Caption"#主标题,纵标题
)#Fig 13
Fig 12
Fig 13
07
ggcoxzph():显示按比例缩放的Schoenfeld残差的图形及曲线
fit <- coxph(Surv(futime, fustat) ~ age + ecog.ps + rx, data=ovarian)
cox.zph.fit <- cox.zph(fit)
ggcoxzph(cox.zph.fit, var = c("ecog.ps", "rx", "age"),
font.main =c(12,"darkred"), caption = "Caption",
submain="Submain",font.y=c(14, "bold.italic", "darkred"),
font.x = c(14, "bold", "red")) #Fig 14
Fig 14
cox回归的等比例风险的验证有多种方法,可以利用Schoenfeld残差图检验。如果Schoenfeld残差与时间t无明显的变化趋势,即Schoenfeld残差与时间t无关,则提示符合等比例风险假设。Fig 14中,三个P值都大于0.05,说明每个变量均满足PH检验。总体P=0.2691,模型整体满足PH检验。图中虚线表示拟合曲线上下2个单位的标准差。
08
ggcompetingrisks():绘制累计关联曲线
require("cmprsk")
#生成随机数
set.seed(2);ss <- rexp(100)
gg <- factor(sample(1:3,100,replace=TRUE),1:3,c('BRCA','LUNG','OV'))
cc <- factor(sample(0:2,100,replace=TRUE),0:2,c('no event', 'death', 'progression'))
strt <- sample(1:2,100,replace=TRUE)
# fit:cmprsk::cuminc类的对象——用cmprsk::cuminc函数创建,或者用survfit函数创建
print(fit <- cmprsk::cuminc(ss,cc,gg,strt)) #Fig 15
ggcompetingrisks(fit)#Fig 16
ggcompetingrisks(fit, multiple_panels = FALSE, conf.int = TRUE)#不分组 Fig 17
fit <- survfit(Surv(time, status, type="mstate") ~ group, data=df)
ggcompetingrisks(fit)#多组 Fig 18
Fig 15
Fig 16
Fig 17
Fig 18
09
ggforest():绘制Cox比例风险模型森林图
model <- coxph( Surv(time, status) ~ sex + rx + adhere, data = colon )
ggforest(model) #Fig 19
colon <- within(colon, {sex <- factor(sex, labels = c("female", "male"))
differ <- factor(differ, labels = c("well", "moderate", "poor"))
extent <- factor(extent, labels = c("submuc.", "muscle", "serosa", "contig."))
})
head(colon)#Fig 20
bigmodel<-coxph(Surv(time, status) ~ sex + rx + adhere + differ + extent + node4,data = colon )
ggforest(bigmodel,
fontsize=1.5,#图中注释的相对大小,默认0.7。
noDigits=4#图中估计值和p值的显示位数
)#Fig 21
Fig 19
Fig 20
Fig 21
10
surv_adjustedcurves() & ggadjustedcurves():计算并绘制coxph模型调整后的生存曲线
fit <- coxph( Surv(stop, event) ~ size + strata(rx), data = bladder )
curve <- surv_adjustedcurves(fit, data = bladder, method = "average", variable = "rx")
ggadjustedcurves(fit, data = bladder,
method = "average",#描述预期生存曲线应如何计算。single(总体平均),average(亚组平均),’marginal’, ’conditional’(再平衡后的亚组平均)。
variable = "rx",
palette=c("red","blue")#曲线颜色
size=0.5#曲线宽度
)# Fig 22
Fig 22
小编总结
今天我们分享了survminer包,生存分析不再是简单的两条“蓝红”曲线了,你还不找出陈年老图给它翻新一下吗?