# Why do I get different results for Tukey HSD compact letter comparisons with different R packages?

I have a dataset that I have been trying to run Tukey HSD tests on and obtain a compact letter display. I originally used the agricolae package with function HSD.Test(). After reviewing my results, I was curious about comparisons between groups that seemed different from what I expected- not that I refuse to accept these results but because it honestly surprised me after I had reviewed raw data. I went and found the multcomp package to check if it would provide the same results from the first package with the functions glht() and cld(). Below is the code that shows the difference in compact letter displays, which I am ultimately wanting to use in a table. My question is: should I trust one over the other? Or is there a fundamental difference in how these are calculated between packages that would aid me in how to reconcile the different results I am getting.

Here is my setup and run of the agricolae.

#import libraries
library (agricolae)
library(multcomp)

#fit analysis of variance model
SPLW<-aov(sepalLengthToWidth~TAXON, data=garden)

#get compact letter comparisons from agricolae
hsd<-HSD.test(SPLW, "TAXON")
hsd

$statistics MSerror Df Mean CV 0.219363 63 3.019053 15.51354$parameters
test name.t ntr StudentizedRange alpha
Tukey  TAXON   6         4.156707  0.05

$means sepalLengthToWidth std r Min Max Q25 Q50 Q75 AGCP 3.003620 0.4542374 24 2.185729 4.048658 2.781238 3.090727 3.254949 CPTN 2.516216 0.1663721 3 2.369989 2.697232 2.425709 2.481429 2.589330 NOAL 3.517557 0.5924274 13 2.439724 4.167752 3.237401 3.745421 3.983766 OUMT 2.633716 0.3395092 6 2.128205 3.056782 2.449174 2.646832 2.866976 TMAL 3.837945 0.7807588 3 3.290098 4.731935 3.390950 3.491803 4.111869 WGCP 2.781738 0.3981362 20 1.854508 3.302083 2.608329 2.882435 3.011886$comparison
NULL

\$groups
sepalLengthToWidth groups
TMAL           3.837945      a
NOAL           3.517557      a
AGCP           3.003620     ab
WGCP           2.781738      b
OUMT           2.633716      b
CPTN           2.516216      b

attr(,"class")
[1] "group"



Then I use the multcomp method on the same aov model.

#get compact letter comparisons from multcomp
tuk <- glht(SPLW, linfct = mcp(TAXON = "Tukey"))
tuk.cld <- cld(tuk)
tuk.cld

AGCP CPTN NOAL OUMT TMAL WGCP
"a"  "a"  "b"  "a"  "b"  "a"