使用自定义函数加速 rollapply 的技术 (r)

Question

我有一个 rollapply 函数，它做一些非常简单的事情，但是超过百万个数据点，这个简单的函数非常慢。我想知道是否可以提供信息来 rollapply 以了解如何进行下一次转换，而不是定义函数本身。

具体来说，我正在为基本的统计异常检测执行滚动窗口。

滚动应用功能：

minmax <- function(x) { max(x) - min(x) }

调用者：

mclapply(data[,eval(vars),with=F], 
         function(x) rollapply(x,width=winSize,FUN=minmax,fill=NA),
         mc.cores=8)

data 是 8 列 data.table，winsize 是 300

此调用在 8 个内核上大约需要 2 分钟。它是整体计算的主要瓶颈之一。但是我可以想象我们可以保持它们排序（按值和索引），然后在每次滑动时进行 Olog(n) 比较。

但是，我经常看到建议不要使用 for 循环并使用 lapply 的帖子。进一步优化的下一个合乎逻辑的步骤是什么？

Answer 1

不确定这是否/如何应用于mclapply 环境，但您可以通过使用zoo 的优化rollmax 函数来获得一点加速。由于他们没有补充 rollmin，因此您需要适应。

minmax <- function(x) max(x) - min(x)
aa <- runif(1e4)
identical(
  zoo::rollapply(aa, width=100, FUN=minmax, fill=NA),
  zoo::rollmax(aa, k=100, fill=NA) + zoo::rollmax(-aa, k=100, fill=NA)
)
# [1] TRUE

microbenchmark::microbenchmark(
  minmax = zoo::rollapply(aa, width=100, FUN=minmax, fill=NA),
  dblmax = zoo::rollmax(aa, k=100, fill=NA) + zoo::rollmax(-aa, k=100, fill=NA)
)
# Unit: milliseconds
#    expr     min      lq     mean   median      uq      max neval
#  minmax 70.7426 76.0469 84.81481 77.99565 81.8047 148.8431   100
#  dblmax 15.6755 17.4501 19.09820 17.93665 18.8650  52.4849   100

（改进将取决于窗口大小，因此您的结果可能会有所不同，但我认为使用优化的函数zoo::rollmax 几乎总是优于每次调用 UDF。）

Answer 2

如果您真的想尽可能提高性能，请使用 Rcpp。自定义循环是 C++ 的一个很好的用例，尤其是当您的函数非常简单时。

先结果再编码：

microbenchmark::microbenchmark(
  minmax = zoo::rollapply(aa, width=100, FUN=minmax, fill=NA),
  dblmax = zoo::rollmax(aa, k=100, fill=NA) + zoo::rollmax(-aa, k=100, fill=NA),
  cminmax = crollapply(aa, width=width), times = 10
)
    Unit: milliseconds
    expr       min         lq       mean    median         uq        max neval cld
  minmax 154.04630 162.728871 188.198416 173.13427 200.928005 298.568673    10   c
  dblmax  37.38127  38.541603  44.818505  41.42796  50.001888  61.024250    10  b 
 cminmax   2.31766   2.363676   2.406835   2.39237   2.438109   2.512162    10 a  

C++/Rcpp 代码：

#include <Rcpp.h>
#include <algorithm>
using namespace Rcpp;

// [[Rcpp::export]]
std::vector<double> crollapply(std::vector<double> aa, int width) {
  if(width > aa.size()) throw exception("width too large :(");
  int start_offset = (width-1) / 2;
  int back_offset = width / 2;
  std::vector<double> results(aa.size());
  int i=0;
  for(; i < start_offset; i++) {
    results[i] = NA_REAL;
  }
  for(; i < results.size() - back_offset; i++) {
    double min = *std::min_element(&aa[i - start_offset], &aa[i + back_offset + 1]);
    double max = *std::max_element(&aa[i - start_offset], &aa[i + back_offset + 1]);
    results[i] = max - min;
  }
  for(; i < results.size(); i++) {
    results[i] = NA_REAL;
  }
  return results;
}

R 代码：

library(dplyr)
library(zoo)
library(microbenchmark)
library(Rcpp)

sourceCpp("~/Desktop/temp.cpp")

minmax <- function(x) max(x) - min(x)
aa <- runif(1e4)
width <- 100
x1 <- zoo::rollapply(aa, width=width, FUN=minmax, fill=NA)
x3 <- crollapply(aa, width=width)
identical(x1,x3)

width <- 101
x1 <- zoo::rollapply(aa, width=width, FUN=minmax, fill=NA)
x3 <- crollapply(aa, width=width)
identical(x1,x3)

microbenchmark::microbenchmark(
  minmax = zoo::rollapply(aa, width=100, FUN=minmax, fill=NA),
  dblmax = zoo::rollmax(aa, k=100, fill=NA) + zoo::rollmax(-aa, k=100, fill=NA),
  cminmax = crollapply(aa, width=width), times = 10
)