time-dependent平滑ROC曲线

前几天把能够实现ROC曲线最佳截点的R包简单汇总了一下：ROC曲线最佳截点

公众号后台回复最佳截点即可获取最佳截点和推文合集；回复ROC即可获取ROC曲线合集推文。

在查看相关资料时发现了risksetROC这个R包，也可以实现time-dependent ROC曲线的绘制，而且可以实现平滑曲线。

平滑曲线这个功能在timeROC和survivalROC中都是不能实现的，我们介绍过一种借助ggplot2实现的方法：生存资料ROC曲线的最佳截点和平滑曲线

下面简单介绍下如何实现平滑的time-dependent ROC曲线。分类资料的平滑ROC曲线可以通过pROC实现，可参考：分类资料的ROC曲线绘制

安装

install.packages("risksetROC")

准备数据

使用survival包中的pbc数据集。

library(risksetROC)

library(survival)
data(pbc)
str(pbc)

## 'data.frame':	418 obs. of  20 variables:
##  $ id      : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ fudays  : int  400 4500 1012 1925 1504 2503 1832 2466 2400 51 ...
##  $ status  : int  2 0 2 2 1 2 0 2 2 2 ...
##  $ drug    : int  1 1 1 1 2 2 2 2 1 2 ...
##  $ age     : int  21464 20617 25594 19994 13918 24201 20284 19379 15526 25772 ...
##  $ sex     : int  1 1 0 1 1 1 1 1 1 1 ...
##  $ ascites : int  1 0 0 0 0 0 0 0 0 1 ...
##  $ hepatom : int  1 1 0 1 1 1 1 0 0 0 ...
##  $ spiders : int  1 1 0 1 1 0 0 0 1 1 ...
##  $ edema   : num  1 0 0.5 0.5 0 0 0 0 0 1 ...
##  $ bili    : num  14.5 1.1 1.4 1.8 3.4 0.8 1 0.3 3.2 12.6 ...
##  $ chol    : int  261 302 176 244 279 248 322 280 562 200 ...
##  $ albumin : num  2.6 4.14 3.48 2.54 3.53 3.98 4.09 4 3.08 2.74 ...
##  $ copper  : int  156 54 210 64 143 50 52 52 79 140 ...
##  $ alkphos : num  1718 7395 516 6122 671 ...
##  $ sgot    : num  137.9 113.5 96.1 60.6 113.2 ...
##  $ trig    : int  172 88 55 92 72 63 213 189 88 143 ...
##  $ platelet: int  190 221 151 183 136 NA 204 373 251 302 ...
##  $ protime : num  12.2 10.6 12 10.3 10.9 11 9.7 11 11 11.5 ...
##  $ stage   : int  4 3 4 4 3 3 3 3 2 4 ...

做一些准备工作。

只使用前312行数据，生存状态用1表示终点事件，0表示删失，然后建立cox模型，计算出线性预测值作为marker：

pbc1 <- pbc[1:312,]
survival.status <- ifelse(pbc1$status==2,1,0)
survival.time <- pbc1$fudays

pbc1$status1 <- survival.status

fit <- coxph(Surv(survival.time, status1) ~ log(bili)+
               log(protime)+edema+albumin+age,
             data = pbc1
             )
eta <- fit$linear.predictors

使用

使用方法和timeROC以及survivalROC基本上是一样的：

nobs <- length(survival.time[survival.status==1])
span <- 1.0*(nobs^(-0.2))

# 3种方法都试一下，然后画在一起
ROC.CC90 <- risksetROC(Stime = survival.time, status = survival.status,
                       marker = eta, predict.time = 90, method = "Cox",
                       main="time-denpendent ROC with riksetROC",
                       lty=2, lwd=2,col="red"
                       )
ROC.SS90 <- risksetROC(Stime = survival.time, status = survival.status,
                       marker = eta, predict.time = 90, method = "Schoenfeld",
                       plot = F, span = span
                       )
ROC.LL90 <- risksetROC(Stime = survival.time, status = survival.status,
                       marker = eta, predict.time = 90, method = "LocalCox",
                       plot = F, span = span
                       )

lines(ROC.SS90$FP, ROC.SS90$TP, lty=3, lwd=2, col="darkblue")
lines(ROC.LL90$FP, ROC.LL90$TP, lty=4, lwd=2, col="green")
legend(0.6,0.25, lty = c(2,3,4),col = c("red","darkblue","green"),
       legend = c("Cox","Schoenfeld","LocalCox"), bty = "n")

完美的平滑曲线，不用自己实现。

对比timeROC

下面对比下timeROC的实现：

library(timeROC)

ROC <- timeROC(T=survival.time,   
               delta=survival.status,   
               marker=eta,   
               cause=1,                
               weighting="cox",   
               times=90)

plot(ROC,time = 90)

不能平滑。

参考资料

https://faculty.washington.edu/heagerty/Software/SurvROC/RisksetROC/risksetROCdiscuss.pdf