爾云間? 一個專門做科研的團隊

今天小云又遇到了一個大麻煩，小云很頭疼，因為小云要做的是獨立性分析，主要包括列線圖和校準曲線，于是小云就去找了個在線工具，結(jié)果畫出來的圖和作為標準的對角線相去甚遠，就像下面這樣。

小云一下子就慌了，嚶嚶嚶，問題到底出在哪兒呢？可能是在線工具默認了一些參數(shù)，無法調(diào)整，小云去查了一下，數(shù)據(jù)在做列線圖之前是要經(jīng)過單因素和多因素cox分析的，刪掉不顯著的特征才能做列線圖。而在線工具把多因素和列線圖連起來了，根本就不給你刪掉不顯著特征的機會，哼，那就只有自己用代碼跑了。

?

代碼如下：

library(survival)

library(survminer)

library(rms)

data=read.table("D:/卵巢癌免疫和內(nèi)質(zhì)網(wǎng)/獨立性分析/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為時間增量，這里可以選擇特征，因為之前小云已經(jīng)做過多因素分析了，階段也就是stage這個特征的多因素結(jié)果是不顯著的，因此小云在畫列線圖時就不帶上它了

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))#執(zhí)行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指數(shù)為1-C，也就是0.666，這種方式獲得的C指數(shù)沒有置信區(qū)間，可以使用boot函數(shù)計算

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)#圖形參數(shù)

plot(cal1,lwd=2,lty=1,subtitles = F,#關(guān)閉副標題

???? 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)#圖形參數(shù)

plot(cal3,lwd=2,lty=1,subtitles = F,#關(guān)閉副標題

???? 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)#圖形參數(shù)

plot(cal5,lwd=2,lty=1,subtitles = F,#關(guān)閉副標題

???? 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")

看，這樣的結(jié)果是不是就好多了，和對角線基本一致了。小云真是太棒了！

?

小伙伴們有沒有看明白呢，如果有不明白的歡迎來和小云討論喲！