[R-lang] Re: Main effects of categorical predictors in lmer

hossein karimi karimihussein@gmail.com
Mon Oct 10 12:41:32 PDT 2011


Hi Florian and Ian,

Sorry about the mistake. It caused a bit of confusion. NP is the same as A,
my first 2-level predictor. Here is the output of the contrasts for B:


> contrasts(B)
                 [,1] [,2]
B2               0    1
B1              -1   -1
B3               1    0

The three levels of B are coded with names rather than numbers and the
baseline condition of B alphabetically comes second; that's why I manually
recoded the matrix so that the baseline receives -1. However, as you can
see, this matrix is different from the one Florian posted. I used the
following code to create the contrasts for B (but I just read Ian's message
which says this is wrong):

> contrasts (B)<- contr.sum(3)
> contrasts(B) <- cbind(c(0,-1,1),c(1,-1,0))

So, if I want to compare the baseline condition of B with the mean of other
two conditions of B, I need to create contrasts like Florian's? When I
helmet code the levels of B with the baseline receiving 2 and the other two
levels of B receiving -1, I get the following confusing contrasts. I have no
idea what this means. I’m just hoping that I’m comparing the baseline with
the mean of the other B levels:
> contrasts (helmertB)
     2
-1  0
2   1


And here is the complete output of the first model (where A is numerical and
centered):

Generalized linear mixed model fit by the Laplace approximation
Formula: response ~ A * B + (A + 1 | sub) + (1 | item)
   Data: mydata
   AIC   BIC logLik deviance
 745.1 794.6 -362.6    725.1
Random effects:
 Groups Name        Variance Std.Dev. Corr
 item   (Intercept) 0.93903  0.96903
 sub    (Intercept) 2.39079  1.54622
        A                0.34163  0.58449  0.551
Number of obs: 1038, groups: item, 42; sub, 36

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)
(Intercept)        1.1661     0.3285   3.550 0.000385 ***
A                     -2.0741     0.1690 -12.270  < 2e-16 ***
B1                   -0.4166     0.1580  -2.636 0.008389 **
B2                    -0.2867     0.1783  -1.608 0.107842
A:B1                 0.9881     0.1594   6.199 5.68e-10 ***
A:B2                 -0.7101     0.1750  -4.059 4.94e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
                  (Intr)     NP     lngt.1 lngt.2 NP:l.1
A                0.121
B1              -0.173  0.251
B2              0.003  0.060 -0.281
A:B1           0.141 -0.304 -0.259 -0.095
A:B2           0.021  0.056 -0.081 -0.300 -0.334

 Florian, you already interpreted this output but it was based on the your
assumption of the contrasts matrix for B…


By the way, I didn't quite get Ian's suggestion for A. Should I keep it
numerical and centered?

Many thanks,

hossein


On Mon, Oct 10, 2011 at 6:11 PM, Finlayson, Ian <IFinlayson@qmu.ac.uk>wrote:

