[R-sig-ME] Predict for glmer

[Marc Girondot]

Hi,

I try to understand exactly how mixed model is working. Let do this 
simple model with fake data:


library(lme4)

invlogit <- function(n) {1/(1 + exp(-n))}

d <- data.frame(A=c(10, 20, 10, 20, 30, 15, 17, 19, 20, 30, 20, 30, 21, 
23),
                 B=c(2, 3, 7, 9, 10, 8, 5, 7, 9, 10, 2, 3, 7, 8),
                 D=c(10, 11, 12, 13, 14, 10, 11, 12, 13, 14, 13, 12, 12, 
13),
                 R1=c("A", "A", "A", "A", "A", "A", "A", "B", "B", "B", 
"B", "B", "B", "B"),
                 R2=c("C", "C", "C", "D", "D", "D", "D", "E", "E", "E", 
"F", "F", "F", "F"))

m <- glmer(formula = cbind(A, B) ~ D + (1 | R1 / R2),
       data=d,
       family = binomial(link = "logit"))

It seems to work well (except the warning for boundary because probably 
I don't introduce enough data).

Let do a predict:

predict(m, type="response")

         1         2         3         4         5 6         7         
8         9        10        11
0.7706127 0.7665437 0.7624248 0.7502806 0.7459699 0.7629174 0.7587547 
0.7572462 0.7530161 0.7487368 0.7696229
        12        13        14
0.7736541 0.7736541 0.7696229

It is ok. Now I try to do the predict by hand to be sure that I 
understand how it works:

ef <- fixef(m)
er <- ranef(m)

fixed <- model.matrix(~ D, data = d) %*% ef
random <- rowSums(model.matrix(~ R1, data = d) * er$R2[d$R2, ])
invlogit(fixed[, 1] + random)

         1         2         3         4         5         6 7         
8         9        10        11

0.7706127 0.7665437 0.7624248 0.7502806 0.7459699 0.7629174 0.7587547 
0.7523144 0.7480270 0.7436906 0.7809063
        12        13        14
0.7847953 0.7847953 0.7809063

The first 7 estimates are ok but not the last 7. So it is related to R1 
factor... but I don't understand why it is badly estimated.


If someone can help me...

Thanks a lot

Marc