Media aritmetica su un array multidimensionale su R e MATLAB: drastica differenza di prestazioni

Question

Sto lavorando con array multidimensionale sia su R che MATLAB, questi array hanno cinque dimensioni (totale di 14,5 M di elementi). Devo rimuovere una dimensione applicando una media aritmetica e ho scoperto una sorprendente differenza di prestazioni usando i due software.

MATLAB:

>> a = rand([144 73 10 6 23]); 
>> tic; b = mean(a,3); toc 
Elapsed time is 0.014454 seconds. 

R:

> a = array(data = runif(144*73*6*23*10), dim = c(144,73,10,6,23)) 
> start <- Sys.time(); b = apply(a, c(1,2,4,5), mean); Sys.time() - start 
Time difference of 1.229083 mins 

So che applicano funzione è lento perché è qualcosa di simile a una funzione di uso generale, ma non so come affrontare questo problema perché questa differenza di prestazioni è davvero un grosso limite per me. Ho provato a cercare una generalizzazione delle funzioni colMeans/rowMeans ma non ci sono riuscito.

EDIT vi mostrerò un po matrice del campione:

> dim(a) 
[1] 2 4 3 
> dput(aa) 
structure(c(7, 8, 5, 8, 10, 11, 9, 9, 6, 12, 9, 10, 12, 10, 14, 
12, 7, 9, 8, 10, 10, 9, 8, 6), .Dim = c(2L, 4L, 3L)) 
a_mean = apply(a, c(2,3), mean) 
> a_mean 
    [,1] [,2] [,3] 
[1,] 7.5 9.0 8.0 
[2,] 6.5 9.5 9.0 
[3,] 10.5 11.0 9.5 
[4,] 9.0 13.0 7.0 

EDIT (2):

ho scoperto che l'applicazione della funzione somma e quindi dividendo per la dimensione del rimosso dimensione è decisamente più veloce:

> start <- Sys.time(); aaout = apply(aa, c(1,2,4,5), sum); Sys.time() - start 
Time difference of 5.528063 secs 

Answer 1

mean è particolarmente lento a causa dell'invio del metodo S3. Questo è più veloce:

set.seed(42) 
a = array(data = runif(144*73*6*23*10), dim = c(144,73,10,6,23)) 

system.time({b = apply(a, c(1,2,4,5), mean.default)}) 
# user system elapsed 
#16.80 0.03 16.94 

Se non è necessario gestire NA s è possibile utilizzare la funzione interna:

system.time({b1 = apply(a, c(1,2,4,5), function(x) .Internal(mean(x)))}) 
# user system elapsed 
# 6.80 0.04 6.86 

Per fare un confronto:

system.time({b2 = apply(a, c(1,2,4,5), function(x) sum(x)/length(x))}) 
# user system elapsed 
# 9.05 0.01 9.08 

system.time({b3 = apply(a, c(1,2,4,5), sum) 
      b3 = b3/dim(a)[[3]]}) 
# user system elapsed 
# 7.44 0.03 7.47 

(Si noti che tutti i tempi sono solo approssimativi. Un corretto benchmarking richiederebbe di eseguire questa operazione in modo ripetitivo, ad esempio, utilizzando uno dei pacchetti bechmarking, ma non sono abbastanza paziente per quello al momento.)

Potrebbe essere possibile velocizzare questa implementazione con un'implementazione Rcpp.

Answer 2

In R, apply non è lo strumento giusto per il compito. Se avevi una matrice e avevi bisogno della riga o della colonna significa che useresti il ​​molto più veloce, vettorizzato rowMeans e colMeans.È comunque possibile utilizzare questi per un array multi-dimensionale, ma è necessario essere un po 'creativo:

Assumendo che l'array ha n dimensioni, e si vuole calcolare mezzi lungo dimensione i:

1. uso aperm a spostare la dimensione i all'ultima posizione n
2. uso rowMeans con dims = n - 1

Allo stesso modo, si potrebbe:

1. uso aperm per spostare la dimensione i alla prima posizione
2. uso colMeans con dims = 1

---

a <- array(data = runif(144*73*6*23*10), dim = c(144,73,10,6,23)) 

means.along <- function(a, i) { 
    n <- length(dim(a)) 
    b <- aperm(a, c(seq_len(n)[-i], i)) 
    rowMeans(b, dims = n - 1) 
} 

system.time(z1 <- apply(a, c(1,2,4,5), mean)) 
# user system elapsed 
# 25.132 0.109 25.239 
system.time(z2 <- means.along(a, 3)) 
# user system elapsed 
# 0.283 0.007 0.289 

identical(z1, z2) 
# [1] TRUE