[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