Skip to content
Prev 67858 / 398506 Next

help on extract variance components from the fitted model by lm

On Sun, 2005-04-17 at 02:38 +0800, wenqing li wrote:
First, the use of as.vector() above is redundant:
rep(5,5)),2))
num [1:50] 1 1 1 1 1 2 2 2 2 2 ...
num [1:50] 1 1 1 1 1 2 2 2 2 2 ...
[1] TRUE


On your question, a couple of options:

# First create a data frame
df <- data.frame(factor(A), factor(B), y)

library(nlme)

# Set your factors as the nested random effects
lme1 <- lme(y ~ 1, random = ~ 1 | A / B, data = df)
Linear mixed-effects model fit by REML
 Data: df
       AIC      BIC    logLik
  147.4539 155.0211 -69.72693

Random effects:
 Formula: ~1 | A
        (Intercept)
StdDev:    3.797784

 Formula: ~1 | B %in% A
        (Intercept)  Residual
StdDev:   0.3768259 0.6957023

Fixed effects: y ~ 1
            Value Std.Error DF t-value p-value
(Intercept)    11  1.702937 25 6.45943       0

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max
-1.7696293 -0.7042397  0.1521976  0.5708701  1.5679386

Number of Observations: 50
Number of Groups:
       A B %in% A
       5       25



Also:

library(ape)
A          B     Within
14.4231663  0.1419978  0.4840016
attr(,"class")
[1] "varcomp"


Note in both cases, lme() is used, not lm().

These are referenced in Section 10.2 (pg 279) of MASS4 by Venables &
Ripley and in Section 4.2.3 (pg 167) of MEMSS by Pinheiro and Bates.

HTH,

Marc Schwartz