post hoc comparison following with multcomp?

classic Classic list List threaded Threaded
2 messages Options
Reply | Threaded
Open this post in threaded view
|

post hoc comparison following with multcomp?

Michaël Coeurdassier
Dear R community,

I would like to check differences between treatment (Trait) following a
glm. From R help list, I tried in the following way using the multcom
library. Please, is it a correct manner to do post hoc comparison with
the data below.

Thank you in advance. Sincerely

*> tabp<-read.delim("ponteM1.txt")*
*> tabp*
   Trait totponte
1      T       10
2      T       11
3      T       11
4      T        9
5      T        7
6      T        7
7      T        9
8      T       12
9      T        9
10     T       10
11     M        9
12     M       10
13     M        8
14     M        8
15     M       10
16     M        8
17     M        8
18     M        9
19     M       11
20     M        7
21     F        8
22     F        8
23     F        5
24     F        9
25     F        5
26     F        7
27     F        5
28     F        7
29     F        6
30     F        3
 
*> attach(tabp)*
*> glm1<-glm(totponte~Trait,family=poisson)* **
*> anova(glm1,test="F")*
Model: poisson, link: log
Response: totponte
Terms added sequentially (first to last)
 
      Df Deviance Resid. Df Resid. Dev      F  Pr(>F)
NULL                     29    16.3879                
Trait  2   7.1770        27     9.2109 3.5885 0.02764 *
* *
*> library(multcomp)*
*> (coefglm1<-coef(glm1))*
 (Intercept)      TraitM      TraitT
  1.8405496   0.3342021   0.4107422
 
*> (vc.trait<-vcov(glm1))*
            (Intercept)      TraitM      TraitT
(Intercept)  0.01587301 -0.01587301 -0.01587301
TraitM      -0.01587301  0.02723665  0.01587301
TraitT      -0.01587301  0.01587301  0.02639933
 
*> (CM<-contrMat(table(tabp$Trait),type="Tukey"))*
     F  M T
M-F -1  1 0
T-F -1  0 1
T-M  0 -1 1
* *
*> csimint(coefglm1,df=27,covm=vc.trait,cmatrix=CM)*
Simultaneous confidence intervals: user-defined contrasts
        95 % confidence intervals
 
    Estimate  2.5 % 97.5 %
M-F   -1.506 -2.177 -0.836
T-F   -1.430 -2.097 -0.763
T-M    0.077 -0.286  0.439


---------------
Michaël COEURDASSIER, PhD
Department of Environmental Biology
UsC INRA EA3184MRT
Institute for Environmental Sciences and Technology

University of Franche-Comte
Place Leclerc
25030 Besançon cedex
FRANCE
Tel : +33 (0)381 665 741
Fax : +33 (0)381 665 797

E-m@il: [hidden email]
http://lbe.univ-fcomte.fr/

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Reply | Threaded
Open this post in threaded view
|

Re: post hoc comparison following with multcomp?

Spencer Graves
          Your use of the multcomp package looks plausible to me (though I
certainly did not check every detail).

          Do you have any specific concerns?  [Also, have you reviewed
vignette("Rmc")?]  And have you considered the plotting, e.g.:

          plot(csimint(coefglm1,df=27,covm=vc.trait,cmatrix=CM))

          hope this helps.
          spencer graves

Michaël Coeurdassier wrote:

> Dear R community,
>
> I would like to check differences between treatment (Trait) following a
> glm. From R help list, I tried in the following way using the multcom
> library. Please, is it a correct manner to do post hoc comparison with
> the data below.
>
> Thank you in advance. Sincerely
>
> *> tabp<-read.delim("ponteM1.txt")*
> *> tabp*
>    Trait totponte
> 1      T       10
> 2      T       11
> 3      T       11
> 4      T        9
> 5      T        7
> 6      T        7
> 7      T        9
> 8      T       12
> 9      T        9
> 10     T       10
> 11     M        9
> 12     M       10
> 13     M        8
> 14     M        8
> 15     M       10
> 16     M        8
> 17     M        8
> 18     M        9
> 19     M       11
> 20     M        7
> 21     F        8
> 22     F        8
> 23     F        5
> 24     F        9
> 25     F        5
> 26     F        7
> 27     F        5
> 28     F        7
> 29     F        6
> 30     F        3
>  
> *> attach(tabp)*
> *> glm1<-glm(totponte~Trait,family=poisson)* **
> *> anova(glm1,test="F")*
> Model: poisson, link: log
> Response: totponte
> Terms added sequentially (first to last)
>  
>       Df Deviance Resid. Df Resid. Dev      F  Pr(>F)
> NULL                     29    16.3879                
> Trait  2   7.1770        27     9.2109 3.5885 0.02764 *
> * *
> *> library(multcomp)*
> *> (coefglm1<-coef(glm1))*
>  (Intercept)      TraitM      TraitT
>   1.8405496   0.3342021   0.4107422
>  
> *> (vc.trait<-vcov(glm1))*
>             (Intercept)      TraitM      TraitT
> (Intercept)  0.01587301 -0.01587301 -0.01587301
> TraitM      -0.01587301  0.02723665  0.01587301
> TraitT      -0.01587301  0.01587301  0.02639933
>  
> *> (CM<-contrMat(table(tabp$Trait),type="Tukey"))*
>      F  M T
> M-F -1  1 0
> T-F -1  0 1
> T-M  0 -1 1
> * *
> *> csimint(coefglm1,df=27,covm=vc.trait,cmatrix=CM)*
> Simultaneous confidence intervals: user-defined contrasts
>         95 % confidence intervals
>  
>     Estimate  2.5 % 97.5 %
> M-F   -1.506 -2.177 -0.836
> T-F   -1.430 -2.097 -0.763
> T-M    0.077 -0.286  0.439
>
>
> ---------------
> Michaël COEURDASSIER, PhD
> Department of Environmental Biology
> UsC INRA EA3184MRT
> Institute for Environmental Sciences and Technology
>
> University of Franche-Comte
> Place Leclerc
> 25030 Besançon cedex
> FRANCE
> Tel : +33 (0)381 665 741
> Fax : +33 (0)381 665 797
>
> E-m@il: [hidden email]
> http://lbe.univ-fcomte.fr/
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html