R Language
Snabba upp svår-att-vektorisera koden

Påskynda för svår att vektorisera för slingor med Rcpp

Tänk på följande tufft att vektorisera för slinga, som skapar en vektor med längden len där det första elementet anges ( first ) och varje element x_i är lika med 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)
}

Denna kod involverar en for-loop med en snabb operation ( cos(x[i-1]+1) ), som ofta har nytta av vektorisering. Det är emellertid inte trivialt att vektorisera denna operation med bas R, eftersom R inte har en "kumulativ cosinus med x + 1" -funktion.

En möjlig metod för att påskynda denna funktion skulle vara att implementera den i C ++ med Rcpp-paketet:

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;
}")

Detta ger ofta betydande hastigheter för stora beräkningar och ger exakt samma resultat:

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 

I detta fall genererar Rcpp-koden en vektor med en längd på 1 miljon på 0,03 sekunder istället för 1,31 sekunder med bas R-metoden.

Påskynda tuff-till-vektorisera för slingor genom byte-sammanställning

Följ Rcpp-exemplet i den här dokumentationsposten, överväg följande tuff-till-vektoriseringsfunktion, som skapar en vektor med längd len där det första elementet anges ( first ) och varje element x_i är lika med 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)
}

Ett enkelt tillvägagångssätt för att påskynda en sådan funktion utan att skriva om en enda kodrader är att byta samman koden med R-kompileringspaketet:

library(compiler)
repeatedCosPlusOneCompiled <- cmpfun(repeatedCosPlusOne)

Den resulterande funktionen är ofta betydligt snabbare medan du fortfarande ger samma resultat:

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 

I detta fall rusade byte-kompilering upp hård-till-vektoriseringsoperationen på en vektor med längden 1 miljon från 1,20 sekunder till 0,34 sekunder.

Anmärka

Kärnan i repeatedCosPlusOne , som den kumulativa tillämpningen av en enda funktion, kan uttryckas mer öppet med 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 kan betraktas som en "vektorisering" av repeatedCosPlusOne . Det kan dock förväntas bli långsammare med en faktor 2:

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