[R-sig-ME] Loglik of model with maxfun = 0 is not the one expected
Marc Girondot
m@rc_grt @end|ng |rom y@hoo@|r
Mon Sep 28 22:25:18 CEST 2020
Thanks a lot to Thierry Onkelinx <thierry.onkelinx using inbo.be> to answer to
my previous demand. It works now !
Another point that I don't understand.
I fit a model. I get -log likelihood
I extract the fitted parameters.
I inject these fitted parameters to estimate the likelihood of the model
and I obtain a different -log likelihood.
So I suppose that I don't extract all the fitted parameters... but which
one ?
Thanks if you can help me to understand what's happened here.
Marc
___________________
library(lme4)
# Ei is simply data that will be used to fit the model
Ei <- structure(list(Ma = c(15, 14, 28, 5, 26, 7, 25, 4, 9, 25,
8, 11, 2, 3, 1, 1, 0, 0, 10, 5, 0, 0, 9, 1, 2,
0, 10, 0, 0, 0,
8, 1, 0, 0, 12, 3, 29, 10, 26, 9, 18, 5, 14,
6, 1, 2, 0, 0),
Fa = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 5, 0, 4, 5, 18, 2,
11, 2, 10, 0, 3, 7, 3, 0, 9, 8, 8, 0, 10, 9,
9, 0, 8, 5,
4, 0, 0, 0, 0, 5, 1, 2, 0, 8, 0, 8, 3, 10, 3),
IT = c(26.4,
26.9, 27.9, 27.9, 28.4, 28.4,
29, 29, 29, 29.5, 29.5, 29.5,
30, 30, 30.3, 30.3, 30.8,
30.8, 28, 29.5, 31, 32.5, 28, 29.5,
31, 32.5, 28, 29.5, 31, 32.5,
28, 29.5, 31, 32.5, 28.15,
28.15, 28.65, 28.65, 29.15,
29.15, 29.55, 29.55, 29.75, 29.75,
30.05, 30.05, 30.65, 30.65),
RMU = structure(c(2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L,
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
3L, 3L, 3L, 3L, 3L, 3L,
3L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L),
.Label = c("ASW", "AW",
"PSW"), class = "factor"),
C = structure(c(4L,
4L, 4L, 3L, 4L, 3L, 4L, 7L, 8L, 4L,
7L, 8L, 7L, 8L, 7L, 8L,
7L, 8L, 5L, 5L, 5L, 5L, 6L, 6L, 6L,
6L, 10L, 10L, 10L, 10L,
9L, 9L, 9L, 9L, 1L, 2L, 1L, 2L, 1L,
2L, 1L, 2L, 1L, 2L, 1L,
2L, 1L, 2L),
.Label = c("A", "B", "C", "D",
"E", "F", "G", "H", "I",
"J"), class = "factor")),
class = "data.frame", row.names = 302:349)
# m1 if the fitted model
m1 <- glmer(formula = cbind(Ma, Fa) ~ 0 + IT + (1 | RMU / C) + (0 + IT |
RMU),
data=Ei,
control = glmerControl(optimizer="bobyqa", optCtrl =
list(maxfun = 1000000)),
family = binomial(link = "logit"))
# I check that all is ok.
logLik(m1)
ranef(m1)
fixef(m1)
theta <- getME(m1, "theta")
beta <- getME(m1, "beta")
# use the fitted values of m1 that I introduce in m2. Note that the
model is not fitted because maxfun = 0. Log likelihood is just calculated
m2 <- glmer(formula = cbind(Ma, Fa) ~ 0 + IT + (1 | RMU / C) + (0 + IT |
RMU),
data=Ei,
control = glmerControl(optimizer="bobyqa", optCtrl =
list(maxfun = 0)),
start = list(theta = theta,
fixef = beta),
family = binomial(link = "logit"))
# Both likelihood should be identical... they are not
logLik(m1)
logLik(m2)
More information about the R-sig-mixed-models
mailing list