[R] model comparison with mixed effects glm
spencer.graves at pdf.com
Wed Apr 5 17:08:57 CEST 2006
Another thought on checking the validity of the suggested
2*log(likelihood ratio) procedure I suggested: If it were my problem, I
think I would do some checking using Monte Carlo, e.g., as described in
sec. 2.6 of the vignette "MlmSoftRev" in the "mlmRev" package. This is
particularly relevant for testing a parameter at a boundary, e.g.,
whether a particular variance component is 0, because the assumptions
for the traditional chi-square approximation to 2*log(LR) do not hold in
that case, as documented in sec. 2.4 of Pinheiro and Bates (2000)
Mixed-Effects Models in S and S-Plus (Springer).
Spencer Graves wrote:
> You are correct on both counts. The exta line is inserted below;
> obviously, I had it but failed to copy it into the email.
> And you are also correct that one needs to be careful that both
> glm and lmer are using comparable definitions for the log(likelihood).
> My crude check on that was just to look compare the lglk0 and lglk.ID1.;
> the numbers seemed too close to be based on different definitions. In
> addition, I think I may have checked this once before, but my memory
> could be faulty on that point.
> Thanks for pointing out both deficiencies in my reply.
> spencer graves
> hadley wickham wrote:
>>> ### To get around that, I computed 2*log(likelihood ratio) manually:
>>> lglk0 <- logLik(fit0)
>>> lglk.ID1. <- logLik(Fit.ID1.)
> chisq.ID. <- 2*(lglk.ID1.-lglk0)
>>> pchisq(as.numeric(chisq.ID.), 1, lower=FALSE)
>>> >  0.008545848
>> (I think you're missing a line in there)
>> But isn't this rather perilous unless you are confident that the two
>> models are using exactly the same formulation of the likelihood? (ie.
>> that they are truly nested)
More information about the R-help