R Language
Accelerare il codice difficile da vectorize

Eccesso di velocità per vettorizzare per cicli con Rcpp

Si consideri il seguente ciclo for to vectorize per il ciclo, che crea un vettore di lunghezza len cui viene specificato il primo elemento ( first ) e ciascun elemento x_i è uguale a 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)
}

Questo codice implica un ciclo for con un'operazione rapida ( cos(x[i-1]+1) ), che spesso beneficiano della vettorizzazione. Tuttavia, non è banale il vettorizzare questa operazione con la base R, poiché R non ha una funzione "coseno cumulativo di x + 1".

Un possibile approccio per accelerare questa funzione sarebbe implementarlo in C ++, usando il pacchetto Rcpp:

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

Questo spesso fornisce significativi aumenti di velocità per i calcoli di grandi dimensioni, fornendo allo stesso tempo gli stessi risultati:

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 

In questo caso, il codice Rcpp genera un vettore di lunghezza 1 milione in 0,03 secondi anziché 1,31 secondi con l'approccio di base R.

Accelerazione difficile da vettorizzare per cicli con la compilazione dei byte

Seguendo l'esempio di Rcpp in questa voce di documentazione, si consideri la seguente funzione tough-to-vectorize, che crea un vettore di lunghezza len cui viene specificato il primo elemento ( first ) e ciascun elemento x_i è uguale a 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)
}

Un semplice approccio per accelerare tale funzione senza riscrivere una singola riga di codice è il byte che compila il codice usando il pacchetto di compilazione R:

library(compiler)
repeatedCosPlusOneCompiled <- cmpfun(repeatedCosPlusOne)

La funzione risultante sarà spesso significativamente più veloce mentre restituisce sempre gli stessi risultati:

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 

In questo caso, la compilazione di byte ha velocizzato l'operazione di difficile interpretazione su un vettore di lunghezza 1 milione da 1,20 secondi a 0,34 secondi.

osservazione

L'essenza di repeatedCosPlusOne , come l'applicazione cumulativa di una singola funzione, può essere espressa in modo più trasparente con 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 può essere considerato come una "vettorizzazione" di repeatedCosPlusOne . Tuttavia, ci si può aspettare che sia più lento di un fattore 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