一步到位：手把手教你R语言竞争风险模型建模-列线图-校准曲线-K折验证-外部验证- 决策曲线

library(QHScrnomo) bc<-read.csv("E:/r/test/qianliexian.csv",sep=',',header=TRUE) ########3 set.seed(123) tr1<- sample(nrow(bc),0.7*nrow(bc))##随机无放抽取 bc_train <- bc[tr1,]#70%数据集 bc_test<- bc[-tr1,]#30%数据集 ########## dd <- datadist(bc_train) options(datadist = "dd") prostate.f <- cph(Surv(TIME_EVENT,EVENT_DOD == 1) ~ TX + rcs(PSA,3) + BX_GLSN_CAT + CLIN_STG + rcs(AGE,3) + RACE_AA, data = bc_train, x = TRUE, y= TRUE, surv=TRUE,time.inc = 144) prostate.crr <- crr.fit(prostate.f,cencode = 0,failcode = 1) summary(prostate.crr) prostate.g <- Newlabels(prostate.crr, c(TX = 'Treatment options', BX_GLSN_CAT = 'Biopsy Gleason Score Sum', CLIN_STG = 'Clinical stage')) nomogram.crr(prostate.g, failtime = 120, lp=FALSE, xfrac=0.65, fun.at = seq(0.2, 0.45, 0.05), funlabel = "Predicted 10-year cumulative incidence") sas.cmprsk(prostate.crr, time = 120) ############## bc_train$preds.tenf <- tenf.crr(prostate.crr, time=120, fold = 10) with(bc_train, cindex(preds.tenf, ftime = TIME_EVENT, fstatus =EVENT_DOD, type = "crr"))["cindex"] with(bc_train, groupci( preds.tenf, ftime = TIME_EVENT, fstatus =EVENT_DOD, g = 5, u = 120, xlab = "Nomogram predicted 10-year cancerspecific mortality", ylab = "Observed predicted 10-year cancerspecific mortality") ) ##################################外部验证，主要是要生成验证集的预测概率 ####建模数据集建立模型 time<- bc_train [, "TIME_EVENT"] status<-bc_train [, "EVENT_DOD"] x<-bc_train [, seq(1, 7)]#x不能包含字符 x$TX<-as.numeric(as.factor(x$TX)) x$CLIN_STG<-as.numeric(as.factor(x$CLIN_STG)) newfit <- crr(time,status,x,cencode = 0,failcode = 1)###建模 ######验证集生成概率 time1<-bc_test[, "TIME_EVENT"] status1<-bc_test[, "EVENT_DOD"] x1<-bc_test[, seq(1, 7)]#x1不能包含字符 x1$TX<-as.numeric(as.factor(x1$TX))##转换位数字 x1$CLIN_STG<-as.numeric(as.factor(x1$CLIN_STG))##转换位数字 #把UNIQIDz从int转成num x3<-apply(x1,2,function(x1) as.numeric(x1)) ### bc_test$preds.tenf1<-QHScrnomo:::pred3.crr(newfit, as.matrix(x3), time = 120, lps = F)###生成验证集概率 ###验证集cindex with(bc_test, cindex(preds.tenf1, ftime = TIME_EVENT, fstatus =EVENT_DOD, type = "crr"))["cindex"] ###验证集校准标准曲线 with(bc_test, groupci( preds.tenf1, ftime = TIME_EVENT, fstatus =EVENT_DOD, g = 5, u = 120, xlab = "Nomogram predicted 10-year cancerspecific mortality", ylab = "Observed predicted 10-year cancerspecific mortality") ) ############决策曲线 library(dcurves) library(gtsummary); library(dplyr); library(tidyr) prostate.f1 <- coxph(Surv(TIME_EVENT,EVENT_DOD == 1) ~ TX + rcs(PSA,3) + BX_GLSN_CAT + CLIN_STG + rcs(AGE,3) + RACE_AA, data = bc_train) newdata1<- broom::augment(prostate.f1, newdata = bc_train %>% mutate(TIME_EVENT = 120), type.predict = "expected" ) %>% mutate(pr_failure120 = 1 - exp(-.fitted) ) %>% select(-.fitted, -.se.fit) dca(Surv(TIME_EVENT, EVENT_DOD) ~ pr_failure120, data = newdata1, time = 120, thresholds = 1:50 / 100) %>% plot(smooth = TRUE)