An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20120625/37a69900/attachment.pl>
sas to R
13 messages · Joshua Wiley, Kevin Wright, Steve Hong +3 more
Hi Steve, That SAS code looks a little strange to me (I feel like it is missing a distribution and some options on the random line). It would probably help if you could also describe in words what you are doing or trying to do. a rough pass at it would be: glmer(y ~ trt + (trt | block) + (trt | trial), data = df, family = "binomial") That will fit a logistic model, with a random intercept and trt effect by block and by trial. If they are nested, that is fine. If they are cross classified, that is fine too. The random intercept and trt effect will be correlated by default. Cheers, Josh
On Mon, Jun 25, 2012 at 7:50 AM, Steve Hong <emptican at gmail.com> wrote:
Dear all, Can anyone help me convert SAS code to R in lme or lmer, especially in random statement part? Here is SAS code: proc glimmix data=df; class trial block trt; model y=trt; random trial block(trial) trt(block*trial); Thanks in advance, Steve ? ? ? ?[[alternative HTML version deleted]]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
Joshua Wiley Ph.D. Student, Health Psychology Programmer Analyst II, Statistical Consulting Group University of California, Los Angeles https://joshuawiley.com/
An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20120625/42fc7ba0/attachment.pl>
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
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
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
Hi Prof. Bates, Thanks for replying. Yes, I loaded both packages together and ran them. I don't understand different/separate R 'session'? Obviously, it seems not different versions (e.g., 2.15.0 vs. 2.14.0). Could you rephrase what you meant by R 'session'? Thanks, Steve
On Mon, Jun 25, 2012 at 12:28 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
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
On Mon, Jun 25, 2012 at 12:41 PM, Steve Hong <emptican at gmail.com> wrote:
Hi Prof. Bates, Thanks for replying. ?Yes, I loaded both packages together and ran them. ?I don't understand different/separate R 'session'? ?Obviously, it seems not different versions (e.g., 2.15.0 vs. 2.14.0). ?Could you rephrase what you meant by R 'session'?
I mean to run R, load lme4 and fit the model. Then quit R and restart it, load nlme and fit that model.
On Mon, Jun 25, 2012 at 12:28 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
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
I restarted R and tried to fit the model. However, I got the same message... So, I checked the session information with sessionInfo(). I found there is still nlme package in "loaded via a namespace (and not attached)". Is it still causing problem? Here is what I did: R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. [Previously saved workspace restored]
df <- bcwrear1 y <- df$day10 trial <- df$trial block <- df$block trt <- df$turf library(lme4)
Loading required package: Matrix
Loading required package: lattice
Attaching package: ?lme4?
The following object(s) are masked from ?package:stats?:
AIC, BIC
fm.lmer <- lmer(log10(day10) ~ 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
sessionInfo()
R version 2.15.0 (2012-03-30) Platform: i386-pc-mingw32/i386 (32-bit) locale: [1] LC_COLLATE=English_United States.1252 [2] LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 [4] LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] lme4_0.999375-42 Matrix_1.0-6 lattice_0.20-6 loaded via a namespace (and not attached): [1] grid_2.15.0 nlme_3.1-103 stats4_2.15.0 Thanks!!! Steve
On Mon, Jun 25, 2012 at 1:01 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
On Mon, Jun 25, 2012 at 12:41 PM, Steve Hong <emptican at gmail.com> wrote:
Hi Prof. Bates, Thanks for replying. ?Yes, I loaded both packages together and ran them. ?I don't understand different/separate R 'session'? ?Obviously, it seems not different versions (e.g., 2.15.0 vs. 2.14.0). ?Could you rephrase what you meant by R 'session'?
I mean to run R, load lme4 and fit the model. ?Then quit R and restart it, load nlme and fit that model.
On Mon, Jun 25, 2012 at 12:28 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
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
An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20120625/4a090771/attachment.pl>
Thank you, all. I restarted R and checked package list using 'search()'. And then, loaded 'lme4' and ran the model. However, the results are same... Below is what I did. Thank you much again!!! Steve R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. [Previously saved workspace restored]
search()
[1] ".GlobalEnv" "package:stats" "package:graphics" [4] "package:grDevices" "package:utils" "package:datasets" [7] "package:methods" "Autoloads" "package:base"
library(lme4)
Loading required package: Matrix
Loading required package: lattice
Attaching package: ?lme4?
The following object(s) are masked from ?package:stats?:
AIC, BIC
search()
[1] ".GlobalEnv" "package:lme4" "package:Matrix" [4] "package:lattice" "package:stats" "package:graphics" [7] "package:grDevices" "package:utils" "package:datasets" [10] "package:methods" "Autoloads" "package:base"
bcwrear1 <- read.table("C:/bcwrear1.txt", header=T)
df <- bcwrear1
df=bcwrear1
df=transform(df, y=day10, trt=turf)
fm.lmer <- lmer(log10(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
On Mon, Jun 25, 2012 at 1:36 PM, Jake Westfall <jake987722 at hotmail.com> wrote:
You can unload nlme by using:
detach("package:nlme")
Jake
Date: Mon, 25 Jun 2012 13:12:52 -0500 From: emptican at gmail.com To: bates at stat.wisc.edu CC: r-sig-mixed-models at r-project.org Subject: Re: [R-sig-ME] sas to R I restarted R and tried to fit the model. ?However, I got the same message... ?So, I checked the session information with sessionInfo(). I found there is still nlme package in "loaded via a namespace (and not attached)". ?Is it still causing problem? Here is what I did: R version 2.15.0 (2012-03-30) Copyright (C) 2012 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i386-pc-mingw32/i386 (32-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. ? Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. [Previously saved workspace restored]
df <- bcwrear1 y <- df$day10 trial <- df$trial block <- df$block trt <- df$turf library(lme4)
Loading required package: Matrix Loading required package: lattice Attaching package: ?lme4? The following object(s) are masked from ?package:stats?: ? ? AIC, BIC
fm.lmer <- lmer(log10(day10) ~ 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
sessionInfo()
R version 2.15.0 (2012-03-30) Platform: i386-pc-mingw32/i386 (32-bit) locale: [1] LC_COLLATE=English_United States.1252 [2] LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 [4] LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] stats ? ? graphics ?grDevices utils ? ? datasets ?methods ? base other attached packages: [1] lme4_0.999375-42 Matrix_1.0-6 ? ? lattice_0.20-6 loaded via a namespace (and not attached): [1] grid_2.15.0 ? nlme_3.1-103 ?stats4_2.15.0 Thanks!!! Steve On Mon, Jun 25, 2012 at 1:01 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
On Mon, Jun 25, 2012 at 12:41 PM, Steve Hong <emptican at gmail.com> wrote:
Hi Prof. Bates, Thanks for replying. ?Yes, I loaded both packages together and ran them. ?I don't understand different/separate R 'session'? ?Obviously, it seems not different versions (e.g., 2.15.0 vs. 2.14.0). ?Could you rephrase what you meant by R 'session'?
I mean to run R, load lme4 and fit the model. ?Then quit R and restart it, load nlme and fit that model.
On Mon, Jun 25, 2012 at 12:28 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
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
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
? ? ? ?[[alternative HTML version deleted]]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
Steve Hong <emptican at ...> writes:
Thank you, all. I restarted R and checked package list using 'search()'. And then, loaded 'lme4' and ran the model. However, the results are same... Below is what I did. Thank you much again!!! Steve
As Kevin Wright said: Make sure that trial, block, trt are factors, [snip snip snip]
bcwrear1 <- read.table("C:/bcwrear1.txt", header=T)
df <- bcwrear1
df=bcwrear1
df=transform(df, y=day10, trt=turf)
df <- transform (trial=factor(trial), block=factor(block), trt=factor(trt))
fm.lmer <- lmer(log10(y) ~ trt + (1|trial/block/trt), data=df,
na.action=na.omit) [snip snip snip]
Thank all of you. Ben and Kevin, you were right! The main reason of my problem that I experienced was not specifying factors (i.e., trial and block)... Once I factorized them, I was able to run the models with both lme and lmer in the SAME R session. I think it was possible since I did not use functions in lmer such as lmList. Also, the results, especially df and F-stats, were same or very close to results from SAS. I greatly appreciate all of your kind/patient suggestions. Thank you again!!! Steve
On Tue, Jun 26, 2012 at 2:10 AM, Ben Bolker <bbolker at gmail.com> wrote:
Steve Hong <emptican at ...> writes:
Thank you, all. I restarted R and checked package list using 'search()'. ?And then, loaded 'lme4' and ran the model. ?However, the results are same... Below is what I did. Thank you much again!!! Steve
As Kevin Wright said: Make sure that trial, block, trt are factors, ?[snip snip snip]
bcwrear1 <- read.table("C:/bcwrear1.txt", header=T)
df <- bcwrear1
df=bcwrear1
df=transform(df, y=day10, trt=turf)
df <- transform (trial=factor(trial), block=factor(block), trt=factor(trt))
fm.lmer <- lmer(log10(y) ~ trt + (1|trial/block/trt), data=df,
?na.action=na.omit) ?[snip snip snip]
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models