> **
>
> You probably are thinking in the ANOVA way (will leave that for others to
> decide if that's good/bad, appropriate/innappropriate) but anova(model) will
> give you an ANOVA-style table for your model (minus the p values).
>
> Not entirely sure what you mean about having added NP, I assume you added
> it as another fixed effect, but it's not listed in your results, so you may
> have meant something else. If B2 becomes significant after including NP then
> that just suggests to me that NP cleans up noise in your data, allowing B2
> to become significant (can't help but notice that the p value in the first
> model is only just over .1, so it's not massively far from significance).
>
> I should have made the point earlier, but I will now: If you're actually
> interested in whether there's a difference between B1 and B2, and between B1
> and B3, then sum coding isn't appropriate. Sum coding will tell you if
> there's a difference between the levels of the factor and the grand mean. If
> you actually want to compare B1 with B2, and B3, then you'll need to use
> treatment coding. Your model suggests this is experimental data (Subject,
> Item etc.) so I'll assume it's a roughly balanced design and collinearity
> might not be a huge concern.
>
> You're correct that you could use subsetting (dividing the data according
> to A makes sense to me, but that is obviously atheoretical as I have no idea
> what A represents) to get a clear idea of your pattern of results. Although,
> as I've just noticed Florian showing, you can get a good picture if you
> understand the coding and work through the coefficients.
>
> HTH,
>
> Ian
>
>
>
> -----Original Message-----
> From: hossein karimi [mailto:karimihussein@gmail.com<karimihussein@gmail.com>
> ]
> Sent: Mon 10/10/2011 16:15
> To: Finlayson, Ian
> Cc: ling-r-lang-l@mailman.ucsd.edu
> Subject: Re: [R-lang] Main effects of categorical predictors in lmer
>
> Hi Ian,
>
> Thank you for your help.I did what you suggested and got the following
> output:
> Fixed effects:
>                  Estimate Std. Error z value Pr(>|z|)
> (Intercept)        1.1661     0.3285   3.550 0.000385 ***
> A                    -2.0741     0.1690 -12.270  < 2e-16 ***
> B1                   -0.4166     0.1580  -2.636 0.008389 **
> B2                   -0.2867     0.1783  -1.608 0.107842
> A:B1                 0.9881     0.1594   6.199 5.68e-10 ***
> A:B2                -0.7101     0.1750  -4.059 4.94e-05 ***
>
> I think this is telling me there is a main effect of A, but my problem is
> that I also need to report a main effect of B, and that's something I don't
> find in this output. I might be thinking in the old ANOVA way, but I think
> I
> should report main effects as well as a main interaction (if there is any)
> and then if there's a significant interaction, I should proceed with
> breaking down the main interaction to see where exactly it lies, am I
> right?
> or, are lmer results reported differently?
>
> By the way, I also ran a model with NP as a factor too and got the
> following
> result:
> Fixed effects:
>                   Estimate Std. Error z value Pr(>|z|)
> (Intercept)         0.9371     0.3313   2.829  0.00468 **
> A1                   2.0857     0.1700  12.269  < 2e-16 ***
> B1                   -0.3075     0.1544  -1.991  0.04646 *
> B2                    -0.3651     0.1735  -2.105  0.03533 *
> A1:B1              -0.9936     0.1603  -6.199 5.68e-10 ***
> A1:B2               0.7139     0.1759   4.058 4.96e-05 ***
>
> The effect of B2 is significant in the first model above but not in
> second!!! why is that?
>
> Many thanks,
>
> Hossein
>
>
>
>
> On Mon, Oct 10, 2011 at 3:27 PM, Finlayson, Ian <IFinlayson@qmu.ac.uk
> >wrote:
>
> > Hi Hossein,****
> >
> > ** **
>
> >
> > As far as I can see, neither model is correct. If B is a predictor with
> > three levels then the model should contain 6 fixed effects: An intercept,
> A,
> > B1 compared to B2, B1 compared to B3, B1 compared to B2 interacting with
> A,
> > and B1 compared to B3 interacting with A.****
> >
> > ** **
>
> >
> > What I suspect you've done is actually enter B as a continuous predictor
> > with a mean of zero (since you used scale()). You'll need to keep B as a
> > factor, but give it appropriate contrasts:****
> >
> > ** **
> >
> > B <- factor(B) # Make sure it's a factor****
>
> >
> > contrasts(B) <- contr.sum(3) # Generate a 3 level contrast matrix using
> sum
> > coding; alternative contr.helmert(3) for Helmhert coding****
> >
> > ** **
> >
> > Then build your models.****
> >
> > ** **
>
> >
> > It's possible that it won't set the right level as a baseline (the
> baseline
> > will be whichever is a row of -1). If so then you'll need to manually set
> up
> > the matrix. For sum coding, something like:****
> >
> > ** **
>
> >
> > contrasts(B) <- cbind(c(-1,0,1),c(-1,1,0)) # Assuming that B1 is the
> first
> > row in contrasts(B) - double check!****
> >
> > ** **
>
> >
> > I'd also do something similar with A so that your estimates are
> comparable
> > (if you've centered 0 and 1, then it'll come out -.5 and .5, rather than
> -1
> > and 1)****
> >
> > ** **
> >
> > HTH,****
> >
> > ** **
> >
> > Ian****
> >
> > ** **
> >
> > ** **
> >
> > ** **
> >
> > ** **
> >
> > *From:* ling-r-lang-l-bounces@mailman.ucsd.edu [mailto:
> > ling-r-lang-l-bounces@mailman.ucsd.edu] *On Behalf Of *hossein karimi
> > *Sent:* 10 October 2011 15:05
> > *To:* ling-r-lang-l@mailman.ucsd.edu
> > *Subject:* [R-lang] Main effects of categorical predictors in lmer****
> >
> > ** **
> >
> > Dear R users,****
>
> >
> > I'm using mixed effects models (lmer) to predict a binary dependent
> > variable as a function of 1.a categorical predictor (A)with 2 levels (A1
> and
> > A2) , 2. another categorical predictor (B) with three levels (B1, B2 and
> B3)
> > and 3. The interaction between these two predictors. I have tried two
> models
> > but they return different results and I'm not sure which one is correct.
> I'm
> > interested in the main effect of B and the interaction between A and B
> > (because A alone has a significant effect in both models). My problem is
> > that there seem to be two sensible ways of examining the main effect of
> B:
> > 1. to helmert code and 2. to center.  But these two methods produce
> opposite
> > results! I don't know which one I should use. Here are the two models
> with
> > some details and their outputs:****
> >
> > ** **
>
> >
> > Model 1: 'A' is centered. 'B' is helmert coded ('B1'(baseline)=2,
> 'B2'=-1,
> > 'B3'=-1) so that I can get a main effect of B by checking to see whether
> > baseline condition in B differs from the mean of B1 and B2 . The lmer
> output
> > returns a significant effect of B and no significant AxB interaction.
> > However, as is highlighted below (in pink), the correlation between B and
> > the 'AxB' interaction is high (-54%).  ****
> >
> >  ****
>
> >
> > > model.1<-lmer(response~A*B+(A+1|sub)+(1|item), mydata,
> family="binomial")
> > ****
> >
> > > print(model.1)****
> >
> > Generalized linear mixed model fit by the Laplace approximation ****
> >
> > Formula: response ~ A * B+ (A + 1 | sub) + (1 | item) ****
> >
> >    Data: mydata****
> >
> >  AIC   BIC logLik deviance****
> >
> > 783 822.6 -383.5      767****
> >
> > Random effects:****
> >
> > Groups Name        Variance Std.Dev. Corr  ****
> >
> >  item   (Intercept) 0.7293   0.85399        ****
> >
> >  sub    (Intercept) 2.0871   1.44468        ****
> >
> >         A          1.3812   1.17524  0.562 ****
> >
> > Number of obs: 1038, groups: item, 42; sub, 36****
> >
> >  ****
> >
> > Fixed effects:****
> >
> >                   Estimate Std. Error z value Pr(>|z|)    ****
> >
> > (Intercept)        1.05261    0.30283   3.476 0.000509 *******
> >
> > A                -3.91080    0.32239 -12.131  < 2e-16 *******
> >
> > B                  0.36128    0.09751   3.705 0.000211 *******
> >
> > A:B            -0.29638    0.18681  -1.586 0.112626    ****
> >
> > ---****
> >
> > Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ****
> >
> >  ****
> >
> > Correlation of Fixed Effects:****
> >
> >                         (Intr)   A      B****
> >
> > A                  0.155              ****
> >
> > B                  0.160 -0.278       ****
> >
> > A:B              -0.156  0.238 -0.540****
> >
> >  ****
>
> >
> > Model 2: 'A' and 'B' are both centered. The lmer output returns no
> > significant effect of B but the A:B interaction is significant. The
> > correlations between predictors are generally lower and the correlation
> > between B and A:B is reduced to -26%. ****
> >
> > ** **
> >
> > Generalized linear mixed model fit by the Laplace approximation ****
> >
> > Formula: resonse ~ A * B + (A + 1 | sub) + (1 | item) ****
> >
> >    Data: mydata ****
> >
> >    AIC   BIC logLik deviance****
> >
> > 756.1 795.7 -370.1    740.1****
> >
> > Random effects:****
> >
> > Groups Name        Variance Std.Dev. Corr  ****
> >
> >  item   (Intercept) 0.87028  0.93289        ****
> >
> >  sub    (Intercept) 2.41707  1.55469        ****
> >
> >         A         1.23669  1.11206  0.533 ****
> >
> > Number of obs: 1038, groups: item, 42; sub, 36****
> >
> >  ****
> >
> > Fixed effects:****
> >
> >                 Estimate Std. Error z value Pr(>|z|)    ****
> >
> > (Intercept)       1.1004     0.3239   3.398 0.000679 *******
> >
> > A               -4.0941     0.3248 -12.605  < 2e-16 *******
> >
> > B                -0.1461     0.1400  -1.043 0.296851    ****
> >
> > A:B             1.7923     0.2818   6.360 2.01e-10 *******
> >
> > ---****
> >
> > Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ****
> >
> >  ****
> >
> > Correlation of Fixed Effects:****
> >
> >                        (Intr)      A        B****
> >
> > A                    0.138              ****
> >
> > B                   -0.148  0.185       ****
> >
> > A:B                0.106 -0.292 -0.265****
> >
> >  ****
>
> >
> > I personally think Model 2 is better but the thing is that I have
> centered
> > a categorical predictor with *three* levels. In my searches in the web, I
>
> > have never seen a three-level predictor to be centered; they were all
> > two-level categorical predictors. ****
>
> >
> > I have used the scale() function to center the predictors (I first
> > converted them to numeric variables and then used the scale () function
> to
> > center them). As I mentioned, my problem is that I don't know how to get
> a
> > main effect of B as well as a *main* A:B interaction. On the one hand, it
>
> > seems logical to compare 'B1' (baseline) with the mean of the other two B
> > conditions to see if the B manipulation has a general effect. On the
> other
> > hand, I hear that one needs to center variables to get a main effect.
> ****
> >
> > ** **
> >
> > I would be grateful of you could please help me****
> >
> > ** **
> >
> > Regards,****
> >
> > ** **
> >
> > Hossein****
> >
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucsd.edu/pipermail/ling-r-lang-l/attachments/20111010/1848bb27/attachment-0001.html 


More information about the ling-r-lang-L mailing list