From: Spencer Graves <spencer.graves_at_pdf.com>

Date: Mon 26 Sep 2005 - 01:03:32 EST

Date: Mon 26 Sep 2005 - 01:03:32 EST

Sorry I couldn't help more. spencer graves

Robert Bagchi wrote:

> Dear R users,

*>
**> I have been having problems getting believable estimates from anova on a
**> model fit from lmer. I get the impression that F is being greatly
**> underestimated, as can be seen by running the example I have given below.
**>
**> First an explanation of what I'm trying to do. I am trying to fit a glmm
**> with binomial errors to some data. The experiment involves 10
**> shadehouses, divided between 2 light treatments (high, low). Within each
**> shadehouse there are 12 seedlings of each of 2 species (hn & sl). 3
**> damage treatments (0, 0.1, 0.25 leaf area removal) were applied to the
**> seedlings (at random) so that there are 4 seedlings of each
**> species*damage treatment in each shadehouse. There maybe a shadehouse
**> effect, so I need to include it as a random effect. Light is applied to
**> a shadehouse, so it is outer to shadehouse. The other 2 factors are
**> inner to shadehouse.
**>
**> We want to assess if light, damage and species affect survival of
**> seedlings. To test this I fitted a binomial mixed effects model with
**> lmer (actually with quasibinomial errors). THe summary function suggests
**> a large effect of both light and species (which agrees with graphical
**> analysis). However, anova produces F values close to 0 and p values
**> close to 1 (see example below).
**>
**> Is this a bug, or am I doing something fundamentally wrong? If anova
**> doesn't work with lmer is there a way to perform hypothesis tests on
**> fixed effects in an lmer model? I was going to just delete terms and
**> then do liklihood ratio tests, but according to Pinheiro & Bates (p. 87)
**> that's very untrustworthy. Any suggestions?
**>
**> I'm using R 2.1.1 on windows XP and lme4 0.98-1
**>
**> Any help will be much appreciated.
**>
**> many thanks
**> Robert
**>
**> ###############################
**> The data are somewhat like this
**>
**> #setting up the dataframe
**>
**> bm.surv<-data.frame(
**> house=rep(1:10, each=6),
**> light=rep(c("h", "l"), each=6, 5),
**> species=rep(c("sl", "hn"), each=3, 10),
**> damage=rep(c(0,.1,.25), 20)
**> )
**>
**> bm.surv$survival<-ifelse(bm.surv$light=="h", rbinom(60, 4, .9),
**> rbinom(60, 4, .6)) # difference in probablility should ensure a
**> light effect
**> bm.surv$death<-4-bm.surv$survival
**>
**> # fitting the model
**> m1<-lmer(cbind(survival, death)~light+species+damage+(1|house),
**> data=bm.surv, family="quasibinomial")
**>
**> summary(m1) # suggests that light is very significant
**> Generalized linear mixed model fit using PQL
**> Formula: cbind(survival, death) ~ light + species + damage + (1 | table)
**> Data: bm.surv
**> Family: quasibinomial(logit link)
**> AIC BIC logLik deviance
**> 227.0558 239.6218 -107.5279 215.0558
**> Random effects:
**> Groups Name Variance Std.Dev.
**> table (Intercept) 1.8158e-09 4.2613e-05
**> Residual 3.6317e+00 1.9057e+00
**> # of obs: 60, groups: table, 10
**>
**> Fixed effects:
**> Estimate Std. Error DF t value Pr(>|t|)
**> (Intercept) 2.35140 0.36832 56 6.3841 3.581e-08 ***
**> lightl -1.71517 0.33281 56 -5.1535 3.447e-06 ***
**> speciessl -0.57418 0.30085 56 -1.9085 0.06145 .
**> damage 1.49963 1.46596 56 1.0230 0.31072
**> ---
**> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
**>
**> Correlation of Fixed Effects:
**> (Intr) lightl spcssl
**> lightl -0.665
**> speciessl -0.494 0.070
**> damage -0.407 -0.038 -0.017
**>
**>
**> anova(m1) # very low F value for light, corresponding to p
**> values approaching 1
**>
**> Analysis of Variance Table
**> Df Sum Sq Mean Sq Denom F value Pr(>F)
**> light 1 0.014 0.014 56.000 0.0018 0.9661
**> species 1 0.002 0.002 56.000 0.0002 0.9887
**> damage 1 0.011 0.011 56.000 0.0014 0.9704
**>
**>
*

-- Spencer Graves, PhD Senior Development Engineer PDF Solutions, Inc. 333 West San Carlos Street Suite 700 San Jose, CA 95110, USA spencer.graves@pdf.com www.pdf.com <http://www.pdf.com> Tel: 408-938-4420 Fax: 408-280-7915 ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.htmlReceived on Mon Sep 26 01:06:14 2005

*
This archive was generated by hypermail 2.1.8
: Fri 03 Mar 2006 - 03:40:26 EST
*