[R] Help please! How to code a mixed-model with 2 within-subject factors using lme or lmer?

[Mark Difford]

Hi Roberto,

The other thing you can do --- if you don't wish to step across to lmer(),
where you will be able to exactly replicate the crossed-factor error
structure --- is stay with aov(... + Error()), but fit the factor you are
interested in last. Assume it is Sex. Then fit your model as

aov.model <- aov(Volume ~ Lobe * Tissue * Sex + Error(Subject/(Lobe *
Tissue))

This should give you a so-called "Type-II" test for Sex. You may verify this
by fitting the model without the Error term and using Anova() from the car
package (which does Type-II/III tests) to look at the SS and F values.

I say should, because the only concern I have is whether this procedure is
affected by the presence of an Error term in the model. Establishing this is
beyond my capabilities.

Regards, Mark.


roberto toro wrote:
> 
> Thanks for answering Mark!
> 
> I tried with the coding of the interaction you suggested:
> 
>> tfac<-with(vlt,interaction(Lobe,Tissue,drop=T))
>> mod<-lme(Volume~Sex*Lobe*Tissue,random=~1|Subject/tfac,data=vlt)
> 
> But is it normal that the DF are 2303? DF is 2303 even for the estimate of
> LobeO that has only 662 values (331 for Tissue=white and 331 for
> Tissue=grey).
> I'm not sure either that Sex, Lobe and Tissue are correctly handled....
> why are
> there different estimates called Sex:LobeO, Sex:LobeP, etc, and not just
> Sex:Lobe as with aov()?. Why there's Tissuew, but not Sex1, for example?
> 
> Thanks again!
> roberto
> 
> ps1. How would you code this with lmer()?
> ps2. this is part of the output of mod<-lme:
>> summary(mod)
> Linear mixed-effects model fit by REML
>  Data: vlt
>        AIC      BIC    logLik
>   57528.35 57639.98 -28745.17
> 
> Random effects:
>  Formula: ~1 | Subject
>         (Intercept)
> StdDev:    11294.65
> 
>  Formula: ~1 | tfac %in% Subject
>         (Intercept) Residual
> StdDev:    10569.03 4587.472
> 
> Fixed effects: Volume ~ Sex * Lobe * Tissue
>                        Value Std.Error   DF    t-value p-value
> (Intercept)        245224.61  1511.124 2303  162.27963  0.0000
> Sex                  2800.01  1866.312  329    1.50029  0.1345
> LobeO             -180794.83  1526.084 2303 -118.46975  0.0000
> LobeP             -131609.27  1526.084 2303  -86.23984  0.0000
> LobeT              -73189.97  1526.084 2303  -47.95932  0.0000
> Tissuew            -72461.05  1526.084 2303  -47.48168  0.0000
> Sex:LobeO            -663.27  1884.789 2303   -0.35191  0.7249
> Sex:LobeP           -2146.08  1884.789 2303   -1.13863  0.2550
> Sex:LobeT            1379.49  1884.789 2303    0.73191  0.4643
> Sex:Tissuew          5387.65  1884.789 2303    2.85849  0.0043
> LobeO:Tissuew       43296.99  2158.209 2303   20.06154  0.0000
> LobeP:Tissuew       50952.21  2158.209 2303   23.60856  0.0000
> LobeT:Tissuew      -15959.31  2158.209 2303   -7.39470  0.0000
> Sex:LobeO:Tissuew   -5228.66  2665.494 2303   -1.96161  0.0499
> Sex:LobeP:Tissuew   -1482.83  2665.494 2303   -0.55631  0.5781
> Sex:LobeT:Tissuew   -6037.49  2665.494 2303   -2.26506  0.0236
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: http://www.nabble.com/Help-please%21-How-to-code-a-mixed-model-with-2-within-subject-factors-using-lme-or-lmer--tp19480815p19489323.html
Sent from the R help mailing list archive at Nabble.com.