[R-sig-ME] testing significance of fixed factors in a mixed model

Julia Chacon Labella juliachacon at gmail.com
Sat Dec 16 19:05:31 CET 2017


Thank you very much for your quick answer. Here is the coding of the model:

m1_y <- lmer(log(total_yield) ~ bench +  dom + MoMix + Divalias + fgalias +
               dom:MoMix + dom:Divalias + dom:fgalias+
         fgcomb + dom:fgcomb + MoMix:fgcomb + Divalias:fgcomb +
fgalias:fgcomb+
         dom:MoMix:fgcomb + dom:Divalias:fgcomb + dom:fgalias:fgcomb+
              (1|spcomb),
            data=data_wo_s)

bench: three levels (a,b,c)
dom: has two levels (c and w)
Now I have the treatments. This is like a tree, where I test three
variables that are included in the previous one.
MoMix: has two levels (mono and mix)
    Divalias: has two levels (2 and 4). This two levels are available
within Mix
       Fgalias: has two levels (1 and 2). This two are only and are
included within level 2 in Divalias
So these three are not orthogonal.
fgcomb has 9 levels.
In the random part I have spcomb that are species combinations and have
many levels.
I copied the summary of the model at the end of the email.

NOW, the anovas!
> anova(m1_y, ddf="Kenward-Roger", type=1)

Analysis of Variance Table of type I  with  Kenward-Roger
approximation for degrees of freedom
                  Sum Sq  Mean Sq NumDF   DenDF F.value    Pr(>F)
bench            0.00037 0.000187     2 143.957  0.0043    0.9957
dom              0.02235 0.022345     1 297.576  0.5147    0.4737
MoMix            0.09241 0.092411     1  32.982  2.1288    0.1540
Divalias         0.01255 0.012553     1  37.111  0.2892    0.5940
fgalias          0.00956 0.009564     1  32.410  0.2203    0.6419
fgcomb           2.15832 0.269790     8  34.559  6.2148 5.608e-05 ***
dom:MoMix        0.11248 0.112478     1 298.040  2.5910    0.1085
dom:Divalias     0.00017 0.000174     1 297.846  0.0040    0.9495
dom:fgalias      0.09687 0.096866     1 298.607  2.2314    0.1363
dom:fgcomb       0.58166 0.072707     8 294.133  1.6749    0.1040
MoMix:fgcomb     0.02578 0.008592     3  32.124  0.1979    0.8970
dom:MoMix:fgcomb 0.21952 0.073172     3 292.107  1.6856    0.1702
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

