又准又好，R语言绘制列线图和校准图
尔云间 一个专门做科研的团队
今天小云又遇到了一个大麻烦，小云很头疼，因为小云要做的是独立性分析，主要包括列线图和校准曲线，于是小云就去找了个在线工具，结果画出来的图和作为标准的对角线相去甚远，就像下面这样。
小云一下子就慌了，嘤嘤嘤，问题到底出在哪儿呢？可能是在线工具默认了一些参数，无法调整，小云去查了一下，数据在做列线图之前是要经过单因素和多因素cox分析的，删掉不显著的特征才能做列线图。而在线工具把多因素和列线图连起来了，根本就不给你删掉不显著特征的机会，哼，那就只有自己用代码跑了。
代码如下：
library(survival)
library(survminer)
library(rms)
data=read.table("D:/卵巢癌免疫和内质网/独立性分析/liexian.txt",header = T,sep="\t")
dd <- datadist(data)
options(datadist="dd")
f <- cph(Surv(Time, Status) ~ chemoresponse+Risk, x=T, y=T, surv=T, data=data,time.inc=1)#time.inc为时间增量，这里可以选择特征，因为之前小云已经做过多因素分析了，阶段也就是stage这个特征的多因素结果是不显著的，因此小云在画列线图时就不带上它了
surv <- Survival(f)
nom <- nomogram(f, fun=list(function(x) surv(1, x), function(x) surv(3, x), function(x) surv(5, x)),
lp=F, funlabel=c("1-year survival", "3-year survival", "5-year survival"), maxscale=100,
fun.at=c(0.95, 0.9, 0.8, 0.7, 0.6, 0.5))#执行Nomogram分析
plot(nom,xfrac=.5,cex.axis=1,cex.var=1)
validate(f, method="boot", B=1000, dxy=T)
rcorrcens(Surv(Time, Status) ~ predict(f), data = data)
#C指数为1-C，也就是0.666，这种方式获得的C指数没有置信区间，可以使用boot函数计算
library(boot)
c_index <- function(data,indices){
dat <- data[indices,]
vames<- c("chemoresponse","Risk")
FML <- as.formula(paste('Surv(Time, Status)~',paste(vames, collapse = "+")))
fit<- coxph(FML,data =dat )
pr1<-predict(fit,newdata=dat)
Cindex=rcorrcens(Surv(time, status) ~ pr1, data =dat)[1]
Cindex=1-Cindex
Cindex
}
c_index(data,1:100)
# [1] 0.666032
results <- boot(data=data, statistic=c_index, R=500)
boot.ci(results,conf = 0.95)
#下面绘制校准曲线
par(mfrow = c(1,3))
f1 <- cph(Surv(Time, Status) ~ chemoresponse+Risk, x=T, y=T, surv=T, data=data, time.inc=1)
cal1 <- calibrate(f1, cmethod="KM", method="boot", u=1,m=43,B=1000)
par(mar=c(6,5,1,2),cex = 1.0)#图形参数
plot(cal1,lwd=2,lty=1,subtitles = F,#关闭副标题
errbar.col=c(rgb(0,0,0,maxColorValue=255)),#误差线的颜色
xlim=c(0,1),ylim=c(0,1),#坐标范围
xlab="Predicted 1-year Overall Survival",ylab="Actual 1-year Overall Survival",
col=c(rgb(255,0,0,maxColorValue=255))#校正曲线的颜色
)
abline(0,1,lty=3,lwd=2,col="blue")#添加参考线
f3 <- cph(Surv(Time, Status) ~ chemoresponse+Risk, x=T, y=T, surv=T, data=data, time.inc=3)
cal3 <- calibrate(f3, cmethod="KM", method="boot", u=3,m=43,B=1000)
par(mar=c(6,5,1,2),cex = 1.0)#图形参数
plot(cal3,lwd=2,lty=1,subtitles = F,#关闭副标题
errbar.col=c(rgb(0,0,0,maxColorValue=255)),
xlim=c(0,1),ylim=c(0,1),
xlab="Predicted 3-year Overall Survival",ylab="Actual 3-year Overall Survival",
col=c(rgb(255,0,0,maxColorValue=255))
)
abline(0,1,lty=3,lwd=2,col="blue")
f5 <- cph(Surv(Time, Status) ~ chemoresponse+Risk, x=T, y=T, surv=T, data=data, time.inc=5)
cal5 <- calibrate(f5, cmethod="KM", method="boot", u=5,m=43,B=1000)
par(mar=c(6,5,1,2),cex = 1.0)#图形参数
plot(cal5,lwd=2,lty=1,subtitles = F,#关闭副标题
errbar.col=c(rgb(0,0,0,maxColorValue=255)),
xlim=c(0,1),ylim=c(0,1),
xlab="Predicted 5-year Overall Survival",ylab="Actual 5-year Overall Survival",
col=c(rgb(255,0,0,maxColorValue=255))
)
abline(0,1,lty=3,lwd=2,col="blue")
看，这样的结果是不是就好多了，和对角线基本一致了。小云真是太棒了！
小伙伴们有没有看明白呢，如果有不明白的欢迎来和小云讨论哟！
