[R] lme: extract variance estimate
Stephen Ellner
spe2 at cornell.edu
Tue Jul 6 22:26:52 CEST 2004
For a Monte Carlo study I need to extract from an lme model
the estimated standard deviation of a random effect
and store it in a vector. If I do a print() or summary()
on the model, the number I need is displayed in the Console
[it's the 0.1590195 in the output below]
>print(fit)
>Linear mixed-effects model fit by maximum likelihood
> Data: datag2
> Log-likelihood: -145.0028
> Fixed: lst1 ~ lst
>(Intercept) lst
> 1.1080629 0.7582595
>
>Random effects:
> Formula: ~1 | yeart
> (Intercept) Residual
>StdDev: 0.1590195 0.3313497
However I have not been able to find that number anywhere in the
model object returned by lme. Can anyone tell me how to pull it
or compute it from a model returned by lme? Ditto for models with
random intercept fitted by glmmPQL with family=binomial.
Details:
The relevant data are the first 3 colums from the
data frame extracted below. Yeart is a factor (year of
study) ranging from 1 to 16. The model is linear regression of
lst1 on lst with constant slope, random intercept varying
across years. I'm using R 1.9.1 in WinXP and the current
version of nlme from CRAN.
Code:
datag2=groupedData(lst1~1|yeart,data=cdata);
mixed.grow=lme(fixed=lst1~lst,random=~1, method="ML",data=datag2);
mixed.surv=glmmPQL(fixed=surv~lst, random=~1|yeart, family=binomial,
data=cdata,niter=50);
Some of the data:
yeart lst lst1 surv flow
1 1 3.158088 3.331216 1 0
2 1 2.472618 2.486410 1 0
3 1 3.582950 3.417807 1 0
4 1 3.554819 3.377117 1 0
5 1 2.830049 2.992426 1 0
6 1 2.779616 3.240449 1 0
7 1 2.580484 2.930154 1 0
Thanks in advance,
Steve
Stephen P. Ellner (spe2 at cornell.edu)
Department of Ecology and Evolutionary Biology
Corson Hall, Cornell University, Ithaca NY 14853-2701
Phone (607) 254-4221 FAX (607) 255-8088
More information about the R-help
mailing list