Data sheet notation and model structure for GLMM with 3 non-factorial factors
Sat, Sep 26, 2009 at 3:11 AM, Raldo Kruger <raldo.kruger at gmail.com> wrote:
Hi Douglas,
Many thanks for the input. I've run two analyses on the same dataset using 1) indicator columns and the 2) a single 'factor / treatment' column for the non-factorial design described in my previous e-mail, and the results were identical (great!).
However, I did the same for a dataset with a factorial design (N, G, N*G, i.e. there were plots with N, plots with G, and plots with both N and G), and the results for the main effects are identical, but the estimates for the interaction effects (N*G) are different between the two analyses (see below). Could you help me make sense of that please (i.e. which one is correct?) !
Generally when you have the possibility of having N and G combined you would treat the design as a two-factor two-level factorial. That is, one factor for presence or absence of G and another factor for presence or absence of N. You could treat it as a single factor with four levels (neither, G only, N only and both N and G) but, as you have seen you need to translate between the representations. In the two-factor, two-level factorial design, let a be the estimate of the main effect for G, b be the estimate of the main effect for N, and c be the interaction estimate. In your example a = 0.14929, b = 0.03766 and c = -0.31633. Then the estimated cell mean for the NG cell is a + b + c =
0.03766 + 0.14929 + (-0.31633)
[1] -0.12938
Thanks, Raldo With expanded treatment notation- Fixed effects: ? ? ? ? ? ? ?Estimate Std. Error z value Pr(>|z|) (Intercept) ? ?2.92060 ? ?0.23834 ?12.254 ?< 2e-16 *** N ? ? ? ? ? ? ?0.03766 ? ?0.03486 ? 1.080 ? 0.2801 G ? ? ? ? ? ? ?0.14929 ? ?0.03395 ? 4.397 1.10e-05 *** Yearthree ? ? -2.85449 ? ?0.10664 -26.768 ?< 2e-16 *** Yeartwo ? ? ? -1.88175 ? ?0.06844 -27.494 ?< 2e-16 *** N:G ? ? ? ? ? -0.31633 ? ?0.04953 ?-6.386 1.70e-10 *** N:Yearthree ? ?0.15710 ? ?0.14428 ? 1.089 ? 0.2762 N:Yeartwo ? ? ?0.14736 ? ?0.09305 ? 1.584 ? 0.1133 G:Yearthree ? -0.25107 ? ?0.15430 ?-1.627 ? 0.1037 G:Yeartwo ? ? ?0.07550 ? ?0.09200 ? 0.821 ? 0.4118 N:G:Yearthree ?0.36353 ? ?0.20810 ? 1.747 ? 0.0807 . N:G:Yeartwo ? -0.01158 ? ?0.12996 ?-0.089 ? 0.9290 With single column treatment notation- Fixed effects: ? ? ? ? ? ? ? ? ? Estimate Std. Error z value Pr(>|z|) (Intercept) ? ? ? ? 2.92057 ? ?0.23836 ?12.253 ?< 2e-16 *** TreatG ? ? ? ? ? ? ?0.14928 ? ?0.03395 ? 4.397 1.10e-05 *** TreatN ? ? ? ? ? ? ?0.03767 ? ?0.03486 ? 1.080 0.279928 TreatNG ? ? ? ? ? ?-0.12938 ? ?0.03639 ?-3.556 0.000377 *** Yearthree ? ? ? ? ?-2.85448 ? ?0.10664 -26.768 ?< 2e-16 *** Yeartwo ? ? ? ? ? ?-1.88175 ? ?0.06844 -27.494 ?< 2e-16 *** TreatG:Yearthree ? -0.25109 ? ?0.15430 ?-1.627 0.103693 TreatN:Yearthree ? ?0.15711 ? ?0.14428 ? 1.089 0.276199 TreatNG :Yearthree ?0.26959 ? ?0.14636 ? 1.842 0.065483 . TreatG:Yeartwo ? ? ?0.07549 ? ?0.09200 ? 0.820 0.411941 TreatN:Yeartwo ? ? ?0.14735 ? ?0.09305 ? 1.583 0.113308 TreatNG :Yeartwo ? ?0.21118 ? ?0.09558 ? 2.210 0.027139 * On Thu, Sep 24, 2009 at 2:10 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
On Thu, Sep 24, 2009 at 1:22 AM, Raldo Kruger <raldo.kruger at gmail.com> wrote:
Hi R users,
I have 3 factors in a non-factorial design (G, K and N), as well as
two time periods (Year) and a random factor (Site), with Plant numbers
as the response variable.
My 1st question relates to the the notation of the treatments in the
data frame. Is it appropriate to use an expanded treatment notation,
such as this, when using glmer{lme4}:
Site ? ?Year ? ?Plant ? G ? ? ? K ? ? ? N
A ? ? ? 1 ? ? ? 5 ? ? ? 0 ? ? ? 0 ? ? ? 0
A ? ? ? 1 ? ? ? 4 ? ? ? 1 ? ? ? 0 ? ? ? 0
A ? ? ? 1 ? ? ? 7 ? ? ? 0 ? ? ? 1 ? ? ? 0
A ? ? ? 1 ? ? ? 10 ? ? ?0 ? ? ? 0 ? ? ? 1
A ? ? ? 2 ? ? ? 3 ? ? ? 0 ? ? ? 0 ? ? ? 0
A ? ? ? 2 ? ? ? 4 ? ? ? 1 ? ? ? 0 ? ? ? 0
A ? ? ? 2 ? ? ? 8 ? ? ? 0 ? ? ? 1 ? ? ? 0
A ? ? ? 2 ? ? ? 12 ? ? ?0 ? ? ? 0 ? ? ? 1
B ? ? ? 1 ? ? ? 7 ? ? ? 0 ? ? ? 0 ? ? ? 0
B ? ? ? 1 ? ? ? 3 ? ? ? 1 ? ? ? 0 ? ? ? 0
B ? ? ? 1 ? ? ? 7 ? ? ? 0 ? ? ? 1 ? ? ? 0
B ? ? ? 1 ? ? ? 12 ? ? ?0 ? ? ? 0 ? ? ? 1
B ? ? ? 2 ? ? ? 4 ? ? ? 0 ? ? ? 0 ? ? ? 0
B ? ? ? 2 ? ? ? 5 ? ? ? 1 ? ? ? 0 ? ? ? 0
B ? ? ? 2 ? ? ? 6 ? ? ? 0 ? ? ? 1 ? ? ? 0
B ? ? ? 2 ? ? ? 11 ? ? ?0 ? ? ? 0 ? ? ? 1
With the model
m1<-glmer(Plant~G+K+N+Year+(1|Site), ...)
Or is it better to use a single column for the treatments, like this:
Site ? ?Year ? ?Plant ? Treatment
A ? ? ? 1 ? ? ? 5 ? ? ? C
A ? ? ? 1 ? ? ? 4 ? ? ? G
A ? ? ? 1 ? ? ? 7 ? ? ? K
A ? ? ? 1 ? ? ? 10 ? ? ?N
A ? ? ? 2 ? ? ? 3 ? ? ? C
A ? ? ? 2 ? ? ? 4 ? ? ? G
A ? ? ? 2 ? ? ? 8 ? ? ? K
A ? ? ? 2 ? ? ? 12 ? ? ?N
B ? ? ? 1 ? ? ? 7 ? ? ? C
B ? ? ? 1 ? ? ? 3 ? ? ? G
B ? ? ? 1 ? ? ? 7 ? ? ? K
B ? ? ? 1 ? ? ? 12 ? ? ?N
B ? ? ? 2 ? ? ? 4 ? ? ? C
B ? ? ? 2 ? ? ? 5 ? ? ? G
B ? ? ? 2 ? ? ? 6 ? ? ? K
B ? ? ? 2 ? ? ? 11 ? ? ?N
With the following model:
m1<-glmer(Plants~Treatment+Year+(1|Site), ...)
The latter is preferred. ?R will generate the indicator columns for the levels of the Treatment factor (the 0/1 columns shown in the first form) and, when appropriate, reduce them to a set of 2 "contrasts" in the model. ?(The reason for quoting the word "contrasts" is that there is a formal mathematical definition of a contrast but the linear combinations generated by R do not always satisfy this definition. The method and results are correct, it is just the name that is inaccurate.) The reason that the latter is preferred is that it is easier to maintain the data in a consistent form (factors maintain consistency and are easy to check in the output from str() or summary(), whereas indicator columns have inter-column dependencies that must be checked separately) and the "when appropriate" clause above. ?Determining a useful parameterization of a linear model incorporating factors is subtle and a lot of code in the R function model.matrix is devoted to a symbolic analysis designed to get this right. ?Also, you can, if you wish, change the parameterization (see ?contrasts).
-- Raldo