Plotting Regression results from lme4 in R using Lattice (or something else)

Question

I have fit a regression using lme4 thanks to a previous answer. Now that I have a regression fit for each state I'd like to use lattice to plot QQ plots for each state. I would also like to plot error plots for each state in a lattice format. How do I make a lattice plot using the results of a lme4 regression?

Below is a simple sample (yeah, I like a good alliteration) using two states. I would like to make a two panel lattice made from the object fits.

library(lme4)
d <- data.frame(state=rep(c('NY', 'CA'), c(10, 10)), year=rep(1:10, 2), response=c(rnorm(10), rnorm(10)))
fits <- lmList(response ~ year | state, data=d)

Answer 1

Instead of using lmList, I'd recommend the more general plyr package.

library(plyr)

d <- data.frame(
 state = rep(c('NY', 'CA'), c(10, 10)), 
 year = rep(1:10, 2), 
 response = c(rnorm(10), rnorm(10))
)

# Create a list of models
# dlply = data frame -> list
models <- dlply(d, ~ state, function(df) { 
  lm(response ~ year, data = df)
})

# Extract the coefficients in a useful form
# ldply = list -> data frame
ldply(models, coef)

# We can get the predictions in a similar way, but we need
# to cast to a data frame so the numbers come out as rows,
# not columns.
predictions <- ldply(models, as.data.frame(predict))

predictions is a regular R data frame and so is easy to plot.

Answer 2

I am not sure you can get this into lattice easily. What you have in fits is a an S4 object containing a .Data slot with a list of standard lm objects:

R> class(fits)
[1] "lmList"
attr(,"package")
[1] "lme4"
R> class(fits@.Data)
[1] "list"
R> class(fits@.Data[[1]])
[1] "lm"
R> op <- par(mfrow=c(2,4))
R> invisible(lapply(fits@.Data, plot))

This last plot simply plots you the standard 2x2 plot for lm objects twice, once for each element of the list of fitted objects. Use the which argument to plot to select subsets of these or for other regression diagnostics.

If you actually want lattice plots of predicted vs actual, you may have to program this.

Answer 3

I've had some trouble with lme4::lmList. For example, summary doesn't seem to work. So you might run into some problems because of that.

So even though I use lmer, instead of lme, I've been explicitly calling nlme::lmList instead. Then summary etc. will work.