This is the output when I do a type I anova in lmerTest. I was doing a type
I anova, because I wanted to "eliminate" or control the effect of bench.
So, the other factors are tested after and there is less variance to be
explained, that is less from the benches of the greenhouse. But, now I am
not sure if this is the best way to proceed.
Because Type I test the terms in the order of appearance and the way it
calculates the sum of squares, you only get the combined effect of the
factor and the below ones. So...with a Type I anova I am not going to see
the effect of each single factor (isn't it?).
Now I think I should go for a type II. But, not completely sure about that
and as you will see the output is completely different.

Here the type 2.

> anova(m1_y, ddf="Kenward-Roger", type=2)
Analysis of Variance Table of type II  with  Kenward-Roger
approximation for degrees of freedom
                  Sum Sq Mean Sq NumDF   DenDF F.value    Pr(>F)
bench            0.14694 0.07347     2 299.947  1.6924  0.185832
dom              0.01050 0.01050     1 292.202  0.2419  0.623175
MoMix            0.01389 0.01389     1  32.110  0.3199  0.575592
Divalias         0.31910 0.31910     1  32.219  7.3508  0.010660 *
fgalias          0.52884 0.52884     1  32.650 12.1823  0.001405 **
fgcomb           2.16056 0.27007     8  34.559  6.2213 5.556e-05 ***
dom:MoMix        0.00350 0.00350     1 292.154  0.0806  0.776621
dom:Divalias     0.24601 0.24601     1 292.329  5.6670  0.017928 *
dom:fgalias      0.26688 0.26688     1 292.277  6.1477  0.013722 *
dom:fgcomb       0.58182 0.07273     8 292.227  1.6753  0.103916
MoMix:fgcomb     0.02578 0.00859     3  32.124  0.1979  0.897033
dom:MoMix:fgcomb 0.21952 0.07317     3 292.107  1.6856  0.170192
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


I am interested in the main effects but also in the interactions, of
course. So, should I test this with a type II or type I anova?
I was also wondering if there is any way of escaping from this, hahaha, was
wondering if going to a Confidence interval approach is a much more elegant
way of solving this...

Thanks so much. I am so sorry for the long email.
Julia

---------------------------------------------------------------------------------------
summary(m1_y)

REML criterion at convergence: 66.4

Scaled residuals:
     Min       1Q   Median       3Q      Max
-2.72351 -0.57177  0.07539  0.58399  2.81317

Random effects:
 Groups   Name        Variance Std.Dev.
 spcomb   (Intercept) 0.07007  0.2647
 Residual             0.04341  0.2084
Number of obs: 361, groups:  spcomb, 52

Fixed effects:
                        Estimate Std. Error        df t value Pr(>|t|)
(Intercept)              2.58512    0.27555  35.22000   9.382 4.12e-11 ***
benchb                   0.05650    0.03070 301.72000   1.840  0.06669 .
benchc                   0.03420    0.03058 301.94000   1.118  0.26441
domw                    -0.25077    0.10469 292.80000  -2.395  0.01723 *
MoMixMon                -0.11301    0.33789  35.40000  -0.334  0.74001
Divalias4                0.61797    0.28504  35.74000   2.168  0.03690 *
fgalias2                 0.88760    0.31157  36.85000   2.849  0.00714 **
fgcombc4                 1.11869    0.38961  35.20000   2.871  0.00688 **
fgcombfb                 0.36312    0.38871  34.88000   0.934  0.35663
fgcomblg                 0.73868    0.38870  34.87000   1.900  0.06567 .
fgcombc4c3              -0.22859    0.20748  45.19000  -1.102  0.27639
fgcombfbc3              -0.55140    0.20666  44.50000  -2.668  0.01060 *
fgcombfbc4              -0.17922    0.20698  44.78000  -0.866  0.39117
fgcombfblg              -0.30621    0.20665  44.50000  -1.482  0.14544
fgcomblgc3              -0.63072    0.20743  45.15000  -3.041  0.00392 **
domw:MoMixMon            0.05878    0.13489 292.76000   0.436  0.66331
domw:Divalias4           0.27973    0.11750 292.99000   2.381  0.01792 *
domw:fgalias2            0.33864    0.13657 292.94000   2.480  0.01372 *
domw:fgcombc4            0.16957    0.14999 292.71000   1.131  0.25917
domw:fgcombfb            0.33597    0.14769 292.75000   2.275  0.02364 *
domw:fgcomblg            0.36547    0.14753 292.73000   2.477  0.01380 *
domw:fgcombc4c3          0.06844    0.12378 293.12000   0.553  0.58072
domw:fgcombfbc3         -0.04844    0.12231 292.95000  -0.396  0.69237
domw:fgcombfbc4          0.04190    0.12282 293.02000   0.341  0.73323
domw:fgcombfblg         -0.07774    0.12201 292.92000  -0.637  0.52449
domw:fgcomblgc3         -0.13964    0.12339 293.12000  -1.132  0.25871
MoMixMon:fgcombc4       -0.02594    0.47926  35.82000  -0.054  0.95713
MoMixMon:fgcombfb        0.23764    0.47830  35.53000   0.497  0.62236
MoMixMon:fgcomblg       -0.10096    0.47803  35.45000  -0.211  0.83395
domw:MoMixMon:fgcombc4   0.14393    0.19339 292.73000   0.744  0.45732
domw:MoMixMon:fgcombfb  -0.21472    0.19170 292.77000  -1.120  0.26360
domw:MoMixMon:fgcomblg  -0.22705    0.19048 292.72000  -1.192  0.23424



2017-12-15 3:51 GMT+01:00 Ben Bolker <bbolker at gmail.com>:

> Can you show us some of your output?
>
> In general, an F-test will be more accurate than an LRT in cases where
> you can figure out how to do it, because an LRT doesn't make any
> finite-size corrections (what df values are you getting from
> KR/Satterthwaite)?
>
> The thing you will have to think hardest about is which variables are
> conditioned on when testing, especially when you have interactions in
> your model and are considering testing main effects.
>
> On Thu, Dec 14, 2017 at 2:25 PM, Julia Chacon Labella
> <juliachacon at gmail.com> wrote:
> > I am sorry if I am asking a naïve question or if that has been already
> > asked many times. I am not an expert in mixed models at all, but want to
> > understand and be confident about what I am doing.
> >
> > Actually, I am analysing a data set of a experiment. I have several
> > treatments, a balance design, randomized by blocks, etc., and generally
> > following the recommendations of the FAQs by Bolker, what are really
> > usefull (Thanks!). http://bbolker.github.io/
> mixedmodels-misc/glmmFAQ.html
> >
> > But, I am pretty confused about how to test the significance of fixed
> > factors in a mixed model.
> > I find huge differences in the significance output when computing 1)
> Anova
> > from the car package or 2) anova from lmerTest with KR or satterthwaite
> > approaches (both, KR or Satterthwaite have similar results), or 3) a LRT
> > via anova. The biggest differences are found for car::Anova!
> >
> > * First, I am not sure if I should employ a LRT or a conditional F-test
> > (lmerTest::anova, KR or satterthwaite approaches). Both results are
> pretty
> > similar in my case, but conditional F-test seem to slightly be more
> > conservative in my case. Can the F-test be somehow problematic in terms
> of
> > statistical power?
> >
> > * Second, in case of using lmerTest::anova, I am not completely sure if I
> > should used a type I or type III anova.  There is any reason for using
> type
> > I anova in balanced designs?
> >
> > Thanks in advance,
> > Julia
> >
> >         [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-mixed-models at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>



-- 
Julia

	[[alternative HTML version deleted]]



More information about the R-sig-mixed-models mailing list