2010-09-10 21 views
35

Tengo una simulación que tiene un gran agregado y un paso combinado en el medio. Creé un prototipo de este proceso utilizando la función ddply() de plyr, que funciona muy bien para un gran porcentaje de mis necesidades. Pero necesito que este paso de agregación sea más rápido ya que tengo que ejecutar simulaciones de 10K. Ya estoy escalando las simulaciones en paralelo, pero si este paso fuera más rápido, podría disminuir en gran medida la cantidad de nodos que necesito.R: acelerar las operaciones de "agrupar por"

Aquí es una simplificación razonable de lo que estoy tratando de hacer:

library(Hmisc) 

# Set up some example data 
year <- sample(1970:2008, 1e6, rep=T) 
state <- sample(1:50, 1e6, rep=T) 
group1 <- sample(1:6, 1e6, rep=T) 
group2 <- sample(1:3, 1e6, rep=T) 
myFact <- rnorm(100, 15, 1e6) 
weights <- rnorm(1e6) 
myDF <- data.frame(year, state, group1, group2, myFact, weights) 

# this is the step I want to make faster 
system.time(aggregateDF <- ddply(myDF, c("year", "state", "group1", "group2"), 
        function(df) wtd.mean(df$myFact, weights=df$weights) 
           ) 
      ) 

Todos los consejos o sugerencias son apreciados!

+1

No relacionado con el rendimiento, pero pago 'weighted.mean' en la base – hadley

+1

Oh, eso es práctico. Puedes ver que aprendí R buscando en Google lo que necesito hacer;) –

Respuesta

37

En lugar de la trama de datos R normal, puede utilizar una trama de datos inmutables que devuelve los punteros a su posición original al subconjunto y puede ser mucho más rápido:

idf <- idata.frame(myDF) 
system.time(aggregateDF <- ddply(idf, c("year", "state", "group1", "group2"), 
    function(df) wtd.mean(df$myFact, weights=df$weights))) 

# user system elapsed 
# 18.032 0.416 19.250 

Si tuviera que escribir una plyr función personalizada exactamente a esta situación, me gustaría hacer algo como esto:

system.time({ 
    ids <- id(myDF[c("year", "state", "group1", "group2")], drop = TRUE) 
    data <- as.matrix(myDF[c("myFact", "weights")]) 
    indices <- plyr:::split_indices(seq_len(nrow(data)), ids, n = attr(ids, "n")) 

    fun <- function(rows) { 
    weighted.mean(data[rows, 1], data[rows, 2]) 
    } 
    values <- vapply(indices, fun, numeric(1)) 

    labels <- myDF[match(seq_len(attr(ids, "n")), ids), 
    c("year", "state", "group1", "group2")] 
    aggregateDF <- cbind(labels, values) 
}) 

# user system elapsed 
# 2.04 0.29 2.33 

es mucho más rápido, ya que evita la copia de los datos, sólo extraer el subconjunto necesario para el cómputo de cada cuando es computado Cambiar los datos a la forma de matriz proporciona otro aumento de velocidad debido a que el subconjunto de matriz es mucho más rápido que el subconjunto de marcos de datos.

+3

'idata.frame' se agregó en plyr 1.0. – hadley

+0

Me había metido con hacer índices y cosas así con data.table y casi había renunciado a esa idea. Esperaba un 50% de mejora. Esto excede por mucho mis expectativas. –

+0

teniendo un pequeño problema para hacer que esto funcione bien ... Pero estoy aprendiendo sobre la marcha ... He cambiado los datos a myDF pero no estoy seguro de dónde está el problema ... –

7

¿Está utilizando la última versión de plyr (nota: esto todavía no ha llegado a todos los espejos CRAN)? Si es así, podrías ejecutar esto en paralelo.

Aquí está el ejemplo llply, pero el mismo debe aplicarse a ddply:

x <- seq_len(20) 
    wait <- function(i) Sys.sleep(0.1) 
    system.time(llply(x, wait)) 
    # user system elapsed 
    # 0.007 0.005 2.005 

    library(doMC) 
    registerDoMC(2) 
    system.time(llply(x, wait, .parallel = TRUE)) 
    # user system elapsed 
    # 0.020 0.011 1.038 

Editar:

Bueno, otros enfoques de bucle son peores, por lo que probablemente requiere ya sea (a) C/Código C++ o (b) un replanteamiento más fundamental de cómo lo está haciendo. Ni siquiera intenté usar by() porque eso es muy lento en mi experiencia.

groups <- unique(myDF[,c("year", "state", "group1", "group2")]) 
system.time(
aggregateDF <- do.call("rbind", lapply(1:nrow(groups), function(i) { 
    df.tmp <- myDF[myDF$year==groups[i,"year"] & myDF$state==groups[i,"state"] & myDF$group1==groups[i,"group1"] & myDF$group2==groups[i,"group2"],] 
    cbind(groups[i,], wtd.mean(df.tmp$myFact, weights=df.tmp$weights)) 
})) 
) 

