[R-sig-ME] LME - varComb and varIdent

[Friso Muijsers]

Am 6/26/2013 5:29 PM, schrieb Friso Muijsers:
> Hello,
>
> I'm quite new to LMEs and have a question to which I did not find an 
> answer in the archives or in the P&B chapter 5 (allthough I have 
> problems undestanding the latter, fully).
> Im trying to fit a linear mixed model to my data (DV = numerical, IV = 
> numerical and factorial).
> I have some issues with variance heterogeneity within two of my factors:
>
> 1) experimental type, two levels: "lab" and "field")
> 2) system, two levels: "marine" and "limnic"
>
> When adding both variance structures to my models, my results differ 
> slightly, depending on whether I use "weights = 
> varIdent(form=~1|system*exp.type)" or "weights = 
> varComb(varIdent(form=~1|exp.type), varIdent(form=~1|system))".
> Which approach would be the better one? What does the " * " exactly 
> do? It somehow uses one additional df.
>
>                             Model     df      AIC BIC 
> logLik               Test       L.Ratio p-value
> lmaicresisa         1             14     127.1082     160.2804 -49.55408
> lmaicresisb         2             15     128.8289     164.3707 
> -49.41447     1 vs 2      0.279231         0.5972
>
>
> And a general question: I've seen (read) many people arguing that one 
> should not use the varFunc too excessively. With only one varIdent 
> Factor, my models have indeed much better p-values (all though I 
> understand that those are less relevant in LME) but my AIC and LogLik 
> increase. Should I use the more parsimonious model with both variance 
> factors (as indicated by AIC) ? Are the low p-values with 2 
> var-functions an indication of a bad model? QQ-Plots indicate a 
> slightly better model fit with 2 var-functions.
>
> It is difficult to add a nice example here, since my data are 
> relatively complex (meta-analysis). If it is necessary, I can try to 
> create a comparable dataset, allthough not sure how to.
>
> This is my model:
>
> lmresisc = 
> lme(resis.log~evenness+exp.type+exp.type:evenness+org.type.merged+org.type.merged:evenness+system+log(exact.duration)+system2,
>                        random =~1|authors.year,
>                        data = 
> data[(!is.na(data$evenness)&!is.na(data$resis.log)),], weights = 
> varIdent(form=~1|system*exp.type), method = "ML")
>         summary(lmresisc)
>         lmaicresisc = stepAIC(lmresisc)
>         summary(lmaicresisc)
>
> This is the last step of the stepAIC:
>
> Step:  AIC=128.83
> resis.log ~ evenness + exp.type + org.type.merged + system +
>     log(exact.duration) + evenness:exp.type + evenness:org.type.merged
>
>                                                    Df    AIC
> <none>                                               128.83
> - evenness:exp.type                1       129.38
> - system                                    1       132.48
> - evenness:org.type.merged  2        133.55
> - log(exact.duration)                1        134.20
>
>
> and this is the summary of the final model
>
>       Value                 Std.Error     DF    t-value p-value
> (Intercept) -1.963651         0.6959103     62     -2.8217014 0.0064
> evenness -2.231725         1.4724964     62     -1.5156064  0.1347
> exp.typelab -0.487326         0.9387122       7     -0.5191437  0.6197
> org.type.mergedheterotroph                   -2.245062 1.3380510       
> 7     -1.6778601  0.1373
> org.type.mergedmixed                            -2.494429 1.1662591   
>     7     -2.1388289  0.0698
> systemmarine 0.505430     0.2204285       7      2.2929430  0.0556
> log(exact.duration) 0.406731     0.1536198     62      2.6476456  0.0103
> evenness:exp.typelab                               2.295373 1.4621236 
>     62      1.5698901  0.1215
> evenness:org.type.mergedheterotroph   3.686948         1.6697807 
> 62      2.2080430  0.0309
> evenness:org.type.mergedmixed            5.930363 2.0642291     
> 62      2.8729191  0.0056
>
>
> I hope I gave enough information, please forgive me, if not. This is 
> my first question here, so i'm not sure about that.
> Thanks in advance!
>
> Friso
>
Sorry for messing up the table structure:

Another try:

Step:  AIC=128.83
resis.log ~ evenness + exp.type + org.type.merged + system +
     log(exact.duration) + evenness:exp.type + evenness:org.type.merged

                            Df    AIC
<none>                        128.83
- evenness:exp.type         1 129.38
- system                    1 132.48
- evenness:org.type.merged  2 133.55
- log(exact.duration)       1 134.20


and this is the summary of the final model

                                         Value Std.Error DF t-value p-value
(Intercept)                         -1.963651 0.6959103 62 -2.8217014  
0.0064
evenness                            -2.231725 1.4724964 62 -1.5156064  
0.1347
exp.typelab                         -0.487326 0.9387122  7 -0.5191437  
0.6197
org.type.mergedheterotroph          -2.245062 1.3380510  7 -1.6778601  
0.1373
org.type.mergedmixed                -2.494429 1.1662591  7 -2.1388289  
0.0698
systemmarine                         0.505430 0.2204285  7 2.2929430  0.0556
log(exact.duration)                  0.406731 0.1536198 62 2.6476456  0.0103
evenness:exp.typelab                 2.295373 1.4621236 62 1.5698901  0.1215
evenness:org.type.mergedheterotroph  3.686948 1.6697807 62 2.2080430  0.0309
evenness:org.type.mergedmixed        5.930363 2.0642291 62 2.8729191  0.0056

-- 
Friso Muijsers

Institute for Chemistry and Biology of the Marine Environment (ICBM)
Carl-von-Ossietzky University Oldenburg
Schleusenstrasse 1
26382 Wilhemshaven