r - 在 iris 数据集中的每个物种之间执行多个成对比较

Question

我想在 iris 数据集中的每个物种之间进行成对比较(使用 t 检验)，以查看哪些物种在哪些变量方面存在显着差异。也就是说，每个成对比较会将给定变量中一个物种的所有测量值与同一变量中另一个物种的所有测量值进行比较。下面列出了所有可能的与 iris 数据集的成对比较。

data(iris)
setosa.only <- iris[iris$Species == "setosa", ]
versicolor.only <- iris[iris$Species == "versicolor", ]
virginica.only <- iris[iris$Species == "virginica", ]

# setosa vs versicolor
t.test(setosa.only$Sepal.Length, versicolor.only$Sepal.Length)
t.test(setosa.only$Sepal.Width, versicolor.only$Sepal.Width)
t.test(setosa.only$Petal.Length, versicolor.only$Petal.Length)
t.test(setosa.only$Petal.Width, versicolor.only$Petal.Width)

# setosa vs virginica
t.test(setosa.only$Sepal.Length, virginica.only$Sepal.Length)
t.test(setosa.only$Sepal.Width, virginica.only$Sepal.Width)
t.test(setosa.only$Petal.Length, virginica.only$Petal.Length)
t.test(setosa.only$Petal.Width, virginica.only$Petal.Width)

# versicolor vs virginica
t.test(versicolor.only$Sepal.Length, virginica.only$Sepal.Length)
t.test(versicolor.only$Sepal.Width, virginica.only$Sepal.Width)
t.test(versicolor.only$Petal.Length, virginica.only$Petal.Length)
t.test(versicolor.only$Petal.Width, virginica.only$Petal.Width)

这样的成对比较很容易用一个小的数据集来一个一个地执行，比如 iris(它只有 12 种可能的比较)，但我想将它应用到有几十个物种和变量(因此有数百个)的更大的数据集可能的比较)。我如何使用单个或几个命令进行上述比较以将它们应用于更大的数据集？由于对 R 语言的了解有限，我一直无法弄清楚如何做到这一点，如果有人有建议，我将不胜感激。

此外，我希望得到一个总结所有成对比较的输出。它可以是一个带有 TRUE 或 FALSE(或类似 1/0 或 Y/N 的等价物)的矩阵，指示哪些物种在哪些变量方面存在显着差异(即，TRUE 指示满足 t 检验的物种对，考虑到 p = 0.05)。如果这样的矩阵同时包含所有物种和所有变量，则可能难以解释，因此它可能是每个变量一个矩阵。例如，由 Sepal.Length 的比较产生的所需输出矩阵类似于:

            setosa   versicolor   virginica
setosa      NA       YES          YES   
versicolor  YES      NA           YES   
virginica   YES      YES          NA   

或者，输出可以是一个数组，就像在调用下面的代码时返回的那样:

tapply(X = iris$Sepal.Length, INDEX = iris$Species, FUN = summary)

Answer 1

我们可以使用 combn 在 Species 列中创建唯一值的组合。对于每个组合，将 t.test 应用于数据集的每一列。

res <- combn(unique(iris$Species), 2, function(x) {
  data1 <- subset(iris, Species == x[1], -Species)
  data2 <- subset(iris, Species == x[2], -Species)
  out <- data.frame(species1 = x[1], species2 = x[2], column = names(data1))
  out$t.test <- Map(t.test, data1, data2)
  out
}, simplify = FALSE)

res 是一个dataframe的列表，合并成一个dataframe

res <- dplyr::bind_rows(res)

然后您可以从每个 t.test 输出中提取 p.value。

unname(sapply(res$t.test, `[[`, "p.value"))