usando predecir con una lista de objetos lm()

Question

Tengo datos en los que regularmente ejecuto regresiones. Cada "porción" de datos se ajusta a una regresión diferente. Cada estado, por ejemplo, podría tener una función diferente que explique el valor dependiente. Esto parece ser un típico problema de "combinación de aplicación dividida", así que estoy usando el paquete plyr. Puedo crear fácilmente una lista de objetos lm() que funciona bien. Sin embargo, no puedo entender cómo uso esos objetos más adelante para predecir valores en un data.frame separado.

Aquí está un ejemplo totalmente artificial que ilustra lo que estoy tratando de hacer:

# setting up some fake data 
set.seed(1) 
funct <- function(myState, myYear){ 
    rnorm(1, 100, 500) + myState + (100 * myYear) 
} 
state <- 50:60 
year <- 10:40 
myData <- expand.grid(year, state) 
names(myData) <- c("year","state") 
myData$value <- apply(myData, 1, function(x) funct(x[2], x[1])) 
## ok, done with the fake data generation. 

require(plyr) 

modelList <- dlply(myData, "state", function(x) lm(value ~ year, data=x)) 
## if you want to see the summaries of the lm() do this: 
    # lapply(modelList, summary) 

state <- 50:60 
year <- 50:60 
newData <- expand.grid(year, state) 
names(newData) <- c("year","state") 
## now how do I predict the values for newData$value 
    # using the regressions in modelList? 

Entonces, ¿cómo puedo usar los objetos contenidos en lm()modelList para predecir los valores utilizando los valores independientes año y el Estado de newData?

Answer 1

Aquí está mi intento:

predNaughty <- ddply(newData, "state", transform, 
    value=predict(modelList[[paste(piece$state[1])]], newdata=piece)) 
head(predNaughty) 
# year state value 
# 1 50 50 5176.326 
# 2 51 50 5274.907 
# 3 52 50 5373.487 
# 4 53 50 5472.068 
# 5 54 50 5570.649 
# 6 55 50 5669.229 
predDiggsApproved <- ddply(newData, "state", function(x) 
    transform(x, value=predict(modelList[[paste(x$state[1])]], newdata=x))) 
head(predDiggsApproved) 
# year state value 
# 1 50 50 5176.326 
# 2 51 50 5274.907 
# 3 52 50 5373.487 
# 4 53 50 5472.068 
# 5 54 50 5570.649 
# 6 55 50 5669.229 

JD largo de edición

Me inspiré suficiente para elaborar una opción adply():

pred3 <- adply(newData, 1, function(x) 
    predict(modelList[[paste(x$state)]], newdata=x)) 
head(pred3) 
# year state  1 
# 1 50 50 5176.326 
# 2 51 50 5274.907 
# 3 52 50 5373.487 
# 4 53 50 5472.068 
# 5 54 50 5570.649 
# 6 55 50 5669.229 

Answer 2

¿Qué hay de malo en

lapply(modelList, predict, newData) 

?

EDIT:

Gracias para explicar lo que está mal con eso. ¿Qué tal:

newData <- data.frame(year) 
ldply(modelList, function(model) { 
    data.frame(newData, predict=predict(model, newData)) 
}) 

iterar sobre los modelos, y aplicar los nuevos datos (que es el mismo para cada estado, ya que acaba de hacer un expand.grid para crearla).

EDIT 2:

Si newData no tiene los mismos valores para year para cada state como en el ejemplo, un enfoque más general se puede utilizar. Tenga en cuenta que esto utiliza la definición original de newData, no la de la primera edición.

ldply(state, function(s) { 
    nd <- newData[newData$state==s,] 
    data.frame(nd, predict=predict(modelList[[as.character(s)]], nd)) 
}) 

Primeras 15 líneas de esta salida:

year state predict 
1 50 50 5176.326 
2 51 50 5274.907 
3 52 50 5373.487 
4 53 50 5472.068 
5 54 50 5570.649 
6 55 50 5669.229 
7 56 50 5767.810 
8 57 50 5866.390 
9 58 50 5964.971 
10 59 50 6063.551 
11 60 50 6162.132 
12 50 51 5514.825 
13 51 51 5626.160 
14 52 51 5737.496 
15 53 51 5848.832 

Answer 3

I toman la parte dura está emparejando cada estado en newData para el modelo correspondiente.

¿Algo así como esto?

predList <- dlply(newData, "state", function(x) { 
    predict(modelList[[as.character(min(x$state))]], x) 
}) 

Aquí he utilizado una forma "hacky" de extraer el correspondiente modelo de estado: as.character(min(x$state))

... Probablemente hay una manera mejor?

Salida:

> predList[1:2] 
$`50` 
     1  2  3  4  5  6  7  8  9  10  11 
5176.326 5274.907 5373.487 5472.068 5570.649 5669.229 5767.810 5866.390 5964.971 6063.551 6162.132 

$`51` 
     12  13  14  15  16  17  18  19  20  21  22 
5514.825 5626.160 5737.496 5848.832 5960.167 6071.503 6182.838 6294.174 6405.510 6516.845 6628.181 

O, si quieres un data.frame como salida:

predData <- ddply(newData, "state", function(x) { 
    y <-predict(modelList[[as.character(min(x$state))]], x) 
    data.frame(id=names(y), value=c(y)) 
}) 

Salida:

head(predData) 
    state id value 
1 50 1 5176.326 
2 50 2 5274.907 
3 50 3 5373.487 
4 50 4 5472.068 
5 50 5 5570.649 
6 50 6 5669.229 

Answer 4

una solución con sólo base R. El formato de la salida es diferente, pero todos los valores son justo ahí.

models <- lapply(split(myData, myData$state), 'lm', formula = value ~ year) 
pred4 <- mapply('predict', models, split(newData, newData$state)) 

Answer 5

Es necesario utilizar mdply para suministrar tanto el modelo como los datos a cada función de llamada:

dataList <- dlply(newData, "state") 

preds <- mdply(cbind(mod = modelList, df = dataList), function(mod, df) { 
    mutate(df, pred = predict(mod, newdata = df)) 
}) 

Answer 6

Tal vez me falta algo, pero creo que lmList es la herramienta ideal aquí,

library(nlme) 
ll = lmList(value ~ year | state, data=myData) 
predict(ll, newData) 


## Or, to show that it produces the same results as the other proposed methods... 
newData[["value"]] <- predict(ll, newData) 
head(newData) 
# year state value 
# 1 50 50 5176.326 
# 2 51 50 5274.907 
# 3 52 50 5373.487 
# 4 53 50 5472.068 
# 5 54 50 5570.649 
# 6 55 50 5669.229