[R-sig-ME] multiple comparison for a LME with interaction

Kay Cecil Cichini Kay.Cichini at uibk.ac.at
Thu Aug 4 18:04:16 CEST 2011


hi,

you will need to provide example data, otherwise it will be hard for  
anyone to offer help.

in any case you could try if the contrasts are set up properly by  
releveling your factors (relevel-function) using the default treatment  
contrasts, then compare with your custom contrasts and see if the  
coefficients are the same between the two models where they should be...

interactionplot() is also nice for checking if your coefficients make sense.

best,
kay



Zitat von Poulin Nicolas <poulin at math.unistra.fr>:

> Dear all,
> I would like to know if you think that what I am doing is correct. I
> have always the same self-questionning when working with (generalized)
> linear mixed models with interactions.
>
> First let me introduce you my data.
>   A group of monkeys have been experimented within 3 different
> conditions (variable Nb with values 1, 2 and 3).
> As a social group of monkeys, the rank among the group of each monkey is
> important in the way of explaining a reaction (variable Rk).
> I want to see the effect of the condition, the rank and the interaction
> of these two variables on the outcome of the experiment (variable Hz).
>
> There are 540 observations among 15 individuals (36 per individual)
>
> The selected model :
> - the variable Nb have been transformed into a factor
> - the square root of the variable Hz was used to be able to use a LME
> instead of a GLMM (valided by a Q-Q plot)
> - Hz seemed to be lineary dependent of the rank so I didn't used a
> factor version of Rk (hence Rk is numeric)
>
> I only tried the random intercept based on the identity (ID) of each
> individual.
>
> So the model I am using is :
> model.sqrt<-lme(sqrtHz~Rk*Nb, random = ~ 1 | ID,data=data1)
>
> Here is the summary of the fixed part of the model :
> Fixed effects: sqrtHz ~ Rk * Nb
>
>                     Value                Std.Error       DF
> t-value         p-value
> (Intercept)  0.5031499      0.04741670    521     10.611238     0.0000
> Rk             -0.0456437      0.00641704     13       -7.112890    0.0000
> Nb2            0.0551886      0.02609983    521       2.114520     0.0349
> Nb3            0.0026095      0.02609983    521      0.099982      0.9204
> Rk:Nb2      -0.0004107     0.00353217    521     -0.116277     0.9075
> Rk:Nb3      0.0144365      0.00353217     521      4.087143     0.0001
>
> To obtain a p-value for the factor Nb and the intercation, I use :
>  > anova(model.sqrt)
>                        numDF denDF      F-value         p-value
> (Intercept)        1          521       113.85619 <.0001
> Rk                    1            13         45.33842 <.0001
> Nb                     2          521        18.49245 <.0001
> Rk:Nb                2         521        11.46233 <.0001
>
>
> Hence both Nb and the interaction have significant effects.
>
> My first question is do I have the right to do such a thing?
>
> I want to compare the different experiment and the interactions.
>
> First, I tried to use summary(glht(model.sqrt,linfct=mcp(Nb="Tukey")))
> that give me the result :
>
>                   Estimate    Std. Error z value   Pr(>|z|)
> 2 - 1 == 0  0.05519    0.02610     2.115      0.0868 .
> 3 - 1 == 0  0.00261    0.02610     0.100      0.9945
> 3 - 2 == 0 -0.05258    0.02610     -2.015     0.1088
>
> This output came with the warning message : covariate interactions found
> -- default contrast might be inappropriate
>
> Moreover I can not apply this to have multiple comparisons among the
> interactions.
>
> Browsing the mailing list, I tried to follow the method given by :
> http://thebiobucket.blogspot.com/2011/06/glmm-with-custom-multiple-comparisons.html
>
> Hence I did :
>  > c <- rbind("Condition 1 - Condition 2" = c(0,0,1,0,0,0),
> +             "Condition 1 - Condition 3" = c(0,0,0,1,0,0),
> +             "Condition 2 - Condition 3" = c(0,0,1,-1,0,0),
> +             "Rank:Condition 1 - Rank:Condition 2" = c(0,0,0,0,1,0),
> +              "Rank:Condition 1 - Rank:Condition 2" = c(0,0,0,0,0,1),
> +              "Rank:Condition 1 - Rank:Condition 2" = c(0,0,0,0,1,-1)
> + )
>
> and
>
>  > summary(glht(model.sqrt,c))
>
>           Simultaneous Tests for General Linear Hypotheses
>
> Fit: lme.formula(fixed = sqrtHz ~ Rk * Nb, data = data1, random = ~1 |
>      ID)
>
> Linear Hypotheses:
>
>    Estimate     Std. Error     z value   Pr(>|z|)
> Condition 1 - Condition 2 == 0                     0.0551886
> 0.0260998      2.115       0.131
> Condition 1 - Condition 3 == 0                     0.0026095
> 0.0260998       0.100      1.000
> Condition 2 - Condition 3 == 0                     0.0525791
> 0.0260998       2.015       0.162
> Rank:Condition 1 - Rank:Condition 2 == 0 -0.0004107  0.0035322
> -0.116      1.000
> Rank:Condition 1 - Rank:Condition 2 == 0  0.0144365  0.0035322
> 4.087 <0.001 ***
> Rank:Condition 1 - Rank:Condition 2 == 0 -0.0148472  0.0035322
> -4.203 <0.001 ***
>
> So my real question is : is this correct?
> If so, it's quite in contradiction with the results of
> summary(model.sqrt) and anova(model.sqrt), for the multiple comparisons
> among the different conditions. I know that the p-values of
> summary(model.sqrt) can be seen as Dunnett comparison with Condition 1.
>
> Here I don't know what is the correction that is used but I think it
> takes into consideration the fact that I have done 6 different tests.
> These corrections seems to be too conservative as there is no
> significant differences among the different conditions.
> Is there a way to select a less conservative correction?
>
> If not I don't know how to explain the differences between the results
> of  both summary(model.sqrt) and anova(model.sqrt) and the last
> displayed one.
>
> Thanks in advance for your answers, and appologies to have written such
> a long message.
> Best,
> Nicolas
> --
> Nicolas Poulin
> Ingénieur de Recherche
> Centre de Statistique de Strasbourg (CeStatS)
> http://www-math.u-strasbg.fr/CeStatS/
> <http://www-math.u-strasbg.fr/CentreStat/>
> IRMA, UMR 7501
> Université de Strasbourg et CNRS
> 7 rue René-Descartes
> 67084 Strasbourg Cedex
>
> Tél. 33 (0)3 68 85 01 89
>
> 	[[alternative HTML version deleted]]
>
>




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