Accelerare un'operazione di matrice di remi
Domanda
Considera la seguente matrice,
nc <- 5000
nr <- 1024
m <- matrix(rnorm(nc*nr), ncol=nc)
Desidero prendere la differenza tra il rowMeans
di due gruppi di dimensioni identiche prese a caso in questa matrice.
n <- 1000 # group size
system.time(replicate(100, {
ind1 <- sample(seq.int(nc), n)
ind2 <- sample(seq.int(nc), n)
rowMeans(m[, ind1]) - rowMeans(m[, ind2])
}))
È abbastanza lento, purtroppo non ho capito l'output di RProf (sembrava che la maggior parte del tempo fosse trascorso is.data.frame
??)
Suggerimenti per qualcosa di più efficiente?
Ho contemplato quanto segue:
Rcpp
: Dalle mie letture online credo che i remi di R siano abbastanza efficienti, quindi non è chiaro che aiuterebbe a questo passaggio. Mi piacerebbe essere convinto di dove il collo di bottiglia è davvero il primo, forse il mio intero design è non ottimale. Se la maggior parte del tempo viene trascorso per creare copie per ciascuna delle matrici più piccole, RCPP funzionerebbe meglio?Aggiornando a R-devel, sembra esserci un nuovo
.rowMeans
funzionare ancora più efficiente. Qualcuno l'ha provato?
Grazie.
Soluzione
A testa rowSums()
chiamare un sottoinsieme di colonne da m
può essere visto come la moltiplicazione della matrice tra m
e un vettore di 0
o 1
indicando le colonne selezionate. Se giustappondi tutti quei vettori, finisci con una moltiplicazione tra due matrici (che è molto più efficiente):
ind1 <- replicate(100, seq.int(nc) %in% sample(seq.int(nc), n))
ind2 <- replicate(100, seq.int(nc) %in% sample(seq.int(nc), n))
output <- m %*% (ind1 - ind2)
Altri suggerimenti
Non hai bisogno di 2 chiamate a rowMeans
. Puoi prima fare la sottrazione e chiamare rowMeans
sul risultato.
x1 <- rowMeans(m[,ind1])-rowMeans(m[,ind2])
x2 <- rowMeans(m[,ind1]-m[,ind2])
all.equal(x1,x2)
# [1] TRUE
is.data.frame
fa parte dei controlli effettuati in rowMeans
.
Aggiornamento: per quanto riguarda .rowMeans
A R-devel, sembra che sia solo una chiamata diretta al codice interno (supponendo do_colsum
non è cambiato). È definito come:
.rowMeans <- function(X, m, n, na.rm = FALSE)
.Internal(rowMeans(X, m, n, na.rm))
Nel tuo caso, m=1024
e n=1000
.