如何计算ggadjustedcurves
survminer 函数的变异指数(SE 或 CI)?我正在使用条件方法。任何人有任何输入或资源?
2 回答
简短的回答:你不能。在文档中ggadjustedcurves
未指定置信区间选项。现在,如果您尝试conf.int
在调用中使用参数,ggadjustedcurves
它不会出错,但它也不会工作。同样在包的github
页面中报告了该问题,并且已请求但尚未添加此功能。
您可以ggsurvplot
通过添加从同一包轻松计算和绘制 CI函数conf.int = TRUE
。它ggsurvplot
是这样的:
library(survminer)
fit <- survfit(Surv(stop, event) ~ rx, data = bladder )
ggsurvplot(fit, data = bladder, conf.int=TRUE)
截至 2019 年底,仍然无法使用该ggadjustedcurves
函数绘制置信区间。survest
然而,可以使用包中的函数从cox 模型及其相应的置信区间估计调整后的生存概率rms
。
使用它,ggplot2
我编写了一个函数reshape2
,pammtools
该函数能够绘制带有置信区间的调整后的生存曲线估计值,由任意数量的变量分层:
library(rms)
library(ggplot2)
library(reshape2)
library(pammtools)
adjusted_surv_curves <- function(model_cph, ..., ci=TRUE, conf.int=0.95,
from=0, to=max(surv_prob$time), step=0.5,
xlab="Time", ylab="Survival Probability", title="",
labels=names(strata), legend.title="Strata",
plot.ylim=c(0, 1), theme=theme_bw(), size=1,
palette="Set1", ci.alpha=0.2, return.data=F) {
# check input
stopifnot(class(model_cph)==c("cph", "rms", "coxph"))
# check if packages are loaded
if(!(all(c("ggplot2", "reshape2", "rms", "pammtools") %in% (.packages())))){
stop("One or more of the following packages are not attached: ggplot2, reshape2, rms, pammtools")
}
# to get strata, this needs to be called without
# specified times
surv_prob <- survest(model_cph)
strata <- surv_prob$strata
# calc adjusted survival probabilities
surv_prob <- survest(model_cph, times = seq(from, to, by=step),
conf.int = conf.int)
# extract estimates from survest object
plotd_surv <- data.frame(time=surv_prob$time)
plotd_lower <- data.frame(time=surv_prob$time)
plotd_upper <- data.frame(time=surv_prob$time)
for (i in 1:length(strata)) {
plotd_surv[,ncol(plotd_surv)+1] <- surv_prob$surv[i,]
plotd_lower[,ncol(plotd_lower)+1] <- surv_prob$lower[i,]
plotd_upper[,ncol(plotd_upper)+1] <- surv_prob$upper[i,]
}
# put together in new data frame
# melt dataframes
plotd_surv <- melt(plotd_surv, id.vars = "time")
colnames(plotd_surv)[3] <- "est"
plotd_lower <- melt(plotd_lower, id.vars = "time")
colnames(plotd_lower)[3] <- "lower"
plotd_upper <- melt(plotd_upper, id.vars = "time")
colnames(plotd_upper)[3] <- "upper"
# merge to one
plotdata <- merge(plotd_surv, plotd_lower, by=c("time", "variable"))
plotdata <- merge(plotdata, plotd_upper, by=c("time", "variable"))
# return data instead, if specified
if (return.data) {
return(plotdata)
}
# plot curves
p <- ggplot(data=plotdata, aes(x=time)) +
geom_step(aes(y=est, color=variable), size=size) +
theme + ylim(plot.ylim) +
scale_colour_brewer(palette = palette, name=legend.title,
labels=labels) +
xlab(xlab) + ylab(ylab) + ggtitle(title)
# add confidence interval if not specified otherwise
if (ci) {
p <- p + pammtools::geom_stepribbon(
aes(ymin=lower, ymax=upper, fill=variable),
alpha=ci.alpha) +
scale_fill_brewer(palette = palette, name=legend.title,
labels=labels)
}
# add additional ggplot parameter, if specified
add_params <- list(...)
if (length(add_params)!=0) {
for (i in 1:length(add_params)) {
p <- p + add_params[[i]]
}
}
# return plot
return(p)
}
作为输入,函数定义的模型cph
是必要的。cph
建议的函数将使用函数在调用中定义的所有层变量的所有值分割曲线strat()
(您还应该指定x=TRUE
和y=TRUE
以进行正确估计)。,to
和参数可用于指定应以什么间隔from
以及step
多久计算一次生存概率。绘图本身的各种参数可以随意修改,因为它还接受ggplot
调用的任何进一步参数。
注意:生存概率的估计及其相应的置信区间将与ggadjustedcurves
函数所做的估计不同,因为它们的计算使用了不同的方法。有关如何计算置信区间的更多信息,您可以查阅官方文档:https ://rdrr.io/cran/rms/man/survest.cph.html
例如,假设函数已经定义:
library(rms)
library(ggplot2)
library(reshape2)
library(pammtools)
# load example dataset
data("lung")
# define model
model <- cph(with(data=lung, Surv(time, status) ~ strat(sex) + age + ph.ecog), x=TRUE, y=TRUE)
# plot adjusted survival curves by sex
adjusted_surv_curves(model)
上面的代码给出了以下输出:
我知道代码很丑陋,而且不是最有效的,但它对我来说是应该的,我希望它对你或任何碰巧遇到这个问题的人有帮助。