Adam Coble
2012-Jul-07 22:51 UTC
[R] Questions about glht() and interpretation of output from Tukey's in multcomp
Hi, I have a few questions about glht() and the interpretation of output from Tukey's in multcomp package for lme() model. The main issue is that I noticed that a plot that I produced with code letters seem to contradict the graph itself. I provide data and code below. I end with my questions. A few things about data set. "LMA.vcp" is continuous response variable. "Canopy.position.vcp" is independent categorical variable and "Leaf.height.vcp" is independent covariate. "plots.vcp" and "tree.vcp" are random effects.> data.vcp ##data-set##LMA.vcp Canopy.position.vcp Leaf.height.vcp plots.vcp tree.vcp 1 96.8 BOTTOM 14.5 2 1 2 104.4 BOTTOM 14.7 2 2 3 95.0 BOTTOM 14.7 2 3 4 105.5 BOTTOM 17.4 4 1 5 105.1 BOTTOM 17.6 4 2 6 94.2 BOTTOM 12.6 11 1 7 85.1 BOTTOM 13.0 11 2 8 101.6 BOTTOM 13.5 11 3 9 95.8 BOTTOM 14.0 12 4 10 91.8 BOTTOM 14.7 12 5 11 101.8 BOTTOM 16.7 17 5 12 108.7 BOTTOM 16.9 17 6 13 83.5 BOTTOM 14.5 19 6 14 86.7 BOTTOM 14.7 19 7 15 88.5 BOTTOM 13.1 21 4 16 103.9 BOTTOM 13.0 21 5 17 98.2 BOTTOM 13.7 21 6 18 121.5 BOTTOM 17.1 23 3 19 112.6 BOTTOM 16.7 23 4 20 93.0 BOTTOM 13.1 28 7 21 84.6 BOTTOM 12.1 28 8 22 86.6 BOTTOM 14.3 29 8 23 108.3 BOTTOM 17.0 31 7 24 78.1 BOTTOM 17.4 31 8 25 100.0 MIDDLE 16.7 2 1 26 99.8 MIDDLE 16.6 2 2 27 82.8 MIDDLE 16.7 2 3 28 103.2 MIDDLE 19.4 4 1 29 104.7 MIDDLE 19.9 4 2 30 80.7 MIDDLE 14.9 11 1 31 87.8 MIDDLE 15.0 11 2 32 88.1 MIDDLE 15.6 11 3 33 85.0 MIDDLE 16.4 12 4 34 98.9 MIDDLE 17.0 12 5 35 69.4 MIDDLE 19.1 17 5 36 106.6 MIDDLE 18.9 17 6 37 96.8 MIDDLE 16.6 19 6 38 96.3 MIDDLE 16.7 19 7 39 83.3 MIDDLE 15.2 21 4 40 96.1 MIDDLE 15.3 21 5 41 91.2 MIDDLE 15.6 21 6 42 118.9 MIDDLE 19.2 23 3 43 111.6 MIDDLE 18.6 23 4 44 88.6 MIDDLE 15.2 28 7 45 96.0 MIDDLE 14.4 28 8 46 92.2 MIDDLE 16.3 29 8 47 114.9 MIDDLE 19.1 31 7 48 106.1 MIDDLE 19.3 31 8 49 95.4 TOP 18.9 2 1 50 92.9 TOP 18.4 2 2 51 93.4 TOP 18.7 2 3 52 90.3 TOP 21.4 4 1 53 89.8 TOP 22.2 4 2 54 82.8 TOP 17.2 11 1 55 75.2 TOP 17.0 11 2 56 85.4 TOP 17.7 11 3 57 99.9 TOP 18.8 12 4 58 90.4 TOP 19.3 12 5 59 93.5 TOP 21.4 17 5 60 98.7 TOP 20.9 17 6 61 85.5 TOP 18.7 19 6 62 97.7 TOP 18.7 19 7 63 72.2 TOP 17.3 21 4 64 93.8 TOP 17.5 21 5 65 79.2 TOP 17.4 21 6 66 118.3 TOP 21.3 23 3 67 101.6 TOP 20.5 23 4 68 83.9 TOP 17.3 28 7 69 85.7 TOP 16.7 28 8 70 94.4 TOP 18.3 29 8 71 117.0 TOP 21.2 31 7 72 101.1 TOP 21.2 31 8 ##Model##>modelVCP1<-lme(LMA.vcp~Canopy.position.vcp+Leaf.height.vcp,random=~1 |plots.vcp/tree.vcp,data=data.vcp,na.action=na.omit) ##lme model##>modelVCP1_glht<-glht(modelVCP1,linfct=mcp(Canopy.position.vcp="Tukey")) >summary(modelVCP1_glht) ##summary output of Tukey's test##*Simultaneous Tests for General Linear Hypotheses* * * *Multiple Comparisons of Means: Tukey Contrasts* * * * * *Fit: lme.formula(fixed = LMA.vcp ~ Canopy.position.vcp + Leaf.height.vcp, * * data = data.vcp, random = ~1 | plots.vcp/tree.vcp, na.action = na.omit) * * * *Linear Hypotheses:* * Estimate Std. Error z value Pr(>|z|) * *MIDDLE - BOTTOM == 0 -8.926 2.777 -3.214 0.00383 ** * *TOP - BOTTOM == 0 -19.817 4.080 -4.857 < 0.001 **** *TOP - MIDDLE == 0 -10.891 2.769 -3.934 < 0.001 **** *---* *Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 * *(Adjusted p values reported -- single-step method)* * * * * Based on the output it appears that the mean (LMA) is greatest at the bottom and lowest at the bottom. However, when I receive the code letters (a,b,c, etc) from line and from a graph, the output suggests otherwise.>modelVCP1_glht_cld<-cld(modelVCP1_glht,level=0.05,decreasing=FALSE) >modelVCP1_glht_cld >old.par<-par(mai=c(1,1,1.25,1)) >plot(modelVCP1_glht_cld,col=c("yellow","red","blue")) >par(old.par)First Question: As you can see, the code letters in the boxplot graph seem to contradict the plot itself as well as the Tukey output that I displayed above. Can anyone explain why this may be happening? Second Question: Are the estimate values in the output the difference between the mean of one level and another level? Third Question: If results from the Tukey's test are reported, does one typically report the difference between levels or the means with code letters? If so, how are the means extracted from the model "modelVCP1_glht"? Thanks! -- Adam P. Coble Ph.D. Student Michigan Technological University School of Forest Resources and Environmental Science Houghton, MI 49931 [[alternative HTML version deleted]]
Maybe Matching Threads
- Multiple comparisons with a mixed effects model
- Odd anova(lm()) order phenomenon, looking for an explanation
- Re : using predict() or fitted() from a model with offset; unsolved, included reproducible code
- Grouped data objects within GLS and Variogram
- Please Help me!