R使用Rcpp加速难以向量化的循环

示例

考虑以下难以向量化的for循环,该循环创建一个长度为矢量的len地方,其中第一个元素被指定为(first),每个元素x_i等于cos(x_{i-1} + 1):

repeatedCosPlusOne <- function(first, len) {
  x <- numeric(len)
  x[1] <- first
  for (i in 2:len) {
    x[i] <- cos(x[i-1] + 1)
  }
  return(x)
}

此代码涉及带有快速操作(cos(x[i-1]+1))的for循环,通常可从向量化中受益。但是,由于R不具有“ x + 1的累积余弦”功能,因此用基数R对这个操作进行矢量化并非易事。

加快此功能的一种可能方法是使用Rcpp软件包在C ++中实现它:

library(Rcpp)
cppFunction("NumericVector repeatedCosPlusOneRcpp(double first, int len) {
  NumericVector x(len);
  x[0] = first;
  for (int i=1; i < len; ++i) {
    x[i] = cos(x[i-1]+1);
  }
  return x;
}")

这通常可为大型计算提供显着的加速,同时产生完全相同的结果:

all.equal(repeatedCosPlusOne(1, 1e6), repeatedCosPlusOneRcpp(1, 1e6))
# [1] TRUE
system.time(repeatedCosPlusOne(1, 1e6))
#    user  system elapsed 
#   1.274   0.015   1.310 
system.time(repeatedCosPlusOneRcpp(1, 1e6))
#    user  system elapsed 
#   0.028   0.001   0.030

在这种情况下,Rcpp代码将在0.03秒内生成长度为100万的向量,而不是使用基本R方法生成的长度为1.31秒。