libreria - dplyr versión de agrupar un marco de datos y luego crear un modelo de regresión en cada grupo

dplyr en r

Question

¿Alguien puede sugerir una respuesta dplyr a la siguiente pregunta? Dividir data.frame por país y crear un modelo de regresión lineal en cada subconjunto

Para completar, la pregunta y la respuesta del enlace se incluyen a continuación.

Pregunta

Como referencia, aquí está la pregunta de Josh:

Tengo un data.frame de datos del Banco Mundial que se parece a esto;

  country date BirthRate     US.   
4   Aruba 2011    10.584 25354.8
5   Aruba 2010    10.804 24289.1
6   Aruba 2009    11.060 24639.9
7   Aruba 2008    11.346 27549.3
8   Aruba 2007    11.653 25921.3
9   Aruba 2006    11.977 24015.4

En general, hay 70 subconjuntos de países en este marco de datos en los que me gustaría ejecutar una regresión lineal. Si utilizo lo siguiente, obtengo un buen lm para un solo país;

andora = subset(high.sub, country == "Andorra")

andora.lm = lm(BirthRate~US., data = andora)

anova(andora.lm)
summary(andora.lm)

Pero cuando trato de usar el mismo tipo de código en un ciclo for, recibo un error que imprimiré debajo del código;

high.sub = subset(highInc, date > 1999 & date < 2012)
high.sub <- na.omit(high.sub)
highnames <- unique(high.sub$country)

for (i in highnames) {
  linmod <- lm(BirthRate~US., data = high.sub, subset = (country == "[i]"))  
}

#Error message:
Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) : 
  0 (non-NA) cases

Si puedo hacer que este bucle se ejecute, me gustaría agregar los coeficientes e incluso mejor los valores r-squared para cada modelo a un data.frame vacío. Cualquier ayuda sería muy apreciada.

Responder

Como referencia, aquí está la respuesta de jlhoward (que incorpora el comentario de BondedDust) haciendo uso de las * funciones de aplicación que se encuentran en esta excelente pregunta: R Funciones de agrupamiento: sapply vs. lapply vs. apply. vs. tapply vs. por vs. agregado

models <- sapply(unique(as.character(df$country)),
                 function(cntry)lm(BirthRate~US.,df,subset=(country==cntry)),
                 simplify=FALSE,USE.NAMES=TRUE)

# to summarize all the models
lapply(models,summary)
# to run anova on all the models
lapply(models,anova)

#This produces a named list of models, so you could extract the model for Aruba as:
models[["Aruba"]]

Answer 1

Aquí hay una forma de usar dplyr . Creé Brazil usando tus datos. Entonces tienes resultados idénticos para los dos países. Verás intersección y pendiente en el marco de datos final.

CODIGO DE DATOS

structure(list(country = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
2L, 2L, 2L, 2L, 2L, 2L), .Label = c("Aruba", "Brazil"), class = "factor"), 
date = c(2011L, 2010L, 2009L, 2008L, 2007L, 2006L, 2011L, 
2010L, 2009L, 2008L, 2007L, 2006L), BirthRate = c(10.584, 
10.804, 11.06, 11.346, 11.653, 11.977, 10.584, 10.804, 11.06, 
11.346, 11.653, 11.977), US. = c(25354.8, 24289.1, 24639.9, 
27549.3, 25921.3, 24015.4, 25354.8, 24289.1, 24639.9, 27549.3, 
25921.3, 24015.4)), .Names = c("country", "date", "BirthRate", 
"US."), class = "data.frame", row.names = c("4", "5", "6", "7", 
"8", "9", "10", "11", "12", "13", "14", "15"))

#   country date BirthRate     US.
#4    Aruba 2011    10.584 25354.8
#5    Aruba 2010    10.804 24289.1
#6    Aruba 2009    11.060 24639.9
#7    Aruba 2008    11.346 27549.3
#8    Aruba 2007    11.653 25921.3
#9    Aruba 2006    11.977 24015.4
#10  Brazil 2011    10.584 25354.8
#11  Brazil 2010    10.804 24289.1
#12  Brazil 2009    11.060 24639.9
#13  Brazil 2008    11.346 27549.3
#14  Brazil 2007    11.653 25921.3
#15  Brazil 2006    11.977 24015.4

group_by(mydf, country) %>%
    do({model = lm(BirthRate ~ US., data = .);
       data.frame(int = coef(model)[1], slope = coef(model)[2])})

#  country      int        slope
#1   Aruba 11.02503 8.393295e-06
#2  Brazil 11.02503 8.393295e-06

Answer 2

Devolver una lista de dplyr aún no es posible. Si solo necesitas el intercepto y la pendiente, la respuesta de @jazzurro es la opción, pero si necesitas todo el modelo, necesitas hacer algo como

library(dplyr)
models <- df %>% group_by(country) %>% do(mod = lm(BirthRate ~ US., data = .))

Entonces, si desea realizar ANOVA en cada modelo ajustado, puede hacerlo utilizando rowwise

models %>% rowwise %>% do(anova(.$mod))

pero una vez más, el resultado es forzado a un marco de datos y no es lo mismo que hacer lapply(models$mod, anova) .

Por ahora (es decir, hasta la próxima versión de dplyr ) si necesita almacenar el resultado completo en una lista, puede usar dlply de plyr , como plyr::dlply(df, "country", function(d) anova(lm(BirthRate ~ US., data = d))) , o por supuesto si no tienes que usar dplyr , puedes ir a la respuesta de @SvenHohenstein que de todos modos parece una mejor manera de hacerlo.

Answer 3

Eche un vistazo a la función nlme paquete nlme :

library(nlme)
lmList(BirthRate ~ US. | country, df)

Aquí, | country | country se usa para crear una regresión para cada país en particular.