[R-sig-ME] Puzzled by simple example.

Rolf Turner r.turner at auckland.ac.nz
Wed Apr 1 01:31:11 CEST 2009


I continue to struggle to understand the syntax of lmer() from the  
lme4 package.

In an attempt to enlighten myself I simulated a rather simple data  
set.  The
imagined scenario is children's heights (in cm.) at ages 4, 5, 6, and  
7, although
no attempt has been made to make the data realistic.

Since the data are balanced, I thought I could check/compare the lmer 
() results
with results produced by aov().

The code I used is as follows:

#
# Script scr.01
#

# Generate the data.  No attempt is made to make these
# data realistic.
set.seed(42)
Beta <- c(40,50,60,70)
B  <- rnorm(80,0,2)
E  <- rnorm(80,0,0.25)
Y  <- Beta + B + E
ht.dat <- data.frame(y=Y,age=factor(rep(4:7,20)),child=factor(rep 
(1:20,each=4)))

# Fit with lmer():
f1 <- lmer(y ~ 0 + age + (1|child),data=ht.dat)

# Fit with aov():
f2 <- aov(y ~ 0 + age + Error(child),data=ht.dat)

# Clean up:
rm(Beta,B,E,Y)

When one does summary(f1) one gets values of the estimates of the  
mean heights
which are very close to the ``true'' values used to generate the  
data.  I am however
puzzled by the random effects estimates.  These are:

Random effects:
  Groups   Name        Variance   Std.Dev.
  child    (Intercept) 1.0362e-07 0.0003219
  Residual             4.6663e+00 2.1601657
Number of obs: 80, groups: child, 20

I thought at the ``child'' variance would be around 4 (the variance
used to simulate B) and that the residual variance would be around  
0.0625
(the variance used to simulate E).  This is clearly not so, so there is
clearly something I am misunderstanding here.  I am not interpreting the
random effects results correctly.  Can someone please point me in the
right direction?

In particular, how can I extract estimates of the variance of B and E?

I am assuming that the model I am fitting in both instances is

	y = X*Beta + Z*B + E

where X = kronecker(One_20,I_4) and Z = kronecker(I_20,One_4), where  
I_k is the
k-x-k identity and One_k is a k-x-1 matrix all of whose entries are 1.
Have I got this much right?  (In the forgoing Beta is a 4-vector of  
fixed effects
and B is a 20-vector of random effects.)

I am even more puzzled by the results of aov() --- I cannot connect  
anything
with anything from these results.  Where/how do I obtain estimates of  
the
variance of B and the variance of E?  Where/how do I obtain estimates  
of the
fixed effects?

If I do coeff(f2) I get

child :
     age4
220.0565

Within :
      age4      age5      age6
-29.86516 -20.91547 -11.09382

The ``Within'' coefficients are the differences age4 - age7, age5 -  
age7, age6 - age7,
respectively.  I don't understand the ``child'' coefficient (age4 =  
220.0565) at all.
Can anyone explain this to me?

One other puzzlement:  If one just types ``f2'' or ``print(f2)'' one  
gets messages,
for stratum 1 (child) ``3 out of 4 effects not estimable / Estimable  
effects are
balanced'', and for stratum 2 (Within) ``1 out of 4 effects not  
estimable / Estimated
effects may be unbalanced''  Huh?  What is going on?  How could a  
design possibly
be more balanced?

I would be grateful for some enlightenment.

	cheers,

		Rolf Turner

######################################################################
Attention:\ This e-mail message is privileged and confid...{{dropped:9}}




More information about the R-sig-mixed-models mailing list