I'm trying to obtain the compact letters associated with the p values I generated with the function pairwise.wilcox.test(). I used the function multcompLetters() which works well, but I want the letters to be ordered so that the highest mean gets the letter "a". In the following example, while using the function multcompLetters() to get the letters ordered, I receive the warning : "Error in tapply(data[, fm[[1]]], data[, fm[[2]]], function(x) do.call(mean, : arguments must have same length".
Is there something I don't understand with the usage of the function ?
Thanks a lot.
Adding data frame
As15 <- structure(list(Community = c("Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat", "Arviat",
"Arviat", "Arviat", "Arviat", "Nain", "Nain", "Nain", "Nain",
"Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain",
"Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain",
"Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain",
"Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain", "Nain",
"Nain", "Nain", "Nain", "Nain", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "Resolute", "Resolute", "Resolute",
"Resolute", "Resolute", "Resolute", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour", "SachsHarbour",
"SachsHarbour", "SachsHarbour", "SachsHarbour"), ArsenicD = c(0.0254,
0.121, 0.0726, 0.021, 0.0782, 0.0538, 0.1078, 0.0454, 0.0368,
0.0618, 0.037, 0.0754, 0.1394, 0.0784, 0.218, 0.1482, 0.0778,
0.0592, 0.0314, 0.0232, 0.0548, 0.0662, 0.0604, 0.0252, 0.0502,
0.0768, 0.036, 0.054, 0.0404, 0.0642, 0.0384, 0.0504, 0.0616,
0.068, 0.0678, 0.06, 0.0526, 0.0454, 0.0574, 0.0462, 0.0558,
0.0506, 0.0764, 0.0466, 0.0518, 0.0544, 0.0362, 0.0472, 0.0374,
0.0564, 0.0512, 0.0442, 0.0526, 0.05, 0.047, 0.0456, 0.073, 0.0798,
0.0544, 0.072, 0.0392, 0.0658, 0.1662, 0.036, 0.0562, 0.0584,
0.0416, 0.037, 0.046, 0.0334, 0.0384, 0.0654, 0.0528, 0.0274,
0.0368, 0.0634, 0.0308, 0.0702, 0.0502, 0.058, 0.037, 0.0456,
0.0404, 0.0516, 0.051, 0.0484, 0.054, 0.0344, 0.064, 0.0548,
0.032, 0.0532, 0.0562, 0.0464, 0.0334, 0.068, 0.0422, 0.298,
0.0344, 0.0338, 0.0508, 0.0356, 0.0446, 0.0484, 0.0408, 0.0148,
0.0374, 0.0244, 0.0644, 0.0574, 0.028, 0.0462, 0.067, 0.0472,
0.053, 0.0418, 0.0324, 0.054, 0.04, 0.0506, 0.0592, 0.0356, 0.049,
0.054, 0.296, 0.276, 0.226, 0.0834, 0.452, 0.306, 0.218, 0.33,
0.208, 0.1628, 0.426, 0.376, 0.0894, 0.438, 0.334, 0.212, 0.1606,
0.082, 0.1178, 0.1128, 0.142, 0.103, 0.0862, 0.1104, 0.0746,
0.0954, 0.202, 0.1362, 0.24, 0.21, 0.172, 0.278, 0.1354, 0.274,
0.228, 0.0854, 0.0924, 0.0992, 0.0648, 0.0548, 0.0768, 0.0996,
0.1008, 0.063, 0.0372, 0.1582, 0.228, 0.1514, 0.218, 0.1612,
0.1608, 0.1652, 0.1474, 0.1904, 0.0396, 0.1396, 0.0816, 0.1132,
0.0968, 1.258, 0.326, 0.236, 0.0854, 0.256, 0.1258, 0.1324, 0.1716,
0.1642, 0.22, 0.0836, 0.18, 0.274, 0.1918, 0.171, 0.356, 0.298,
0.1084, 0.1436, 0.238, 0.348, 0.1496, 0.342, 0.434, 0.146, 0.174,
0.288, 0.396, 0.1834, 0.64, 0.342, 0.686, 0.43, 0.1808, 0.454,
0.482, 0.352, 0.476, 0.422, 0.286, 0.81, 0.656, 0.564, 0.119,
0.296, 0.1778, 0.512, 0.278, 0.22, 0.208, 0.1474, 0.55, 0.1526,
0.1858, 0.128, 0.224, 0.1752, 0.0724, 0.1302, 0.1724, 0.1354,
0.104, 0.1132, 0.292, 0.238, 0.348, 0.1626, 0.312, 0.238, 0.236,
0.1276, 0.1228, 0.0978, 0.376, 0.1968, 0.1164, 0.448, 0.2, 0.1322,
0.117)), row.names = c(NA, -263L), class = c("tbl_df", "tbl",
"data.frame"))
Wilcoxon Test
As15Wilcox <- suppressWarnings(pairwise.wilcox.test(As15$ArsenicD, As15$Community, p.adjust.method = "BH"))
As15pvals <- c(na.omit(setNames(c(As15Wilcox$p.value),
do.call("paste", c(as.list(expand.grid(rownames(As15Wilcox$p.value),
colnames(As15Wilcox$p.value))), sep = "-")))))
Obtaining Letters
multcompLetters2(ArsenicD ~ Community, As15pvals, As15)
I verified with length(), and both ArsenicD and Community possess the same length.
It's not clear why you think you can pass the output of pairwise.wilcox.test to the second argument of multcompLetters2 (the parameter x). According to the docs:
x
- One of the following: (1) A square, symmetric matrix with row names. (2) A vector with hyphenated names, which identify individual items or factor levels after "strsplit". (3) An object of class "dist". If x (or x[1]) is not already of class "logical", it is replaced with do.call(compare, list(x, threshold)), which by default converts numbers (typically p-values) less than 0.05 to TRUE and everything else to FALSE. If x is a matrix, its diagonal must be or must convert to FALSE.
...whereas wilcox1 is an object of class "pairwise.htest".
However, it does contain a matrix of p values of the comparisons, which we can convert into one of the above formats. Let's demonstrate with a dataset similar to yours.
library(multcompView)
wilcox1 <- suppressWarnings(
pairwise.wilcox.test(site21$ArsenicD, site21$Community,
p.adjust.method = "BH"))
Using your code, we get an error:
multcompLetters2(ArsenicD ~ Community, wilcox1, site21)
#> Error in vec2mat2(namx, sep): Names must contain exactly one '-' each;
#> instead got method, data.name, p.value, p.adjust.method
Let us take the p values from our wilcox1 object and convert to a named vector with the correct pair names:
pvals <- c(na.omit(setNames(c(wilcox1$p.value),
do.call("paste", c(as.list(expand.grid(rownames(wilcox1$p.value),
colnames(wilcox1$p.value))), sep = "-")))))
pvals
#> North-East South-East West-East South-North West-North West-South
#> 4.746437e-14 1.077333e-04 3.164291e-14 3.164291e-14 3.164291e-14 7.175058e-08
And now we can do:
multcompLetters2(ArsenicD ~ Community, pvals, site21)
#> West South East North
#> "a" "b" "c" "d"
Data used
Obviously, we don't have your data, so here is a reproducible example with the same names and basic structure that I used for this example:
set.seed(1)
site21 <- data.frame(Community = rep(c("North", "East",
"South", "West"), each = 25),
ArsenicD = c(rnorm(25, 5), rnorm(25, 8),
rnorm(25, 9), rnorm(25, 11)))
Created on 2023-07-26 with reprex v2.0.2