使用场景

我们在做数据分析的时候经常需要用到ROC曲线，我们一般都在下面的两种情况用到ROC曲线：

1. 我们在对一个连续性的指标进行cutoff分析的时候，如果结局变量是二分类的话，可以通过ROC曲线来计算最佳的cut-off值。
2. 我们在做构建了模型的时候，可以通过ROC曲线来评价构建的模型的诊断效能。如果我们要评价多个模型的同样的也是可以通过ROC曲线来比较

ROC曲线有一些特定的值都是需要在做ROC的时候给出的。例如：灵敏度；特意度；阈值；曲线下面积以及约登指数等等。

R语言实现

R语言可以实现绘制ROC曲线的有很多，例如：pROC;RROC;cutpointr都可以对一个变量进行ROC计算以及相关数据的统计。这次我们就通过pROC包的相关函数来批量计算ROC曲线的结果。

1. pROC包当中提供了一个函数叫roc，我们可以通过这个函数来计算两组之间的ROC结果
2. 提供了ci.auc函数了计算ROC曲线的曲线下面积以及95%CI
3. 提供了coords来计算ROC曲线的相关信息。这个函数可以接受的参数包括：“threshold”, “specificity”, “sensitivity”, “accuracy”, “tn” (true negative count), “tp” (true positive count), “fn” (false negative count), “fp” (false positive count), “npv” (negative predictive value), “ppv” (positive predictive value), “precision”, “recall”. “1-specificity”, “1-sensitivity”, “1-accuracy”, “1-npv” and “1-ppv”`。我们可以选择不同的参数来返回不同的结果
4. 这个函数包好像没有通过约登指数的计算，不过越登指数就是灵敏度+特意度-1。我们可以自己计算。

通过上面的需求，我们可以创建一个函数来返回相关的ROC分析的结果

### 加载使用的包
library(tidyverse)
library(pROC)

ROCStatFunc <- function(dat, group, var,retype = c("threshold", "specificity", "sensitivity"),
                        auc = T,youden = T, digit = 3){
     subgroup <- levels(as.factor(dat[[group]]))
     subgroup1 <- paste0(subgroup[2], " vs ", subgroup[1])
     rocmodel <- roc(dat[[group]], dat[[var]])
     other <- coords(rocmodel, "b", ret = retype)
     other <- round(other, digit)
     if(auc == T){
         auc <- round(ci.auc(rocmodel),digit)
         auc <- paste0(auc[2],"(",auc[1],"-",auc[3],")")
         if(youden == T){
             abc <- coords(rocmodel, "b", ret = c("specificity", "sensitivity"))
             youdenres <- abc[1] + abc[2] - 1
             youdenres <- round(youdenres, digit)
             result <- c(group, subgroup1, auc, other, youdenres)
             names(result) <- c("group", "subgroup","auc(95%CI)", retype, "youden")
         }else{
             result <- c(group, subgroup1, auc, other)
             names(result) <- c("group", "subgroup", "auc(95%CI)", retype)
         }
     }else{
         if(youden == T){
             abc <- coords(rocmodel, "b", ret = c("specificity", "sensitivity"))
             youdenres <- abc[1] + abc[2] - 1
             youdenres <- round(youdenres, digit)
             result <- c(group, subgroup1, other, youdenres)
             names(result) <- c("group","subgroup", retype, "youden")
         }else{
             result <- c(group, subgroup1,other)
             names(result) <- c("group", "subgroup",retype)
         }
     }
     return(result)
}

函数解读

我们构建的函数一共包括7个参数。分别是：

• dat: 我们需要输入我们的数据集
• group: 想要进行ROC的结局变量
• var： 需要分析的变量
• retype: 想要计算的具体信息: 默认的返回的包括：阈值，灵敏度， 特异度。
• auc： 是否返回曲线下面积以及95%CI，逻辑值，默认为TRUE。
• youden: 是否返回约登指数，逻辑值，默认为TURE。
• digit: 结果保留几位小数点

由于不是自己写的函数，我们在使用pROC包的时候会返回很多信息。我们可以去掉这些结果

quiteROCFunc <- quietly(ROCStatFunc)

我们使用示例数据来查看结果

data("aSAH")
head(aSAH)

### 计算outcome为结局变量的age的相关信息
quiteROCFunc(aSAH, group = "s100b", var = "age")$result

# 批量计算变量的ROC结果
## 定义group
multigroup <- c("age", "s100b", "ndka")
rocRes <- lapply(multigroup, function(x) quiteROCFunc(aSAH, "outcome", x)$result)
rocResDat <- do.call(rbind, rocRes)
rocResDat

PS:

1. 这个函数依赖的是tidyververse以及pROC包，运行之前需要提前加载
2. 函数对于分组的话必须是两组。
3. 结果当中的subgroup是case vs control。