[R-sig-ME] Predict for glmer

[Thierry Onkelinx]

Dear Marc,

You probably want random <- er$R1[d$R1, ] + er$R2[d$R2, ]

Best regards,

ir. Thierry Onkelinx
Statisticus / Statistician

Vlaamse Overheid / Government of Flanders
INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE AND
FOREST
Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
thierry.onkelinx using inbo.be
Havenlaan 88 bus 73, 1000 Brussel
www.inbo.be

///////////////////////////////////////////////////////////////////////////////////////////
To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to say
what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey
///////////////////////////////////////////////////////////////////////////////////////////

<https://www.inbo.be>


Op zo 27 sep. 2020 om 10:53 schreef Marc Girondot via R-sig-mixed-models <
r-sig-mixed-models using r-project.org>:

> 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
>
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>

	[[alternative HTML version deleted]]