Hi, all,

 

I have a question about random effects model. I am dealing with a
three-factor experiment dataset. The response variable y is modeled
against three factors: Samples, Operators, and Runs. The experimental
design is as follow:

 

4 samples were randomly chosen from a large pool of test samples. Each
of the 4 samples was analyzed by 4 operators, randomly selected from a
group of operators. Each operator independently analyzed same samples
over 5 runs (runs nested in operator). I would like to know the
following things:

 

(1)                     the standard deviation within each run;

(2)                     the standard deviation between runs;

(3)                     the standard deviation within operator

(4)                     the standard deviation between operator.

 

With this data, I assumed the three factors are all random effects. So
the model I am looking for is

 

Model:  y  = Samples(random) + Operator(random) + Operator:Run(random) +
Error(Operator) + Error(Operator:Run)  + Residuals

 

I am using lme function in nlme package. Here is the R code I have

 

1.	 using lme:

First I created a grouped data using

gx <- groupedData(y ~ 1 | Sample, data=x)

gx$dummy <- factor(rep(1,nrow(gx)))

 

then I run the lme

 

fm<- lme(y ~ 1, data=gx,
random=list(dummy=pdBlocked(list(pdIdent(~Sample-1),

            pdIdent(~Operator-1), 

            pdIdent(~Operator:Run-1)))))

 

    finally, I use VarCorr to extract the variance components

 

           vc <- VarCorr(fm)

 

                     Variance           StdDev  

Operator:Run 1.595713e-10:20   1.263215e-05:20  

Sample       5.035235e+00: 4   2.243933e+00: 4  

Operator     5.483145e-04: 4   2.341612e-02: 4  

Residuals    8.543601e-02: 1   2.922944e-01: 1  

 

 

2.	Using lmer in Matrix package:

 

fm <- lmer(y ~ (1 | Sample) + (1 | Operator) + 

           (1|Operator:Run), data=x)

     summary(fm)

 

Linear mixed-effects model fit by REML 

Formula: H.I.Index ~ (1 | Sample.Name) + (1 | Operator) + (1 |
Operator:Run) 

          Data: x 

      AIC      BIC    logLik MLdeviance REMLdeviance

 96.73522 109.0108 -44.36761   90.80064     88.73522

Random effects:

 Groups       Name        Variance   Std.Dev.  

 Operator:Run (Intercept) 4.2718e-11 6.5359e-06

 Operator     (Intercept) 5.4821e-04 2.3414e-02

 Sample       (Intercept) 5.0352e+00 2.2439e+00

 Residual                 8.5436e-02 2.9229e-01

number of obs: 159, groups: Operator:Run, 20; Operator, 4; Sample.Name,
4

 

Fixed effects:

             Estimate Std. Error  t value

(Intercept) 0.0020818  1.1222683 0.001855

 

 

There is a difference between lmer and lme is for the factor
Operator:Run.  I cannot find where the problem is. Could anyone point me
out if my model specification is correct for the problem I am dealing
with. I am pretty new user to lme and lmer. Thanks for your help!

 

 

Wilson Wang

 

 


	[[alternative HTML version deleted]]

