4
votes

Is there a way to get significance codes after a pairwise comparisons to a Kruskall wallis test? With significance codes I mean letter codes that are assigned to populations to indicate where differences are significant.

With a traditional anova such a test can be performed using HSD.test from the agricolae library but for non parametric counterparts of anova I have not been able to find anything.

A small toy example:

dv  <-  c(runif(100, 5.0, 10))
iv  <-  as.factor( c(rep("I", 10),  rep("II", 10),  rep("III", 10),  rep("IV", 10), rep("V", 10),
                    rep("VI", 10), rep("VII", 10), rep("VIII", 10), rep("IX", 10), rep("X", 10)))

df  <-  data.frame(dv, iv)

# with anova
library(agricolae)
aov.000  <-  aov(dv ~ iv,  data=df)
HSD.test(aov.000, "iv")

# after KW test: 
(kt  <-  kruskal.test(dv ~ iv,  data=df))

library(coin)
library(multcomp)
NDWD <- oneway_test(dv ~ iv, data = df,
        ytrafo = function(data) trafo(data, numeric_trafo = rank),
        xtrafo = function(data) trafo(data, factor_trafo = function(x)
            model.matrix(~x - 1) %*% t(contrMat(table(x), "Tukey"))),
        teststat = "max", distribution = approximate(B=1000))

### global p-value
print(pvalue(NDWD))

### sites (I = II) != (III = IV) at alpha = 0.01 (page 244)
print(pvalue(NDWD, method = "single-step"))
4

4 Answers

5
votes

because it might be of use to others, the following seems to work (using the multcompView library):

library(multcompView)
mat  <-  data.frame(print(pvalue(NDWD, method = "single-step")))
(a   <-  c(mat[, 1]));  names(a)  <-  rownames(mat)
multcompLetters(a)

Alternatively the following will work:

test  <-  pairwise.wilcox.test(dv, iv, p.adj="bonferroni", exact=FALSE)
# test  <-  pairwise.wilcox.test(et.ef, s.t, p.adj="holm", exact=FALSE)

library(multcompView)
test$p.value
library(reshape)
(a <- melt(test$p.value))
a.cc  <-  na.omit(a)
a.pvals  <-  a.cc[, 3]
names(a.pvals)  <-  paste(a.cc[, 1], a.cc[, 2], sep="-")
a.pvals
multcompLetters(a.pvals)
3
votes

You can also use the cldList function from the rcompanion package (see https://rcompanion.org/rcompanion/d_06.html). Example:

k_test <- k_test$res

library(rcompanion)

cldList(comparison = k_test$Comparison,
    p.value    = PT$P.adj,
    threshold  = 0.05)


Error: No significant differences.

I used it in combination with the Dunn post-hoc and it worked perfectly.

2
votes

You can at least do it graphically using the multicomp package:

dv  <-  c(runif(100, 5.0, 10))
iv  <-  as.factor( c(rep("I", 10),  rep("II", 10),  rep("III", 10),  rep("IV", 10), rep("V", 10),
                rep("VI", 10), rep("VII", 10), rep("VIII", 10), rep("IX", 10), rep("X", 10)))
df  <-  data.frame(dv, iv)
anova_results  <-  aov(dv ~ iv,  data=df)
library(multcomp)
tuk <- glht(anova_results, linfct = mcp(iv = "Tukey"))
summary(tuk)          # standard display
tuk.cld <- cld(tuk)   # letter-based display
opar <- par(mai=c(1,1,1.5,1))
plot(tuk.cld)
par(opar)

Of course given your randomly generated data, the resulting plot is not very interesting, but will give you the groupings-

enter image description here

This is one of my plots, using the same approach:

enter image description here Lastly, if you do not want the graphics, you can dig into the package and easily find the string that stores the grouping information to be used elsewhere.

0
votes

If you want the compact letter display for the Kruskal test, the same library agricolae seems to allow it with the function kruskal. Using your own data:

library(agricolae)
kruskal(df$dv, df$iv, group=TRUE, p.adj="bonferroni")$groups
####    trt  means M
#### 1  VI    59.2 a
#### 2  VII   57.0 a
#### 3  IX    56.4 a
#### 4  II    55.0 a
#### ...

(well, in this example the groups are not considered different...)