R: acelerar las operaciones de “agrupar por”

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.

Aquí hay 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!

En lugar del dataframe R normal, puede usar un dataframe inmutables que devuelve punteros al original cuando se subconjunta 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 función plyr personalizada para esta situación, haría 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 porque evita copiar los datos, solo extrae el subconjunto necesario para cada cómputo cuando se calcula. 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.

Más velocidad 2x y un 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 primera publicación, así que por favor sé amable;)


De data.table v1.9.2, se setDT función setDT que convertirá data.frame a data.table por referencia (de acuerdo con el lenguaje de 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 opuesto a 1.264 segundos con la solución de OP anterior, donde data.table(.) Se usa para crear dtb .

Me gustaría perfil con la base 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, toma 5 segundos en comparación con 67 segundos con el código original.

EDITAR Acaba de encontrar otra velocidad con la función de rowsum :

 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 

¡Toma 3 segundos!

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

Aquí está el ejemplo llply, pero lo mismo debería 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 requiera (a) código C / 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)))) } ) 

Usualmente uso un vector de índice con tapply cuando la función que se está aplicando tiene múltiples args vectoriales:

 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 

Utilizo una envoltura simple que es equivalente pero oculta el desastre:

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