[R-sig-ME] post-hoc comparison on interaction term in lme, using contrasts

Lenth, Russell V russell-lenth at uiowa.edu
Wed Jan 15 15:30:35 CET 2014


The lsmeans package might be helpful here.

You can visualize the predictions via an interaction-plot:

    library(lsmeans)
    lsmip(MY_MODEL, treat ~ species)

It looks like you want the treatment comparisons A-D, B-D, C-D for each
species -- is that right? If so, it can be done using

    lsmeans(MY_MODEL, trt.vs.ctrl ~ treat | species, ref = 4)

Russ

Russell V. Lenth  -  Professor Emeritus
Department of Statistics and Actuarial Science   
The University of Iowa  -  Iowa City, IA 52242  USA   
Voice (319)335-0712 (Dept. office)  -  FAX (319)335-3017

On 1/15/2014 5:02 AM, r-sig-mixed-models-request at r-project.org wrote:
> Dear R experts,
>
> I have a significant interaction term in my lme model and I was searching for a way, how to perform post hoc test.
> I was told that direct Tukey test of interaction using glht within lme gives unreliable results and should be avoided.
> I have searched for solutions and found out only a recomendation to build up a contrast matrix using function contrast.
>
> I have 7 plant species and 4 treatments which significantly interact and I used folowing syntax to build the matrix
>
>> > cm<-contrast(MY MODEL, a=list(species=c("A","B","C","D","E","F","G"), treat=c("A","B","C")),
> b=list( species=c("A","B","C","D","E","F","G"), treat=c("D","D","D")))
>
>> > cmtrx <- cm$X
>> > ttgl<-glht( MY MODEL,lin=cmtrx)
>> > confint(ttgl)
> the outcome is 21 rows numbered from 1 to 21 (corresponds to 7 species * 3 treatments), and it is not clear to me, how the combinations are ordered?
> i.e. 1 == 0  is for species A : treat A against species A : treat D
> 2 == 0   is for species A : treat B against species A : treat D
>
> the outcome is estimated values plus confidence intervals such as
>                Estimate       lwr            upr       
> 1 == 0  -0.2935212 -0.4847410 -0.1023014
> 2 == 0  -0.4448065 -0.6360263 -0.2535867
> ..........
>
> I suppose this needs to be further digested by some function from multcomp package to estimate significance of these tests.
> I tried to extend ttgl to ttgl<-glht( MY MODEL,lin=cmtrx(tension = "Tukey")) but that does not work.
> summary(ttgl)$test$pvalues also has not yield what I need.
> Can you please advise me how to go on?
>
> Or is there easier way how to deal with interactions in lme?
>
> Thank you very much for any help.
> Best from Jana




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