又准又好，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-CindexCindex}c_index(data,1:100)# [1] 0.666032results <- 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")
看，这样的结果是不是就好多了，和对角线基本一致了。小果真是太棒了！
小伙伴们有没有看明白呢，如果有不明白的欢迎来和小果讨论哟！
扫码关注我们
shengxinguoer
生信果
生信硬核知识解答
和小果一起学生信