如何在R中将wilcox.test应用于整个数据框?

11

我有一个数据框,其中有一个分组因子(第一列)具有多个级别(超过两个),以及几列数据。我想将 wilcox.test 应用于整个数据框,以比较每个组变量之间的差异。我该怎么办?

更新: 我知道 wilcox.test 只能测试两个组之间的差异,而我的数据框包含三个组。但我更关心的是如何做到这一点,而不是使用哪种测试方法。很可能会删除一个组,但我还没有决定,因此我想测试所有变体。

以下是示例:

structure(list(group = c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), var1 = c(9.3, 
9.05, 7.78, 7.11, 7.14, 8.12, 7.5, 7.84, 7.8, 7.52, 8.84, 6.98, 
6.1, 6.89, 6.5, 7.5, 7.8, 5.5, 6.61, 7.65, 7.68), var2 = c(11L, 
11L, 10L, 1L, 3L, 7L, 11L, 11L, 11L, 11L, 4L, 1L, 1L, 1L, 2L, 
2L, 1L, 4L, 8L, 8L, 1L), var3 = c(7L, 11L, 3L, 7L, 11L, 2L, 11L, 
5L, 11L, 11L, 5L, 11L, 11L, 2L, 9L, 9L, 3L, 8L, 11L, 11L, 2L), 
    var4 = c(11L, 11L, 11L, 11L, 6L, 11L, 11L, 11L, 10L, 7L, 
    11L, 2L, 11L, 3L, 11L, 11L, 6L, 11L, 1L, 11L, 11L), var5 = c(11L, 
    1L, 2L, 2L, 11L, 11L, 1L, 10L, 2L, 11L, 1L, 3L, 11L, 11L, 
    8L, 8L, 11L, 11L, 11L, 2L, 9L)), .Names = c("group", "var1", 
"var2", "var3", "var4", "var5"), class = "data.frame", row.names = c(NA, 
-21L))

更新

感谢大家提供的所有答案!


wilcox.test 只能测试两组之间的差异。而你的数据框包含了三组。你确定这是你想要的测试吗?如果是,你是否需要进行所有可能的成对比较? - RoyalTS
@RoyalTS,我知道这个。但我更感兴趣的是如何做到这一点,而不是使用什么测试。假设有一组将被移除,但我还没有决定,因此我想测试所有变体。 - Yurié
3个回答

14

pairwise.wilcox.test 函数似乎在这里会很有用;也许像这样使用?

out <- lapply(2:6, function(x) pairwise.wilcox.test(d[[x]], d$group))
names(out) <- names(d)[2:6]
out

如果您只想获取p值,您可以逐个提取这些值并制作成矩阵。

sapply(out, function(x) {
    p <- x$p.value
    n <- outer(rownames(p), colnames(p), paste, sep='v')
    p <- as.vector(p)
    names(p) <- n
    p
})
##         var1      var2      var3 var4      var5
## 2v1 0.5414627 0.8205958 0.4851572    1 1.0000000
## 3v1 0.1778222 0.3479835 1.0000000    1 1.0000000
## 2v2        NA        NA        NA   NA        NA
## 3v2 0.5414627 0.3479835 0.3784941    1 0.6919826

请注意,pairwise.wilcox.test使用Holm方法对多重比较进行了调整;如果您想要尝试其他方法,请查看 p.adjust 参数。


哎呀,这里有一个pairwise.wilcox.test!这比单独对每一组进行wilcox.test要好得多。 - RoyalTS
矩阵出了问题 - 请查看第三行和 var4 列。 - Yurié
我认为这是正确的;第三行正在对自身进行测试,这没有意义,因此是NA;对于var4,未校正的p值足够大,使得Holm校正将它们全部等于1。 - Aaron left Stack Overflow
非常感谢,也感谢@royalts的回答,我从你们每个人那里都学到了一些东西! - Yurié

6

更新我的答案,使其跨列工作


