-1

我习惯在 minitab 中进行 Tukey 事后测试。当我这样做时，我通常会得到因变量/预测变量的家庭分组。

Minitab 输出 - Tukey 成对

在 R 中，TukeyHSD()不显示（或计算？）使用族分组。它仅显示每个因/预测变量之间的关系。是否可以像在 minitab 中一样显示家庭分组？

使用diamonds数据集：

av <- aov(price ~ cut, data = diamonds)
tk <- TukeyHSD(av, ordered = T, which = "cut")
plot(tk)

输出：

Fit: aov(formula = price ~ cut, data = diamonds)

$cut
                    diff        lwr       upr     p adj
Good-Ideal         471.32248  300.28228  642.3627 0.0000000
Very Good-Ideal    524.21792  401.33117  647.1047 0.0000000
Fair-Ideal         901.21579  621.86019 1180.5714 0.0000000
Premium-Ideal     1126.71573 1008.80880 1244.6227 0.0000000
Very Good-Good      52.89544 -130.15186  235.9427 0.9341158
Fair-Good          429.89331  119.33783  740.4488 0.0014980
Premium-Good       655.39325  475.65120  835.1353 0.0000000
Fair-Very Good     376.99787   90.13360  663.8622 0.0031094
Premium-Very Good  602.49781  467.76249  737.2331 0.0000000
Premium-Fair       225.49994  -59.26664  510.2665 0.1950425

添加图片以帮助澄清我对 Maruits 评论的回应： 在此处输入图像描述

4

1 回答 1

0

这是有关如何为数据集重现 minitab 表的分步示例ggplot2::diamonds。我已经尽可能地包含了细节/解释。

请注意，据我所知，minitab 表中显示的结果与 Tukey 事后检验的结果无关。它们基于方差分析的结果。Tukey 的诚实显着性差异 (HSD) 检验是一种事后检验，根据 ANOVA 结果确定哪些比较（在所有可能的组均值成对比较中）是（诚实地）具有统计显着性的。

为了重现 minitabs“均值分组”汇总表（参见 minitab Express Support 的“解释结果：步骤 3”的第一个表），我建议（重新）运行线性模型以提取均值和置信区间。请注意，这正是aov适合每个组的方差分析模型的方式。

拟合线性模型

我们指定一个0偏移量来获得每个组的绝对估计值（而不是相对于偏移量的变化的估计值）。

fit <- lm(price ~ 0 + cut, data = diamonds)
coef <- summary(fit)$coef;
coef;
#    Estimate Std. Error   t value Pr(>|t|)
#cutFair      4358.758   98.78795  44.12236        0
#cutGood      3928.864   56.59175  69.42468        0
#cutVery Good 3981.760   36.06181 110.41487        0
#cutPremium   4584.258   33.75352 135.81570        0
#cutIdeal     3457.542   27.00121 128.05137        0

确定家庭分组

为了获得类似于 minitab 的“家庭分组”的东西，我们采用以下方法：

  1. 计算所有参数的置信区间
  2. 对所有参数的置信区间数据执行层次聚类分析
  3. 在与 CI 的标准偏差相对应的高度切割生成的树。这将为我们提供一组基于置信区间的参数估计值。这是一种有点经验的方法，但也是合理的，因为树测量了置信区间之间的成对距离，标准差可以解释为欧几里得距离。

我们首先计算置信区间，然后使用完全链接的层次聚类对得到的距离矩阵进行聚类。 

CI <- confint(fit);
hc <- hclust(dist(CI));

我们检查集群树状图

plot(hc);

在此处输入图像描述

我们现在在与所有参数估计中所有 CI 的标准差相对应的高度切割树，以获得“家庭分组”

grps <- cutree(hc, h = sd(CI))

总结结果

最后，我们整理所有数量并将结果存储在一个类似于 minitab 的“均值分组”表的表中。

library(tidyverse)
bind_cols(
    cut = rownames(coef),
    N = as.integer(table(fit$model$cut)),
    Mean = coef[, 1],
    Groupings = grps) %>%
    as.data.frame()
#           cut     N     Mean Groupings
#1      cutFair  1610 4358.758         1
#2      cutGood  4906 3928.864         2
#3 cutVery Good 12082 3981.760         2
#4   cutPremium 13791 4584.258         1
#5     cutIdeal 21551 3457.542         3

请注意，我们的结果与 minitab“平均分组”表中的结果几乎完全一致：单独cut = Ideal属于组3（Cminitab 表中的组），而Fair+Premium共享组1(minitab: group A) 和Good+Very Good共享组2(minitab:组B）。

于 2018-06-27T13:43:00.670 回答