2013-09-05 2 views
9

Ich arbeite mit mehrdimensionalen Array sowohl auf R und MATLAB, diese Arrays haben fünf Dimensionen (insgesamt 14,5 M Elemente). Ich muss eine Dimension entfernen, die ein arithmetisches Mittel anwendet, und ich entdeckte einen erstaunlichen Unterschied der Leistungen, die die zwei Software verwenden.Arithmetisches Mittel auf einem mehrdimensionalen Array auf R und MATLAB: drastische Differenz der Leistungen

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 

Ich weiß, dass anwenden Funktion langsam ist, da so etwas wie eine Allzweck-Funktion ist, aber ich weiß nicht, wie mit diesem Problem beschäftigen weil dieser Leistungsunterschied wirklich eine große Grenze für mich darstellt. Ich habe versucht, nach einer Verallgemeinerung von colMeans/rowMeans Funktionen zu suchen, aber es ist mir nicht gelungen.

EDIT ich eine kleine Probenmatrix zeigen werden:

> 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):

ich entdeckte, dass Summenfunktion anwendet und dann durch die Größe des entfernten Dividieren Dimension ist definitiv schneller:

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

Können Sie reduzieren die Eingabe/gewünschte Ausgabe bewegen eine kleine 3-dimensionale Anordnung für illustrative Zwecke, z eine 3 * 3 * 2 Matrix? –

+0

@Matteodefelice siehe http://StackOverflow.com/Questions/18604406/Why-is-mean-so-Slow vor allem Joshua die Antwort in Bezug auf Präzision. –

Antwort

5

mean ist wegen des S3-Methodenversands besonders langsam. Dies ist schneller:

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

Zum Vergleich:

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 

Wenn Sie nicht NA s behandeln müssen Sie die interne Funktion können

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 

(Beachten Sie, dass alle Timings sind nur ungefähre Werte, und ein korrektes Benchmarking würde es erfordern, dies wiederholt auszuführen, zB mit einem der Bechmarking-Pakete. Aber dafür bin ich nicht geduldig genug.)

Es könnte möglich sein, dies mit einer Rcpp-Implementierung zu beschleunigen.

+1

[** Siehe hier **] (http://Stackoverflow.com/a/18604487/1478381) für weitere Informationen. –

+0

Ich habe auch versucht 'Bibliothek (data.table); system.time ({b3 = apply (a, c (1,2,4,5), Funktion (x) .External ("Cfastmean", x, FALSE))}) ', aber es war nicht schneller. – Roland

+0

Danke, du hast mir definitiv geholfen! Ich habe auch etwas sehr Interessantes über R interne Mechanismen gelernt ... –

20

In R, apply ist nicht das richtige Werkzeug für die Aufgabe. Wenn Sie eine Matrix hätten und die Zeilen- oder Spaltenmittel benötigen, würden Sie die viel viel schneller, vektorisierten rowMeans und colMeans verwenden.Sie können immer noch diese für ein mehrdimensionales Array verwenden, aber Sie müssen ein wenig kreativ sein:

Array Unter der Annahme hat n Dimensionen, und Sie wollen Mittel entlang Dimension berechnen i:

  1. Verwendung aperm zu bewegen, um die Dimension i auf die letzte Position n
  2. Verwendung rowMeans mit dims = n - 1

, Sie Ähnlich könnte:

  1. Verwendung aperm die Dimension i in die erste Position
  2. Verwendung colMeans mit 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 
+0

Genau. Verwenden Sie vektorisierte Funktionen immer über Schleifen oder * wenden Sie sie nach Möglichkeit an. –

+0

Dies sollte * definitiv * die akzeptierte Antwort sein. +1 für eine großartige Erklärung des verallgemeinerten Falles. Ich sah "aperm" aus, konnte es aber einfach nicht richtig machen. Vielen Dank! –

+0

Der Vollständigkeit halber verwendet 'rowMeans' nicht den gleichen Algorithmus' mean'; der erste ist naive Single-Pass-Akkumulation und -Teilung; Letzteres hat einen zweiten Durchgang, um die numerische Stabilität zu verbessern. –

Verwandte Themen