aggregateDF <- data.frame() 
system.time(
for(i in 1:nrow(groups)) { 
    df.tmp <- myDF[myDF$year==groups[i,"year"] & myDF$state==groups[i,"state"] & myDF$group1==groups[i,"group1"] & myDF$group2==groups[i,"group2"],] 
    aggregateDF <- rbind(aggregateDF, data.frame(cbind(groups[i,], wtd.mean(df.tmp$myFact, weights=df.tmp$weights)))) 
} 
) 
+0

que me ayuda en el caso de una sola máquina, pero ya estoy explotando en Hadoop y sobrescribiendo cada nodo (más procesos que procesadores). ¡Pero estoy muy contento de ver la paralelización convirtiéndolo en plyr! –

8

me habría perfil con base de R

g <- with(myDF, paste(year, state, group1, group2)) 
x <- with(myDF, c(tapply(weights * myFact, g, sum)/tapply(weights, g, sum))) 
aggregateDF <- myDF[match(names(x), g), c("year", "state", "group1", "group2")] 
aggregateDF$V1 <- x 

En mi máquina se tarda 5 segundos para comparar con 67sec código original.

EDITAR acaba de encontrar otra velocidad con rowsum función:

g <- with(myDF, paste(year, state, group1, group2)) 
X <- with(myDF, rowsum(data.frame(a=weights*myFact, b=weights), g)) 
x <- X$a/X$b 
aggregateDF2 <- myDF[match(rownames(X), g), c("year", "state", "group1", "group2")] 
aggregateDF2$V1 <- x 

Se tarda 3 segundos!

+2

El segundo tarda 5 segundos en mi computadora, por lo que plyr sigue golpeando la base;) (Además ordena las filas correctamente) – hadley

+2

Pero gracias por el puntero a 'rowsum' - es tan difícil mantenerse al día con la plétora de funciones de agregación en la base R. – hadley

+0

Sabía que tenía que ser una forma tapply de hacer esto también, pero estaba luchando para resolverlo. En general, tengo esta lucha con la familia de aplicar. –

5

generalmente utilizo un vector de índice con tapply cuando se aplica la función tiene múltiples args vector:

system.time(tapply(1:nrow(myDF), myDF[c('year', 'state', 'group1', 'group2')], function(s) weighted.mean(myDF$myFact[s], myDF$weights[s]))) 
# user system elapsed 
# 1.36 0.08 1.44 

I utilizan una envoltura simple que es equivalente pero esconde el desorden:

tmapply(list(myDF$myFact, myDF$weights), myDF[c('year', 'state', 'group1', 'group2')], weighted.mean) 

Editado para incluir tmapply para comentarios a continuación:

tmapply = function(XS, INDEX, FUN, ..., simplify=T) { 
    FUN = match.fun(FUN) 
    if (!is.list(XS)) 
    XS = list(XS) 
    tapply(1:length(XS[[1L]]), INDEX, function(s, ...) 
    do.call(FUN, c(lapply(XS, `[`, s), list(...))), ..., simplify=simplify) 
} 
+0

que es muy ingenioso para ver eso hecho en la base R. ¡Gracias! –

+1

Solo para agregar: 'as.data.frame (as.table (RESULTS))' es una manera fácil de crear 'data.frame' desde la salida. – Marek

+0

¿Es este el 'tmapply' que estás usando? https://stat.ethz.ch/pipermail/r-help/2002-October/025773.html – Shane

25

Además aceleración 2x ​​y el código más conciso:

library(data.table) 
dtb <- data.table(myDF, key="year,state,group1,group2") 
system.time( 
    res <- dtb[, weighted.mean(myFact, weights), by=list(year, state, group1, group2)] 
) 
# user system elapsed 
# 0.950 0.050 1.007 

Mi primer post, así que por favor ser agradable;)


De data.table v1.9.2, setDT función se exporta que va a convertirlos en data.frame a data.tablepor referencia (según el lenguaje data.table - todas las funciones set* modifican el objeto por referencia). Esto significa que no hay copia innecesaria y, por lo tanto, es rápido. Puedes cronometrarlo, pero será negligente.

require(data.table) 
system.time({ 
    setDT(myDF) 
    res <- myDF[, weighted.mean(myFact, weights), 
      by=list(year, state, group1, group2)] 
}) 
# user system elapsed 
# 0.970 0.024 1.015 

Esto es tan opuesto a 1.264 segundos con la solución de OP de arriba, donde data.table(.) se utiliza para crear dtb.

+0

¡Buena publicación! Gracias por la respuesta. Sin embargo, para ser coherente con los otros métodos, el paso que crea la tabla de datos y el índice debe estar dentro del paso system.time(). –

+2

De hecho, pero sigue siendo el más rápido. Sería bueno tener una opción en ddply para operar en data.tables o usar data.tables bajo el capó (Acabo de descubrir data.table buscando soluciones para el mismo problema, pero preferiría una versión más ddply-like sintaxis para este caso). – datasmurf