R Language
厳しいからベクトル化するコードの高速化
サーチ…
Rcppを使用したループのタフなベクトル化の高速化
最初の要素が指定され( first
)、各要素x_i
がcos(x_{i-1} + 1)
等しい長さlen
ベクトルを作成する次のtough-to-vectorize forループを考えてみましょう。
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コードは、基底Rのアプローチでは1.31秒ではなく、0.03秒で100万の長さのベクトルを生成します。
バイトコンパイルによるループの難易度の高いベクトル化の高速化
最初の要素が指定され( first
)、各要素x_i
がcos(x_{i-1} + 1)
等しい長さlen
ベクトルを作成する次の難しいtoベクトル化関数を考えてみましょう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)
}
1行のコードを書き直さずにこのような関数を高速化する簡単な方法の1つは、Rコンパイルパッケージを使用してコードをバイトコンパイルすることです。
library(compiler)
repeatedCosPlusOneCompiled <- cmpfun(repeatedCosPlusOne)
結果の関数は、しばしば同じ結果を返しながら、かなり速くなります:
all.equal(repeatedCosPlusOne(1, 1e6), repeatedCosPlusOneCompiled(1, 1e6))
# [1] TRUE
system.time(repeatedCosPlusOne(1, 1e6))
# user system elapsed
# 1.175 0.014 1.201
system.time(repeatedCosPlusOneCompiled(1, 1e6))
# user system elapsed
# 0.339 0.002 0.341
この場合、バイトコンパイルは、1.20秒から0.34秒までの長さ1百万のベクトル上の厳しいto-vectorize操作をスピードアップしました。
リマーク
repeatedCosPlusOne
の本質は、単一の関数の累積的なアプリケーションとして、 Reduce
より透過的に表現できます。
iterFunc <- function(init, n, func) {
funcs <- replicate(n, func)
Reduce(function(., f) f(.), funcs, init = init, accumulate = TRUE)
}
repeatedCosPlusOne_vec <- function(first, len) {
iterFunc(first, len - 1, function(.) cos(. + 1))
}
repeatedCosPlusOne_vec
は、 repeatedCosPlusOne
「ベクトル化」とみなすことができます。しかし、それは2分の1で遅くなることが予想される:
library(microbenchmark)
microbenchmark(
repeatedCosPlusOne(1, 1e4),
repeatedCosPlusOne_vec(1, 1e4)
)
#> Unit: milliseconds
#> expr min lq mean median uq max neval cld
#> repeatedCosPlusOne(1, 10000) 8.349261 9.216724 10.22715 10.23095 11.10817 14.33763 100 a
#> repeatedCosPlusOne_vec(1, 10000) 14.406291 16.236153 17.55571 17.22295 18.59085 24.37059 100 b