[R-sig-ME] testing random effect

Sundar Dorai-Raj sundar.dorai-raj at pdf.com
Mon Feb 18 17:31:09 CET 2008



Douglas Bates said the following on 2/16/2008 7:08 AM:
> 2008/2/16 Maria Silvina Fenoglio <msfeno at yahoo.com.ar>:
>> Hi all,
>>   I am trying to test the significance of the random effect in the following model but I didn't find the way to do that.
>> I think that I can do it with LRT, removing the random effect and comparing both models, but an error appears...
>>   Can anybody help me with that?
>>
>>   mod<-lmer(Riqsp~Distpza+logtam+logaisl+Densmin+(1|Sitio),sil,family=poisson,method="ML")
> 
> Could you give us a bit more detail, please?  How did you try to
> perform the likelihood ratio test?
> 
> If possible, could you make your data available so we can check the
> fits against different versions of the lme4 package?
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


Hi, Prof. Bates.

Here's an example that may help. "sessionInfo" at bottom.

library(lme4)
set.seed(1)
b <- rnorm(10, sd = 0.5)
z <- factor(rep(1:10, each = 10))
eta <- 1 + b[z]
y <- rpois(100, exp(eta))

## fit glmm model
(fit <- lmer(y ~ 1 + (1 | z), family = poisson, method = "ML"))
## fit null model without random effect
(fit0 <- glm(y ~ 1, family = poisson))

## this fails
anova(fit, fit0)
# Error in FUN(X[[1L]], ...) :
#   no slot of name "call" for this object of class "glm"
# In addition: Warning message:
# In logLik.glm(X[[2L]], ...) : extra arguments discarded

## compare log-likelihoods
(ll1.a <- logLik(fit))
## 'log Lik.' -65.58766 (df=2)
## the latter appears to be missing a constant
## making it impossible to compare to the following
(ll0.a <- logLik(fit0))
## 'log Lik.' -1164.184 (df=1)
## compute manually
(ll1.b <- sum(dpois(y, exp(fitted(fit)), log = TRUE)))
## [1] -187.6246
(ll0.b <- sum(dpois(y, fitted(fit0), log = TRUE)))
## [1] -1164.184

(LR <- -2 * (ll0.b - ll1.b))
## [1] 1953.119

## sorry haven't upgraded to 2.6.2 yet
 > sessionInfo()
R version 2.6.1 (2007-11-26)
i386-pc-mingw32

locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United 
States.1252;LC_MONETARY=English_United 
States.1252;LC_NUMERIC=C;LC_TIME=English_United States.1252

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] lme4_0.99875-9    Matrix_0.999375-4 lattice_0.17-4

loaded via a namespace (and not attached):
[1] grid_2.6.1  tools_2.6.1
 >




More information about the R-sig-mixed-models mailing list