46w数据的auc,pROC的auc函数算的非常非常非常慢。
同样的数据量在python中scikit learn包大概只要10s。
求问是否有高效的包

    这个恐怕得深入研究源代码以及用代码分析工具(如 profvis)测速之后才能回答了。感觉这是一个潜在的铜矿可以挖挖,研究明白之后没准可以大幅提高计算速度。

    Hi 热心西红柿啰嗦两句,
    pROC慢的原因是,用auc()的时候,它调用了roc类及其下面的方法,所以你会看到算了sensitivity,specificity等一些指标。从profvis也可以看出其层层调用,这是这个包的面向对象设计的问题。有人已经提出这个issue了:https://github.com/xrobin/pROC/issues/17

    ROCR很快是因为performance()只调用了auc方法,所以只算了AUC这个指标,代码如下:

    .performance.auc <-
      function(predictions, labels, cutoffs, fp, tp, fn, tn,
               n.pos, n.neg, n.pos.pred, n.neg.pred, fpr.stop) {
          
          x <- fp / n.neg
          y <- tp / n.pos
    
          finite.bool <- is.finite(x) & is.finite(y)
          x <- x[ finite.bool ]
          y <- y[ finite.bool ]
          if (length(x) < 2) {
              stop(paste("Not enough distinct predictions to compute area",
                         "under the ROC curve."))
          }
    
          if (fpr.stop < 1) {
            ind <- max(which( x <= fpr.stop ))
            tpr.stop <- approxfun( x[ind:(ind+1)], y[ind:(ind+1)] )(fpr.stop)
            x <- c(x[1:ind], fpr.stop)
            y <- c(y[1:ind], tpr.stop)
          }
          
          ans <- list()
          auc <- 0
          for (i in 2:length(x)) {
              auc <- auc + 0.5 * (x[i] - x[i-1]) * (y[i] + y[i-1])
          }
          ans <- list( c(), auc)
          names(ans) <- c("x.values","y.values")
          return(ans)
      }

    以上。嘻嘻~

      tomatoiscoding 赞热心以及用科学的方法谈速度。

      所以我说这些包里面大有铜矿可挖嘛(算不上金矿,但多少有些价值),我就随便瞄了一眼,显然这个循环是可以很轻易向量化的:

      for (i in 2:length(x)) {
        auc <- auc + 0.5 * (x[i] - x[i-1]) * (y[i] + y[i-1])
      }

      循个啥子环嘛,不就是算一系列梯形的面积之和咩,用一点几何想象力,把向量分别掐头、去尾、横坐标作差、纵坐标求和,最后点乘、除以 2 就是曲线下的面积了噻。

      n = length(x)
      auc = (x[-1] - x[-n]) * (y[-1] + y[-n]) / 2

      我不知道横坐标有没有等间距的假设,要是有的话,连两两作差都不用了。假设越多,代码就有越快的可能性。

      期待某位少侠来把这事儿提速二百倍,没准儿将来可以作为在 R 江湖的第一件成名作哦。

        5 天 后
        9 天 后

        lijian 对 x 来说是 diff(),对 y 来说是相邻元素两两求和。我之所以没用 diff() 是因为我猜想:

        1. 对 y 两两求和的时候需要用到向量长度,所以我提前算了 n = length(x),然后这个 n 可以在 x 和 y 的计算中二度利用;如果用 diff(),其内部还是需要求 x 的长度。

        2. 下标运算和减法应该都足够快,所以用 diff()应该优势不大。

        刚简单测了一下速,对长度较短的向量,我的办法会比 diff() 明显快;随着长度加大,两个办法的速度会趋于相同。下面测速只考虑了我的第二点猜想,没考虑第一点可能的优势(但长度计算应该可以快到忽略不计)。

        N = 100
        x = runif(N)
        microbenchmark::microbenchmark(
          {n = length(x); x[-1] - x[-n]},
          diff(x),
          times = 1000
        )
        Unit: microseconds
                                            expr   min    lq     mean median     uq     max neval
         {     n = length(x)     x[-1] - x[-n] } 2.616 2.894 4.474535  3.165 4.9205 159.119  1000
                                         diff(x) 5.581 6.058 8.613288  6.706 8.7680 150.211  1000

        yihui 估计这个差不多就是极限了。分别在200个cutoff上的fpr和tpr,然后就是这200个小矩形的和。

        3 个月 后

        我写的一个用于评分卡建模的包scorecard提供了计算auc的函数 perf_eva(y, pred, show_plot = F, type = "roc"),计算效率还蛮高的,1.3万的样本0.04秒,大家可以试试看