Re: [R] Components of variance

From: Renaud Lancelot <renaud.lancelot_at_cirad.fr>
Date: Sun 26 Jun 2005 - 16:23:16 EST

Spencer Graves a écrit :
> As far as I know, the best available is lme in library(nlme). For
> more information, see the the following:
>
> Pinheiro and Bates (2000) Mixed-Effects Models in S and S-Plus
> (Springer)
>
> Consider the following example:
>
> > set.seed(2)
> > lot <- rep(LETTERS[1:3], each=9)
> > lot.e <- rep(rnorm(3), each=9)
> > wf <- paste(lot, rep(1:9, each=3))
> > wf.e <- rep(rnorm(9), each=3)
> > DF <- data.frame(lot=lot, wafer=wf,

> + Vt=lot.e+wf.e+rnorm(27))
> > (fit <- lme(Vt~1, random=~1|lot/wafer, DF))
> Linear mixed-effects model fit by REML
> Data: DF
> Log-restricted-likelihood: -48.44022
> Fixed: Vt ~ 1
> (Intercept)
> 0.6083933
>
> Random effects:
> Formula: ~1 | lot
> (Intercept)
> StdDev: 1.230572
>
> Formula: ~1 | wafer %in% lot
> (Intercept) Residual
> StdDev: 0.9801403 1.161218
>
> Number of Observations: 27
> Number of Groups:
> lot wafer %in% lot
> 3 9
> > (CI.fit <- intervals(fit))
> Approximate 95% confidence intervals
>
> Fixed effects:
> lower est. upper
> (Intercept) -1.100281 0.6083933 2.317068
> attr(,"label")
> [1] "Fixed effects:"
>
> Random Effects:
> Level: lot
> lower est. upper
> sd((Intercept)) 0.3368174 1.230572 4.495931
> Level: wafer
> lower est. upper
> sd((Intercept)) 0.426171 0.9801403 2.254201
>
> Within-group standard error:
> lower est. upper
> 0.8378296 1.1612179 1.6094289
> > str(CI.fit)
> List of 3
> $ fixed : num [1, 1:3] -1.100 0.608 2.317
> ..- attr(*, "dimnames")=List of 2
> .. ..$ : chr "(Intercept)"
> .. ..$ : chr [1:3] "lower" "est." "upper"
> ..- attr(*, "label")= chr "Fixed effects:"
> $ reStruct:List of 2
> ..$ lot :`data.frame': 1 obs. of 3 variables:
> .. ..$ lower: num 0.337
> .. ..$ est. : num 1.23
> .. ..$ upper: num 4.5
> ..$ wafer:`data.frame': 1 obs. of 3 variables:
> .. ..$ lower: num 0.426
> .. ..$ est. : num 0.98
> .. ..$ upper: num 2.25
> ..- attr(*, "label")= chr "Random Effects:"
> $ sigma : atomic [1:3] 0.838 1.161 1.609
> ..- attr(*, "label")= chr "Within-group standard error:"
> - attr(*, "level")= num 0.95
> - attr(*, "class")= chr "intervals.lme"
> > diff(log(CI.fit$sigma))
> est. upper
> 0.32641 0.32641
>
> The last line combined with help for intervals.lme shows that the
> confidence interval for sigma (and doubtless also for lot and wafer
> variance components is based on a normal approximation for the
> distribution of log(sigma).
>
> The state of the art is reflected in "lmer" in library(lme4),
> described in the following:
>
> Doug Bates (2005) "Fitting linear mixed models in R" in R News 5/1
> available from "www.r-project.org" -> Newsletter
>
> However, an "intervals" function is not yet available for "lmer"
> objects.

The issue of variance components in lme4 was recently by D. Bates:

http://tolstoy.newcastle.edu.au/R/help/05/06/7291.html

Also, a colleague wrote (in French) a nice report - with data and R code on how to compute variance components, their bias and distribution with lme (normal approximation , Monte Carlo simulation and delta method). Let me know if you want it.

Best,

Renaud

>
> spencer graves
>
> John Sorkin wrote:
>
>

>>Could someone identify a function that I might use to perform a
>>components of variance analysis? In addition to the variance
>>attributable to each factor, I would also like to obtain the SE of the
>>variances.
>>Thank you,
>>John
>> 
>>John Sorkin M.D., Ph.D.
>>Chief, Biostatistics and Informatics
>>Baltimore VA Medical Center GRECC and
>>University of Maryland School of Medicine Claude Pepper OAIC
>> 
>>University of Maryland School of Medicine
>>Division of Gerontology
>>Baltimore VA Medical Center
>>10 North Greene Street
>>GRECC (BT/18/GR)
>>Baltimore, MD 21201-1524
>> 
>>410-605-7119 
>>--- NOTE NEW EMAIL ADDRESS:
>>jsorkin@grecc.umaryland.edu
>>
>>	[[alternative HTML version deleted]]
>>
>>______________________________________________
>>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.html

>
>
-- 
Dr Renaud Lancelot, vétérinaire
Projet FSP régional épidémiologie vétérinaire
C/0 Ambassade de France - SCAC
BP 834 Antananarivo 101 - Madagascar

e-mail: renaud.lancelot@cirad.fr
tel.:   +261 32 40 165 53 (cell)
         +261 20 22 665 36 ext. 225 (work)
         +261 20 22 494 37 (home)

______________________________________________
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.html
Received on Sun Jun 26 16:29:24 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:33:02 EST