[R-sig-ME] No estimated scale value given
bolker at ufl.edu
Tue Feb 17 21:23:28 CET 2009
The not-necessarily-obvious thing is that none of the
parameter estimates change when you go from poisson,
binomial, etc. to their quasi- equivalents -- the variance
on all points is inflated by the same amount, so the
point estimates don't change. So you can go ahead and
use the quasi- variant to get the estimated scale.
I'm still struggling with exactly what "sigma" is in
the quasi- case (as is Doug Bates) -- it seems
NOT equal to the Pearson residuals^2/(resid df) ?
It would be worth exploring this some more ...
ntot <- 1000
x <- runif(ntot)
f <- rep(1:20,each=50)
reff <- rnorm(20,mean=0,sd=0.5)
y <- rnbinom(ntot,mu=exp(2*x-1+reff[f]),size=1)
m1 <- glmer(y~x+(1|f),family=poisson)
m2 <- update(m1,family=quasipoisson)
lme4:::sigma(m1) ## 1
lme4:::sigma(m2) ## 1.04 (??)
## residuals returned are Pearson residuals
## (uncorrected for overdispersion)
myres <- (y-fitted(m1))/sqrt(fitted(m1))
## not quite identical, but nearly
all.equal(residuals(m1),residuals(m2)) ## TRUE
sum(residuals(m1)^2)/ntot ## NOT sigmaML
m1 at deviance["sigmaML"]
m2 at deviance["sigmaML"]
Amy Wade wrote:
> When I fit a quasi- variant it gives what seems to be a sensible value for
> the sigma. The trouble is I would like to know if my models are
> overdispersed without the quasi-variant. I'm using data with poisson
> distribution from which I know what the estimated scale is and it still
> gives '1' when it should be 0.97. There was a discussion about this
> previously and it was suggested to use 'lme4:::sigma(model)', that just
> gives me the same result as using 'summary(model)@sigma'.
> Thanks very much for your help.
> Should I have put this on the public thread? Not really sure how to do that!
> -----Original Message-----
> From: Ben Bolker [mailto:bolker at ufl.edu]
> Sent: 16 February 2009 14:11
> To: Amy Wade
> Subject: Re: [R-sig-ME] No estimated scale value given
> What happens if you fit the same model with a quasi- variant and then try
> to access sigma?
> Ben Bolker
> Amy Wade wrote:
>> I'm trying to run generalized linear mixed models using the lmer()
>> function in the lme4 package. The problem is that the output does not
>> give me a value for the Estimated scale. The rest of the output is as
>> it should be. This makes it very difficult to assess whether the model
>> is overdispersed. My colleague tried the same code on his computer and
>> it did give an estimated scale value. I tried unistalling R and
>> reinstalling the latest edition
>> (2.8.1) with the latest lme4, but this made no difference. I also
>> tried 'summary(model)@sigma' which people suggested on the CRAN forums
>> to extract the scale parameter. This always gave me an answer of '1'
>> even when I used data where I knew this was not the case. When I load
>> up the lme4 library it does warn me about several objects being masked
> from 'stats' and 'base'.
>> Any idea why the output is not giving me the estimate scale? Or any
>> idea how to extract it somehow?
>> Many thanks,
>> [[alternative HTML version deleted]]
>> R-sig-mixed-models at r-project.org mailing list
> Ben Bolker
> Associate professor, Biology Dep't, Univ. of Florida bolker at ufl.edu /
> www.zoology.ufl.edu/bolker GPG key:
Associate professor, Biology Dep't, Univ. of Florida
bolker at ufl.edu / www.zoology.ufl.edu/bolker
GPG key: www.zoology.ufl.edu/bolker/benbolker-publickey.asc
More information about the R-sig-mixed-models