Not on my computer:

lmer(score~trt+(1|id/eye),dat)
Linear mixed model fit by REML
Formula: score ~ trt + (1 | id/eye)
   Data: dat
   AIC   BIC logLik deviance REMLdev
 425.2 474.2 -201.6    412.7   403.2
Random effects:
 Groups   Name        Variance Std.Dev.
 eye:id   (Intercept) 3.59532  1.89613
 id       (Intercept) 3.51024  1.87356
 Residual             0.01875  0.13693
Number of obs: 640, groups: eye:id, 160; id, 80

Fixed effects:
            Estimate Std. Error t value
(Intercept)   2.5500     0.7287   3.499
trtB          0.9875     1.0305   0.958
trtC          0.8625     1.0305   0.837
trtControl    1.1500     1.0305   1.116
trtD          2.4125     1.0305   2.341
trtE          1.4375     1.0305   1.395
trtF          2.1750     1.0305   2.111
trtG         -2.5500     1.0305  -2.474

Correlation of Fixed Effects:
           (Intr) trtB   trtC   trtCnt trtD   trtE   trtF
trtB       -0.707
trtC       -0.707  0.500
trtControl -0.707  0.500  0.500
trtD       -0.707  0.500  0.500  0.500
trtE       -0.707  0.500  0.500  0.500  0.500
trtF       -0.707  0.500  0.500  0.500  0.500  0.500
trtG       -0.707  0.500  0.500  0.500  0.500  0.500  0.500
>


On Fri, Sep 17, 2010 at 3:39 PM, array chip <arrayprofile@yahoo.com> wrote:

> Hi, I have a dataset of animals receiving some eye treatments. There are 8
> treatments, each animal's right and left eye was measured with some scores
> (ranging from 0 to 7) 4 times after treatment. So there are nesting groups
> eyes
> within animal. Dataset attached
>
> > dat<-read.table("dat.txt",sep='\t',header=T,row.names=1)
> > dat$id<-factor(dat$id)
> > str(dat)
> 'data.frame':   640 obs. of  5 variables:
>  $ score: int  0 2 0 7 4 7 0 2 0 7 ...
>  $ id   : Factor w/ 80 levels "1","3","6","10",..: 7 48 66 54 18 26 38 52
> 39 63
> ...
>  $ rep  : int  1 1 1 1 1 1 1 1 1 1 ...
>  $ eye  : Factor w/ 2 levels "L","R": 2 2 2 2 2 2 2 2 2 2 ...
>  $ trt  : Factor w/ 8 levels "A","B","C","Control",..: 1 1 1 1 1 1 1 1 1 1
> ...
>
> I fit a mixed model using both lmer() from lme4 package and lme() from nlme
> package:
>
> > lmer(score~trt+(1|id/eye),dat)
>
> Linear mixed model fit by REML
> Formula: score ~ trt + (1 | id/eye)
>   Data: dat
>   AIC   BIC logLik deviance REMLdev
>  446.7 495.8 -212.4    430.9   424.7
> Random effects:
>  Groups   Name        Variance   Std.Dev.
>  eye:id   (Intercept) 6.9208e+00 2.630742315798
>  id       (Intercept) 1.4471e-16 0.000000012030
>  Residual             1.8750e-02 0.136930641909
> Number of obs: 640, groups: eye:id, 160; id, 80
>
> > summary(lme(score~trt, random=(~1|id/eye), dat))
>
> Linear mixed-effects model fit by REML
>  Data: dat
>       AIC      BIC    logLik
>  425.1569 474.0947 -201.5785
>
> Random effects:
>  Formula: ~1 | id
>        (Intercept)
> StdDev:    1.873576
>
>  Formula: ~1 | eye %in% id
>        (Intercept)  Residual
> StdDev:    1.896126 0.1369306
>
> As you can see, the variance components estimates of random effects are
> quite
> different between the 2 model fits. From the data, I know that the variance
> component for "id" can't be near 0, which is what lmer() fit produced, so I
> think the lme() fit is correct while lmer() fit is off. This can also be
> seen
> from AIC, BIC etc. lme() fit has better values than lmer() fit.
>
>
> I guess this might be due to lmer() didn't converge very well, is there
> anyway
> to adjust to make lmer() converge better to get similar results as lme()?
>
> Thanks
>
> John
>
>
>
> _______________________________________________
> R-sig-mixed-models@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>


-- 
Daniel Marcelino
Skype: dmsilv
http://bit.ly/pol4vc

	[[alternative HTML version deleted]]


