Wie kann ich gruppierte Daten mit broom und dplyr auf gruppierte Modelle anwenden?

Question

Ich würde gerne das Äquivalent der Anpassung eines Modells von gpm (Gallonen pro Meile = 1/mpg) an wt im mtcars-Datensatz tun. Das scheint einfach:

data(mtcars) 
library(dplyr) 
library(tidyr) 
library(broom) 
library(ggplot2) 
library(scales) 

mtcars2 <- 
    mtcars %>% 
    mutate(gpm = 1/mpg) %>% 
    group_by(cyl, am) 

lm1 <- 
    mtcars2 %>% 
    do(fit = lm(gpm ~ wt, data = .)) 

Das bringt mir einen reihenweisen Datenrahmen mit 6 Zeilen, wie erwartet.

Dieser Graph bestätigt, dass es sechs Gruppen:

p1 <- 
    qplot(wt, gpm, data = mtcars2) + 
    facet_grid(cyl ~ am) + 
    stat_smooth(method='lm',se=FALSE, fullrange = TRUE) + 
    scale_x_continuous(limits = c(0,NA)) 

kann ich Augmentation() verwenden, um den Einbau Ausgänge zu erhalten:

lm1 %>% augment(fit) 

Das ist mir 32 Zeilen gibt, eine für jede Zeile in mtcars2, wie erwartet.

Nun ist die Herausforderung: Ich möchte Einbau Ausgänge mit newdata bekommen, wo ich wt von Zyl erhöht haben/4:

newdata <- 
    mtcars2 %>% 
    mutate(
     wt = wt + cyl/4) 

Ich gehe davon aus, dass dies einen Datenrahmen der gleichen Größe produzieren als lm1%>% augment (fit): eine Zeile für jede Zeile in newdata, da broom Modelle und newdata mit den Gruppierungsvariablen cyl und am abgleicht.

Leider

pred1 <- 
    lm1 %>% 
    augment(
     fit, 
     newdata = newdata) 

gibt mir einen Datenrahmen mit 192 Zeilen (= 6 x 32), scheinbar jedes Modell newdata an jede Zeile des Einpassen.

Aus dem Lesen von woanders, versichere ich, dass group_by und zeilenweise Datenrahmen nicht kompatibel sind, so lm1 ist ungruppiert, und Augment kann keine Modelle und newdata zuordnen. Gibt es ein anderes Designmuster, das mir dies ermöglicht? Es wäre schön, wenn es so einfach und transparent wäre wie der obige Versuch, aber es ist wichtiger, dass es funktioniert.

Hier ist meine Session():

> sessionInfo() 
R version 3.3.1 (2016-06-21) 
Platform: x86_64-w64-mingw32/x64 (64-bit) 
Running under: Windows 7 x64 (build 7601) Service Pack 1 

locale: 
[1] LC_COLLATE=English_United States.1252 
[2] LC_CTYPE=English_United States.1252 
[3] LC_MONETARY=English_United States.1252 
[4] LC_NUMERIC=C       
[5] LC_TIME=English_United States.1252  

attached base packages: 
[1] stats  graphics grDevices utils  datasets methods base  

other attached packages: 
[1] scales_0.4.0 ggplot2_2.1.0 broom_0.4.1 tidyr_0.6.0 dplyr_0.5.0 

loaded via a namespace (and not attached): 
[1] Rcpp_0.12.7  magrittr_1.5  mnormt_1.5-4  munsell_0.4.3 
[5] colorspace_1.2-6 lattice_0.20-34 R6_2.1.3   stringr_1.1.0 
[9] plyr_1.8.4  tools_3.3.1  parallel_3.3.1 grid_3.3.1  
[13] nlme_3.1-128  gtable_0.2.0  psych_1.6.9  DBI_0.5-1  
[17] lazyeval_0.2.0 assertthat_0.1 tibble_1.2  reshape2_1.4.1 
[21] labeling_0.3  stringi_1.1.1 compiler_3.3.1 foreign_0.8-67 

EDIT:

@aosmith: Ich habe Ihre zweite Option erforscht, und ich mag es. Wenn ich es jedoch auf meinen realen Daten versuche, habe ich ein Problem im muate-Befehl: es gibt "Fehler: augment weiß nicht, wie man mit Daten der Klassenliste umgeht" zurück.

Mein richtige Code ist mehr wie:

newdata %>% 
dplyr::select(cyl, am, wt) %>% # wt holds new predictor values 
group_by(cyl, am) %>% 
nest() %>% 
inner_join(regressions, .) %>% 
## looks like yours at this point 
mutate(pred = list(augment(fit, newdata = data))) %>% # Error here 
unnest(pred) 

Wo ich sagen, es ist wie das Ihre aussieht, ich meine, ich habe die folgenden Spalten (hier für Konsistenz umbenannt): ID (chr), attr1 (dbl), Zyl (dbl), am (chr), passen (Liste) und Daten (Liste). Sie haben Zyl, Am (dbl), Fit und Daten. Ich änderte mein bin zu dbl, aber das half nicht.

Ich denke, der Unterschied ist, dass ich 3 (ID ... ähnlich den rownames in mtcars) x 2 (cyl) x 2 (am) Einheiten in diesem Beispiel (mit jeder Probe mit 12 Messungen) haben, während die mtcars Beispiel hat 3 (cyl) x 2 (am) Zellen x eine zufällige Anzahl von Autotypen pro Zelle. In meiner Analyse muss ich die ID-Werte sehen, aber newdata gilt gleichermaßen für alle Einheiten. Wenn es hilft, denke an es als die Geschwindigkeit eines Gegenwindes, der auf jedes Auto im Test angewendet wird. Ist dies ein Grund für die Beschwerde von augment, dass sie nicht mit Daten der Klassenliste umgehen kann?

EDIT: Das Zusammenführen der ID mit den newdata (mit Full = TRUE) löste das letzte Problem.Ich verwende derzeit Ihre erste vorgeschlagene Lösung.

Answer 1

Ich habe map2 aus Paket purrr für diese Art von Situation verwendet. map2 durchläuft die Elemente von zwei Listen gleichzeitig. Die Listen müssen die gleiche Länge haben und in der gleichen Reihenfolge sein.

Die Elemente der Listen werden als Argumente für einige Funktionen verwendet, die Sie anwenden möchten (augment, in Ihrem Fall). Hier wären Ihre zwei Listen eine Liste von Modellen und eine Liste von Datensätzen (eine Liste für jede cyl/ Kombination).

Mit map2_df werden die Ergebnisse als data.frame anstelle einer Liste zurückgegeben.

library(purrr) 

Ich habe die Liste der data.frames mit der Verwendung von split vorherzusagen. Die Reihenfolge der zu teilenden Faktoren bestimmte die Reihenfolge der Liste, so dass ich sicherstellte, dass sie in der gleichen Reihenfolge wie lm1 lag.

test_split = split(newdata, list(newdata$am, newdata$cyl) 

map2_df(lm1$fit, test_split, ~augment(.x, newdata = .y)) 

Um zu vermeiden, so viel über, um sich Gedanken, um die Vorhersagedaten von Gruppen nest könnte, kommen diese zu lm1, und die Ergebnisse der augment als Liste zurück zur Entschachtelung.

newdata %>% 
    group_by(cyl, am) %>% 
    nest() %>% 
    inner_join(lm1, .) %>% 
    mutate(pred = list(augment(fit, newdata = data))) %>% 
    unnest(pred)