Message-ID: <CAO7JsnQtmKLXoNHsf26x4i+6G7hJjyozJ+7MPnEF=mCsfNr1EA@mail.gmail.com>
Date: 2012-06-25T17:28:23Z
From: Douglas Bates
Subject: sas to R
In-Reply-To: <CALSKosAzMbibgn6nUKTc7swddvio_af0BQLpOygD_SPhPJb+MA@mail.gmail.com>
Are you trying to load both the nlme and the lme4 packages at the same
time? That can cause problems. You are better off fitting the lmer
model in one R session and the lme model in another.
On Mon, Jun 25, 2012 at 11:50 AM, Steve Hong <emptican at gmail.com> wrote:
> Thank all of you for replying to me.
>
> I tried lmer, lme, and SAS. ?I was able to get outputs when I use
> 'lme' whereas no results from 'lmer'. ?I don't know why. ?Does anyone
> know what the warning message mean? ?Outputs from ?'lme' were similar
> with those from SAS. ?Below is selected outputs from lmer, lme, and
> SAS, FYI.
>
> Thanks again,
>
> Steve Hong
>
>> fm.lmer <- lmer(y ~ trt + (1|trial/block/trt), data=df, na.action=na.omit)
> Error: length(f1) == length(f2) is not TRUE
> In addition: Warning messages:
> 1: In block:trial :
> ? numerical expression has 92 elements: only the first used
> 2: In block:trial :
> ? numerical expression has 92 elements: only the first used
> 3: In trt:(block:trial) :
> ? numerical expression has 92 elements: only the first used
> 4: In block:trial :
> ? numerical expression has 92 elements: only the first used
> 5: In block:trial :
> ? numerical expression has 92 elements: only the first used
>> fm.lme <- lme(y ~ trt, random=(~1|trial/block/trt), data = df, na.action=na.omit)
>> summary(fm.lme)
> Linear mixed-effects model fit by REML
> ?Data: df
> ? ? ? ? AIC ? ? ? BIC ? logLik
> ? -85.22388 -60.68041 52.61194
>
> Random effects:
> ?Formula: ~1 | trial
> ? ? ? ? (Intercept)
> StdDev: ? 0.1112442
>
> ?Formula: ~1 | block %in% trial
> ? ? ? ? ?(Intercept)
> StdDev: 1.449228e-06
>
> ?Formula: ~1 | trt %in% block %in% trial
> ? ? ? ? (Intercept) ?Residual
> StdDev: ?0.07081356 0.1020226
>
> Fixed effects: y ~ trt
> ? ? ? ? ? ? ? ? ? Value ?Std.Error DF ? ?t-value p-value
> (Intercept) ?0.24428523 0.08793775 56 ?2.7779337 ?0.0074
> trtau2 ? ? ?-0.00996643 0.05605221 25 -0.1778063 ?0.8603
> trtberm ? ? -0.12786905 0.05686903 25 -2.2484830 ?0.0336
> trtls44 ? ? ?0.12326637 0.05478364 25 ?2.2500582 ?0.0335
> trtsr10y5 ? ?0.02513355 0.05517460 25 ?0.4555275 ?0.6527
> trtsr10y6 ? ?0.01932992 0.05478364 25 ?0.3528410 ?0.7272
> ?Correlation:
> ? ? ? ? ? (Intr) trtau2 trtbrm trtl44 trt105
> trtau2 ? ?-0.314
> trtberm ? -0.309 ?0.486
> trtls44 ? -0.321 ?0.504 ?0.497
> trtsr10y5 -0.319 ?0.500 ?0.493 ?0.511
> trtsr10y6 -0.321 ?0.504 ?0.497 ?0.515 ?0.511
>
> Standardized Within-Group Residuals:
> ? ? ? ? ? Min ? ? ? ? ? ?Q1 ? ? ? ? ? Med ? ? ? ? ? ?Q3 ? ? ? ? ? Max
> -2.614096e+00 -5.666986e-01 -9.727356e-05 ?4.692685e-01 ?2.410879e+00
>
> Number of Observations: 92
> Number of Groups:
> ? ? ? ? ? ? ? ? ? ? trial ? ? ? ? ?block %in% trial trt %in% block %in% trial
> ? ? ? ? ? ? ? ? ? ? ? ? 2 ? ? ? ? ? ? ? ? ? ? ? ? 6 ? ? ? ? ? ? ? ? ? ? ? ?36
>> anova(fm.lme)
> ? ? ? ? ? ? numDF denDF ?F-value p-value
> (Intercept) ? ? 1 ? ?56 9.907983 ?0.0026
> trt ? ? ? ? ? ? 5 ? ?25 4.122070 ?0.0072
>
>
> SAS code and outputs:
> proc glimmix data=df;
> model y=trt;
> random trial block(trial) turf(block*turf);
> run;
>
> ? ? Covariance Parameter Estimates
>
> ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? Standard
> Cov Parm ? ? ? ? ? ? Estimate ? ? ? Error
>
> trial ? ? ? ? ? ? ? ? 0.01237 ? ? 0.01823
> block(trial) ? ? ? ? ? ? ? ?0 ? ? ? ? ? .
> trt(trial*block) ? ?0.005015 ? ?0.002546
> Residual ? ? ? ? ? ? ?0.01041 ? ?0.001963
>
>
> ? ? ? ?Type III Tests of Fixed Effects
>
> ? ? ? ? ? ? ?Num ? ? ?Den
> Effect ? ? ? ? DF ? ? ? DF ? ?F Value ? ?Pr > F
>
> trt ? ? ? ? ? ?5 ? ? ? 25 ? ? ? 4.12 ? ?0.0072
>
>
>
> On Mon, Jun 25, 2012 at 10:25 AM, Kevin Wright <kw.stat at gmail.com> wrote:
>>
>> This could be similar to a multi-location RCB design were "trial" is
>> location. ?Since no distribution is specified, the distribution is
>> assumed to be Gaussian. ?Make sure that trial, block, trt are factors,
>> this should be similar to SAS:
>>
>> lmer(y ~ trt + (1|trial/block/trt), data=df)
>>
>> > proc glimmix data=df;
>> > class trial block trt;
>> > model y=trt;
>> > random trial block(trial) trt(block*trial);
>>
>> Kevin Wright
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models