test.fun <- function(dat, col) { 

 c1 <- combn(unique(dat$group),2)
 sigs <- list()
 for(i in 1:ncol(c1)) {
    sigs[[i]] <- wilcox.test(
                   dat[dat$group == c1[1,i],col],
                   dat[dat$group == c1[2,i],col]
                 )
    }
    names(sigs) <- paste("Group",c1[1,],"by Group",c1[2,])

 tests <- data.frame(Test=names(sigs),
                    W=unlist(lapply(sigs,function(x) x$statistic)),
                    p=unlist(lapply(sigs,function(x) x$p.value)),row.names=NULL)

 return(tests)
}


tests <- lapply(colnames(dat)[-1],function(x) test.fun(dat,x))
names(tests) <- colnames(dat)[-1]
# tests <- do.call(rbind, tests) reprints as data.frame

# This solution is not "slow" and outperforms the other answers significantly: 
system.time(
  rep(
   tests <- lapply(colnames(dat)[-1],function(x) test.fun(dat,x)),10000
  )
)

#   user  system elapsed 
#  0.056   0.000   0.053 

并且结果是:

tests

$var1
                Test  W          p
1 Group 1 by Group 2 28 0.36596737
2 Group 1 by Group 3 39 0.05927406
3 Group 2 by Group 3 38 0.27073136

$var2
                Test    W         p
1 Group 1 by Group 2 19.0 0.8205958
2 Group 1 by Group 3 36.5 0.1159945
3 Group 2 by Group 3 40.5 0.1522726

$var3
                Test    W         p
1 Group 1 by Group 2 13.0 0.2425786
2 Group 1 by Group 3 23.5 1.0000000
3 Group 2 by Group 3 41.0 0.1261647

$var4
                Test  W         p
1 Group 1 by Group 2 26 0.4323470
2 Group 1 by Group 3 30 0.3729664
3 Group 2 by Group 3 29 0.9479518

$var5
                Test    W         p
1 Group 1 by Group 2 24.0 0.7100968
2 Group 1 by Group 3 19.0 0.5324295
3 Group 2 by Group 3 17.5 0.2306609

谢谢你,Brandon!为了分析所有变量,我在你的代码中添加了第二个循环(请参见我的问题中的最后一次更新),但我不知道如何进一步调整它。现在只打印了最后一个变量的结果,而不是全部。你能帮忙吗? - Yurié
我喜欢你的回答,但我的意思是你只分析了一个变量,而数据集包含5个变量。 - Yurié
当我第一次阅读您的问题时,我认为您只想针对一个列进行测试。我认为其他两个解决方案更完整,您可以以与我上面描述的方式相同探索结果。虽然我不想评论使用替代测试的情况。 - Brandon Bertelsen
现在没问题了,谢谢!但我遇到了另一个问题。分析我的真实数据后得到的输出结果是一个包含116个元素(大表格)的列表。这个列表能否以某种方式进行格式化或导出,以便我可以将其复制/粘贴到Excel电子表格中? - Yurié
@lurie,我认为最好你提出一个新问题来解决你的新问题。 - Brandon Bertelsen
显示剩余4条评论

6
你可以使用 apply 来遍历列,然后使用匿名函数将列传递给你想要使用的任何测试,如下所示(假设数据框命名为df):
apply(df[-1],2,function(x) kruskal.test(x,df$group))

注意:我使用了Kruskal-Wallis检验,因为它适用于多个组。 如果只有两个组,则上述内容同样适用于使用Wilcoxon检验。
如果您确实希望对所有变量进行成对的Wilcoxon检验,则以下是一个两行代码,将循环遍历所有列和所有成对,并将结果作为列表返回。
group.pairs <- combn(unique(df$group),2,simplify=FALSE)
# this loops over the 2nd margin - the columns - of df and makes each column
# available as x
apply(df[-1], 2, function(x)
             # this loops over the list of group pairs and makes each such pair
             # available as an integer vector y
             lapply(group.pairs, function(y)
                    wilcox.test(x[df$group %in% y],df$group[df$group %in% y])))

1
你可以写成 apply(df[-1],2,kruskal.test,df$group) - ziggystar

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接