[R] model comparison with mixed effects glm

[Spencer Graves]

	  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)
>>> > [1] 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)
>>
>> Hadley
> 
>