[R-sig-ME] Interpreting lmer() interactions with Helmert contrasts
Becky Gilbert
beckyannegilbert at gmail.com
Mon Aug 24 13:54:01 CEST 2015
Thanks very much everyone for the responses.
@Dan: Thank you for the recommendation about my factor contrast
coefficients. I hadn't given much thought to the sign/level association,
but now that you point it out, it seems obvious that I should do it the way
you describe. Here are the model coefficients with recoded contrasts:
> contrasts(rtData$Time)
[,1]
-1 -0.5 # pre-test
1 0.5 # post-test
> contrasts(rtData$WordType)
[,1] [,2]
0 -0.6666667 0.0 # untrained
1 0.3333333 0.5 # trained-related
2 0.3333333 -0.5 # trained-unrelated
Estimate Std. Error t value
(Intercept) 2.8765116 0.0177527 162.03
WordType1 -0.0111628 0.0110852 -1.01
WordType2 -0.0007306 0.0071519 -0.10
Time1 0.0268310 0.0195248 1.37
WordType1:Time1 0.0301627 0.0115349 2.61
WordType2:Time1 -0.0089123 0.0141624 -0.63
My interpretations of the interaction coefficients are:
1) log RT increases (i.e. RTs slow down) for the two trained (vs untrained)
Word Types at post-test (Time = 1)
2) log RT decreases (i.e. RTs speed up) for the trained-related (vs
trained-unrelated) Word Type at post-test (Time = 1)..
However, this doesn't really answer my original question about how to
assess (and report) the contribution of these two interactions to the model
fit. Obviously the t statistic is larger for the Time1:WordType1 compared
to the Time1:WordType2 interaction coefficients, but that only tells me
their relative contributions - I would need to know degrees of freedom to
get p-values, which I understand is not straightforward. Also, I've read
that the t statistics for coefficients that are output by summary() for an
lmer model are sequential tests and thus not the appropriate/desired
statistics for assessing the contribution of factors (someone please
correct me if I'm wrong!). Hence the reason for using LRT to assess this.
This still leaves me with the problem of not being able to test the
interactions between Time and the two contrasts for WordType - I can test
the whole WordType factor and Time:WordType interaction via LRTs, but not
each contrast within WordType.
@Steven: thanks for your explanation re interpreting main effects in the
presence of an interaction, and of the Chi-square LRTs for assessing the
contribution of factors/terms.
However I'm confused by this:
An omnibus test for the statistical significance of a variable of interest
> (say variable A), when that variable is in a model involving an interaction
> with another variable (say variable B) will test the interaction term A:B
> and the main effect A. The full model has A + B + A:B and the reduced
> model has only B. Thus a proper omnibus test for the usefulness of A in
> the model will involve the interaction A:B and the main effect A. This
> test really should be done before testing A:B for proper multiple
> comparisons control.
Is this what you're saying?
1. test A: (A + B + A:B) vs (B)
2. test B: (A + B + A:B) vs (A)
then, if either of the above are significant:
3. test A:B: (A + B + A:B) vs (A + B)
Which I think is the procedure described here:
https://mailman.ucsd.edu/pipermail/ling-r-lang-l/2011-October/000305.html
Assuming this is what you meant, will this procedure always get you to step
3 (assessing the interaction) in the case of a significant interaction
without main effects (as in a cross-over interaction). Sorry if I've
completely misunderstood!
Becky
____________________________________________
Dr Becky Gilbert (nee Prince)
http://www.york.ac.uk/psychology/staff/postgrads/becky.gilbert/
http://www.researchgate.net/profile/Becky_Gilbert2
http://twitter.com/BeckyAGilbert
On 22 August 2015 at 02:53, Steven McKinney <smckinney at bccrc.ca> wrote:
>
>
> > -----Original Message-----
> > From: R-sig-mixed-models [mailto:
> r-sig-mixed-models-bounces at r-project.org]
> > On Behalf Of Ken Beath
> > Sent: August-21-15 5:19 PM
> > To: Dan McCloy
> > Cc: r-sig-mixed-models at r-project.org
> > Subject: Re: [R-sig-ME] Interpreting lmer() interactions with Helmert
> > contrasts
> >
> > The first thing to do is to decide if an interaction is really needed.
> > Applying anova() to the model with interaction should give this. Note
> that
> > if there is an interaction it is not possible to conclude anything about
> > the main effects.
>
> In the case of a significant interaction, the conclusion is that the
> change in one of the main effects depends on the level of the other main
> effect. Both main effects are statistically significant, but a single
> value for each main effect will not adequately summarize the degree of
> response. One can make conclusions about the main effects, but the
> conclusions are more complex than a single estimate for each main effect.
>
>
>
> > The tests that you have performed use Chi-sq statistics
> > which is not correct for a linear model, it should be F.
>
> For linear mixed effects models, F tests are not possible. The Chi-square
> tests presented are likelihood ratio chi-square tests comparing the
> likelihoods of the full model (under the alternative hypothesis) and the
> reduced model (under the null hypothesis). The anova() method for lmer
> objects automagically refits the two models using maximum likelihood (as
> opposed to restricted maximum likelihood) and performs the likelihood ratio
> chi-square test. The degrees of freedom of that test equal the number of
> parameters posited to be zero under the reduced (null hypothesis) model as
> compared to the full (alternative hypothesis) model.
>
>
> > Playing around
> > with lme4 I found that wasn't possible except when using a single model.
> > Testing for removal of combined interaction and main effect is not
> correct.
>
> An omnibus test for the statistical significance of a variable of interest
> (say variable A), when that variable is in a model involving an interaction
> with another variable (say variable B) will test the interaction term A:B
> and the main effect A. The full model has A + B + A:B and the reduced
> model has only B. Thus a proper omnibus test for the usefulness of A in
> the model will involve the interaction A:B and the main effect A. This
> test really should be done before testing A:B for proper multiple
> comparisons control.
>
>
>
> Steven McKinney, Ph.D.
>
> Statistician
> Molecular Oncology and Breast Cancer Program
> British Columbia Cancer Research Centre
>
>
>
>
>
>
> >
> > If there is an interaction I would suggest using the standard contrasts.
> > Then you will need to summarise this and there are a few ways in R to do
> > this. One is svycontrast in the survey package. The result should be a
> > different effect of time for each word type or vice-versa. This is
> > something that should be talked over with a statistician.
> >
> > On 22 August 2015 at 03:39, Dan McCloy <drmccloy at uw.edu> wrote:
> >
> > > As a word of caution, you seem to have set up your factor coding to
> make
> > > interpretation especially tricky. The coding of your "Time1" variable
> is
> > > set up so that your factor level of "-1" has a positive coefficient,
> and
> > > your factor level of "1" has a negative coefficient. Before doing
> > anything
> > > else, I recommend you re-run the model after re-setting the contrasts
> for
> > > "Time" so that your textual levels have the same sign as their
> > coefficients
> > > in the model (personally I would go further and re-code the factor as
> > "Pos"
> > > and "Neg" or some other textual shorthand that cannot be confused with
> > row
> > > or column numbers of the contrast matrix). I also usually set the row
> > > names of contrast matrices to be actual names, so that the lmer output
> > > names the coefficients in a way that is harder for me to mis-interpret
> > > (e.g., as "TimePos" or "TimeNeg" instead of "Time1"). While you're at
> > it,
> > > if you're interested in "treatment" vs "no treatment" you might
> consider
> > > re-setting the contrasts for the WordType factor as well. You have
> this:
> > >
> > > [,1] [,2]
> > > 0 0.6666667 0.0 # untrained
> > > 1 -0.3333333 -0.5 # trained-related
> > > 2 -0.3333333 0.5 # trained-unrelated
> > >
> > > which means that *positive* coefficient estimates for factor 1 mean
> that
> > > "untrained" increases RT. Similar comment for related vs. unrelated.
> I
> > > would recommend swapping the signs on both factors so that anything
> that
> > is
> > > "un-" is negative, like this:
> > >
> > > [,1] [,2]
> > > 0 -0.6666667 0.0 # untrained
> > > 1 0.3333333 0.5 # trained-related
> > > 2 0.3333333 -0.5 # trained-unrelated
> > >
> > > As far as interpreting the model coefficients for the interactions:
> > >
> > > WordType1:Time1 0.0301627 0.0115349 2.61
> > > WordType2:Time1 -0.0089123 0.0141624 -0.63
> > >
> > > This says that comparing cases of "WordType1" (which curently means
> > > "untrained minus trained" in your experiment) combined with "Time1"
> > (which
> > > I think means Time=1 or what I'm calling "Pos") has a positive
> > coefficient
> > > (the combination increases log reaction time, or slows people down)
> > > relative to what you would expect if "WordType" and "Time" contributed
> > > independently to reaction time. In other words, I think this means
> that
> > > lack of training slows people down more when Time=1 than when Time=-1
> > > (though the mismatch between signs of the factor levels and contrast
> > > coefficients for the Time variable make me hesitate as to whether I
> said
> > > that last bit backwards).
> > >
> > > Hope it helps, and good luck.
> > > -- dan
> > >
> > > Daniel McCloy
> > > http://dan.mccloy.info/
> > > Postdoctoral Research Fellow
> > > Institute for Learning and Brain Sciences
> > > University of Washington
> > >
> > >
> > > On Fri, Aug 21, 2015 at 6:23 AM, Becky Gilbert
> > <beckyannegilbert at gmail.com
> > > >
> > > wrote:
> > >
> > > > Hi Paul/List,
> > > >
> > > > After thinking about this a bit more, I don't think planned
> comparisons
> > > > gives me what I'm looking for. I want to know whether the effect of
> > Time
> > > > is different for WordType = 0 (level 1) vs the other two levels
> > combined
> > > -
> > > > this is WordType contrast 1. I also want to know whether the effect
> of
> > > > Time is different for WordType = 1 vs WordType = 2 (i.e. level 2 vs
> > level
> > > > 3) - this is WordType contrast 2.
> > > >
> > > > I think the use of planned comparisons here would defeat the purpose
> of
> > > my
> > > > contrasts, but maybe I'm missing something?
> > > >
> > > > Thanks!
> > > > Becky
> > > >
> > > > ____________________________________________
> > > >
> > > > Dr Becky Gilbert
> > > >
> > > > On 21 August 2015 at 13:46, Becky Gilbert <
> beckyannegilbert at gmail.com>
> > > > wrote:
> > > >
> > > > > Hi Paul,
> > > > >
> > > > > Thanks very much for the suggestion! I tried using lsmeans() to
> get
> > > the
> > > > > pairwise comparisons as you suggested, and the results are below.
> > > > >
> > > > > I'm a little confused by the results because the pairwise
> comparison
> > > > tests
> > > > > all show p > .05, but the WordType x Time interaction was
> significant
> > > > when
> > > > > tested via model comparisons...? I think this might be due to the
> > > Tukey
> > > > > adjustment for multiple comparisons, but I'm not sure.
> Specifically
> > > the
> > > > > contrast for the two levels of Time at WordType = 2 looks like it
> > might
> > > > > have been significant before the multiple comparisons correction,
> > thus
> > > > > accounting for the significance of the interaction term in model
> > > > > comparisons. Any thoughts?
> > > > >
> > > > > Thanks again!
> > > > > Becky
> > > > >
> > > > > $lsmeans
> > > > > WordType = 0:
> > > > > Time lsmean SE df lower.CL upper.CL
> > > > > -1 2.880592 0.02209390 21.58 2.834721 2.926464
> > > > > 1 2.887315 0.02144245 22.13 2.842860 2.931769
> > > > >
> > > > > WordType = 1:
> > > > > Time lsmean SE df lower.CL upper.CL
> > > > > -1 2.856211 0.02156603 19.78 2.811193 2.901229
> > > > > 1 2.888640 0.02089339 20.17 2.845080 2.932200
> > > > >
> > > > > WordType = 2:
> > > > > Time lsmean SE df lower.CL upper.CL
> > > > > -1 2.852485 0.02181905 20.72 2.807072 2.897898
> > > > > 1 2.893827 0.02113775 21.12 2.849883 2.937770
> > > > >
> > > > > Confidence level used: 0.95
> > > > >
> > > > > $contrasts
> > > > > WordType = 0:
> > > > > contrast estimate SE df t.ratio p.value
> > > > > -1 - 1 -0.00672255 0.02078469 19.31 -0.323 0.7498
> > > > >
> > > > > WordType = 1:
> > > > > contrast estimate SE df t.ratio p.value
> > > > > -1 - 1 -0.03242907 0.02097452 20.02 -1.546 0.1377
> > > > >
> > > > > WordType = 2:
> > > > > contrast estimate SE df t.ratio p.value
> > > > > -1 - 1 -0.04134141 0.02146707 21.93 -1.926 0.0672
> > > > >
> > > > >
> > > > > ____________________________________________
> > > > >
> > > > > Dr Becky Gilbert
> > > > >
> > > > > On 21 August 2015 at 12:19, paul debes <paul.debes at utu.fi> wrote:
> > > > >
> > > > >> Hi Becky,
> > > > >>
> > > > >> Maybe you are interested in pairwise comparisons? The "lsmeans"
> > > package
> > > > >> comes in handy.
> > > > >>
> > > > >> Try something like this:
> > > > >>
> > > > >> library("pbkrtest") # gives you KW-adjusted denDF for tests, but
> > must
> > > be
> > > > >> installed
> > > > >> library("lsmeans")
> > > > >>
> > > > >> Model.lmer.means = lsmeans(Model, spec = pairwise ~ WordType|Time)
> > > > >> Model.lmer.means = summary(Model.lmer.means)
> > > > >> Model.lmer.means
> > > > >>
> > > > >> Maybe you want the contrast conditional on WordType, not Time?
> Swap
> > it
> > > > to:
> > > > >> "spec = pairwise ~ Time|WordType"
> > > > >>
> > > > >> Best,
> > > > >> Paul
> > > > >>
> > > > >>
> > > > >> On Fri, 21 Aug 2015 14:04:07 +0300, Becky Gilbert <
> > > > >> beckyannegilbert at gmail.com> wrote:
> > > > >>
> > > > >> Dear list,
> > > > >>>
> > > > >>> I'm wondering if someone could help me interpret an interaction
> > > between
> > > > >>> two
> > > > >>> factors, when one of the factors uses Helmert contrasts?
> > > > >>>
> > > > >>> I ran a linear mixed effects model (lmer) with reaction times as
> > the
> > > > DV,
> > > > >>> 2
> > > > >>> fixed factors: Time (2 levels) and Word Type (3 levels), and 2
> > random
> > > > >>> factors: Subjects and Items. I used Helmert contrasts for the
> Word
> > > > Type
> > > > >>> factor:
> > > > >>> - Contrast 1 = level 1 (Untrained) vs levels 2 & 3
> (Trained-related
> > > and
> > > > >>> Trained-unrelated)
> > > > >>> - Contrast 2 = level 2 vs. level 3 (Trained-related vs
> > > > Trained-unrelated)
> > > > >>> The data, contrasts, model, summary and model comparisons are
> > listed
> > > at
> > > > >>> the
> > > > >>> end of the message.
> > > > >>>
> > > > >>> Model comparisons with anova() showed a significant interaction
> > > between
> > > > >>> Time and Word Type. However, I don't know how to get the
> > statistics
> > > > for
> > > > >>> the interactions between Time and each Word Type contrast.
> > > > >>>
> > > > >>> Based on the t-values for coefficients in the model summary, it
> > looks
> > > > >>> like
> > > > >>> the significant Word Type x Time interaction is driven by the
> > > > interaction
> > > > >>> with the 1st contrast for Word Type (t = 2.61). However I don't
> > > think
> > > > >>> that
> > > > >>> the statistics for the fixed effects coefficients are exactly
> what
> > > I'm
> > > > >>> looking forward (they are sequential tests, right?). And if
> these
> > > are
> > > > >>> the
> > > > >>> appropriate statistics, I'm aware of the problems with trying to
> > get
> > > > >>> p-values from these estimates. So is there a way to do
> likelihood
> > > > ratio
> > > > >>> tests for each Word Type contrast, or some other way of
> > interpreting
> > > > the
> > > > >>> Word Type x Time interaction?
> > > > >>>
> > > > >>> Data structure:
> > > > >>>
> > > > >>>> str(rtData)
> > > > >>>>
> > > > >>> 'data.frame': 1244 obs. of 11 variables:
> > > > >>> $ Subject : Factor w/ 16 levels "AB","AS","AW",..: 1 1 1 1
> 1
> > 1
> > > 1
> > > > 1
> > > > >>> 1
> > > > >>> 1 ...
> > > > >>> $ Item : Factor w/ 48 levels "ANT","BANDAGE",..: 3 4 6 12
> 13
> > > 14
> > > > 22
> > > > >>> 29 30 34 ...
> > > > >>> $ Response : int 960 1255 651 1043 671 643 743 695 965 589
> > ...
> > > > >>> $ Time : Factor w/ 2 levels "-1","1": 1 1 1 1 1 1 1 1 1 1
> ...
> > > > >>> $ WordType : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1
> 1 1
> > 1
> > > > ...
> > > > >>> $ logRT : num 2.98 3.1 2.81 3.02 2.83 ...
> > > > >>>
> > > > >>> contrasts(rtData$Time)
> > > > >>>>
> > > > >>> [,1]
> > > > >>> -1 0.5
> > > > >>> 1 -0.5
> > > > >>>
> > > > >>> contrasts(rtData$WordType)
> > > > >>>>
> > > > >>> [,1] [,2]
> > > > >>> 0 0.6666667 0.0
> > > > >>> 1 -0.3333333 -0.5
> > > > >>> 2 -0.3333333 0.5
> > > > >>>
> > > > >>> Model:
> > > > >>> lmer(logRT ~ 1 + WordType + Time + WordType:Time +
> > > > >>> (1 + Time|Subject) +
> > > > >>> (1|Item),
> > > > >>> data = rtData)
> > > > >>>
> > > > >>> REML criterion at convergence: -2061.2
> > > > >>> Scaled residuals:
> > > > >>> Min 1Q Median 3Q Max
> > > > >>> -2.7228 -0.6588 -0.0872 0.5712 3.7790
> > > > >>> Random effects:
> > > > >>> Groups Name Variance Std.Dev. Corr
> > > > >>> Item (Intercept) 0.000933 0.03054
> > > > >>> Subject (Intercept) 0.004590 0.06775
> > > > >>> Time1 0.005591 0.07478 0.05
> > > > >>> Residual 0.009575 0.09785
> > > > >>> Number of obs: 1244, groups: Target, 46; Subject, 16
> > > > >>> Fixed effects:
> > > > >>> Estimate Std. Error t value
> > > > >>> (Intercept) 2.8765116 0.0177527 162.03
> > > > >>> WordType1 0.0111628 0.0110852 1.01
> > > > >>> WordType2 0.0007306 0.0071519 0.10
> > > > >>> Time1 -0.0268310 0.0195248 -1.37
> > > > >>> WordType1:Time1 0.0301627 0.0115349 2.61
> > > > >>> WordType2:Time1 -0.0089123 0.0141624 -0.63
> > > > >>>
> > > > >>> Model comparisons with anova() for main effects and interaction:
> > > > >>>
> > > > >>> -full model vs no Word Type x Time interaction
> > > > >>> Df AIC BIC logLik deviance
> > > > Chisq
> > > > >>> Chi Df Pr(>Chisq)
> > > > >>> rtModelNoInteraction 9 -2077.5 -2031.3 1047.7 -2095.5
> > > > >>>
> > > > >>> rtModelFull 11 -2080.5 -2024.1 1051.2 -2102.5
> 7.0388
> > > > >>> 2
> > > > >>> 0.02962 *
> > > > >>>
> > > > >>> -full model vs model without Time and interaction
> > > > >>> Df AIC BIC logLik deviance Chisq
> > Chi
> > > Df
> > > > >>> Pr(>Chisq)
> > > > >>> rtModelNoTime 8 -2077.8 -2036.7 1046.9 -2093.8
> > > > >>> rtModelFull 11 -2080.5 -2024.1 1051.2 -2102.5 8.7424 3
> > > > >>> 0.03292 *
> > > > >>>
> > > > >>> -full model vs model without Word Type and interaction
> > > > >>> Df AIC BIC logLik deviance Chisq
> Chi
> > > Df
> > > > >>> Pr(>Chisq)
> > > > >>> rtModelNoWT 7 -2080.4 -2044.5 1047.2 -2094.4
> > > > >>> rtModelFull 11 -2080.5 -2024.1 1051.2 -2102.5 8.0875 4
> > > > >>> 0.08842
> > > > >>> .
> > > > >>>
> > > > >>> Thanks in advance for any advice!
> > > > >>> Becky
> > > > >>> ____________________________________________
> > > > >>>
> > > > >>> Dr Becky Gilbert
> > > > >>>
> > > > >>> [[alternative HTML version deleted]]
> > > > >>>
> > > > >>> _______________________________________________
> > > > >>> R-sig-mixed-models at r-project.org mailing list
> > > > >>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> > > > >>>
> > > > >>
> > > > >>
> > > > >> --
> > > > >> Paul Debes
> > > > >> DFG Research Fellow
> > > > >> University of Turku
> > > > >> Department of Biology
> > > > >> Itäinen Pitkäkatu 4
> > > > >> 20520 Turku
> > > > >> Finland
> > > > >>
> > > > >
> > > > >
> > > >
> > > > [[alternative HTML version deleted]]
> > > >
> > > > _______________________________________________
> > > > R-sig-mixed-models at r-project.org mailing list
> > > > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> > > >
> > >
> > > [[alternative HTML version deleted]]
> > >
> > > _______________________________________________
> > > R-sig-mixed-models at r-project.org mailing list
> > > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> > >
> >
> >
> >
> > --
> >
> > *Ken Beath*
> > Lecturer
> > Statistics Department
> > MACQUARIE UNIVERSITY NSW 2109, Australia
> >
> > Phone: +61 (0)2 9850 8516
> >
> > Level 2, AHH
> > http://stat.mq.edu.au/our_staff/staff_-_alphabetical/staff/beath,_ken/
> >
> > CRICOS Provider No 00002J
> > This message is intended for the addressee named and may...{{dropped:9}}
> >
> > _______________________________________________
> > R-sig-mixed-models at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list