Untitled
shao wang
2019年9月17日
1.独立性检验
R提供了多种检验类别型变量独立性的方法。本节中描述的三种检验分别为卡方独立性检验、Fisher精确检验和Cochran-Mantel-Haenszel检验。
1.1卡方独立性检验
你可以使用chisq.test()函数对二维表的行变量和列变量进行卡方独立性检验。示例参见代码清单:
library(vcd) ## Warning: package 'vcd' was built under R version 3.5.2## Loading required package: gridmytable <-xtabs(~Treatment+Improved, data=Arthritis)
chisq.test(mytable)##
## Pearson's Chi-squared test
##
## data: mytable
## X-squared = 13.055, df = 2, p-value = 0.001463mytable <- xtabs(~Improved+Sex, data=Arthritis)
chisq.test(mytable)## Warning in chisq.test(mytable): Chi-squared approximation may be incorrect##
## Pearson's Chi-squared test
##
## data: mytable
## X-squared = 4.8407, df = 2, p-value = 0.08889在结果1中，患者接受的治疗和改善的水平看上去存在着某种关系（p<0.01）。而患者性别和改善情况之间却不存在关系（p>0.05）。这里的p值表示从总体中抽取的样本行变量与列变量是相互独立的概率。由于的概率值很小，所以你拒绝了治疗类型和治疗结果相互独立的原假设。由于的概率不够小，故没有足够的理由说明治疗结果和性别之间是不独立的。代码清单中产生警告信息的原因是，表中的6个单元格之一（男性－一定程度上的改善）有一个小于5的值，这可能会使卡方近似无效。
1.2Fisher精确检验
可以使用fisher.test()函数进行Fisher精确检验。Fisher精确检验的原假设是：边界固定的列联表中行和列是相互独立的。其调用格式为fisher.test(mytable)，其中的mytable是一个二维列联表。示例如下：
mytable <- xtabs(~Treatment+Improved, data=Arthritis)
fisher.test(mytable)##
## Fisher's Exact Test for Count Data
##
## data: mytable
## p-value = 0.001393
## alternative hypothesis: two.sided与许多统计软件不同的是，这里的fisher.test()函数可以在任意行列数大于等于2的二维列联表上使用，但不能用于2×2的列联表。
1.3.Cochran-Mantel-Haenszel检验
mantelhaen.test()函数可用来进行Cochran-Mantel-Haenszel卡方检验，其原假设是，两个名义变量在第三个变量的每一层中都是条件独立的。下列代码可以检验治疗情况和改善情况在性别的每一水平下是否独立。此检验假设不存在三阶交互作用（治疗情况×改善情况×性别）。
mytable <- xtabs(~Treatment+Improved+Sex, data=Arthritis)
mantelhaen.test(mytable)##
## Cochran-Mantel-Haenszel test
##
## data: mytable
## Cochran-Mantel-Haenszel M^2 = 14.632, df = 2, p-value = 0.0006647结果表明，患者接受的治疗与得到的改善在性别的每一水平下并不独立（分性别来看，用药治疗的患者较接受安慰剂的患者有了更多的改善）。
2.相关性的度量
上一节中的显著性检验评估了是否存在充分的证据以拒绝变量间相互独立的原假设。如果可以拒绝原假设，那么你的兴趣就会自然而然地转向用以衡量相关性强弱的相关性度量。vcd包中的assocstats()函数可以用来计算二维列联表的phi系数、列联系数和Cramer’s V系数。以下码清单给出了一个示例。 ## 2.1 二维列联表的相关性度量
library(vcd)
mytable <- xtabs(~Treatment+Improved, data=Arthritis)
assocstats(mytable)## X^2 df P(> X^2)
## Likelihood Ratio 13.530 2 0.0011536
## Pearson 13.055 2 0.0014626
##
## Phi-Coefficient : NA
## Contingency Coeff.: 0.367
## Cramer's V : 0.394总体来说，较大的值意味着较强的相关性。vcd包也提供了一个kappa()函数，可以计算混淆矩阵的Cohen’s kappa值以及加权的kappa值。（举例来说，混淆矩阵可以表示两位评判者对于一系列对象进行分类所得结果的一致程度。）
3.相关
相关系数可以用来描述定量变量之间的关系。相关系数的符号（±）表明关系的方向（正相关或负相关），其值的大小表示关系的强弱程度（完全不相关时为0，完全相关时为1）。 >本节中，我们将关注多种相关系数和相关性的显著性检验。我们将使用R基础安装中的state.x77数据集，它提供了美国50个州在1977年的人口、收入、文盲率、预期寿命、谋杀率和高中毕业率数据。数据集中还收录了气温和土地面积数据，但为了节约空间，这里将其丢弃。你可以使用help(state.x77)了解数据集的更多信息。除了基础安装以外，我们还将使用psych和ggm包。
3.1 相关的类型
R可以计算多种相关系数，包括Pearson相关系数、Spearman相关系数、Kendall相关系数、偏相关系数、多分格（polychoric）相关系数和多系列（polyserial）相关系数。下面让我们依次理解这些相关系数。
3.1.1 Pearson、Spearman和Kendall相关
- Pearson积差相关系数衡量了两个定量变量之间的线性相关程度。
- Spearman等级相关系数则衡量分级定序变量之间的相关程度。
- Kendall’s Tau相关系数也是一种非参数的等级相关度量。
cor()函数可以计算这三种相关系数，而cov()函数可用来计算协方差。两个函数的参数有很多，其中与相关系数的计算有关的参数可以简化为： cor(x, use= , method= ) 这些参数详述于表7-2中。
>默认参数为use=“everything”和method=“pearson”。你可以在代码清单7-14中看到一个示例。
#协方差和相关系数
states<- state.x77[,1:6]
cov(states)## Population Income Illiteracy Life Exp Murder
## Population 19931683.7588 571229.7796 292.8679592 -407.8424612 5663.523714
## Income 571229.7796 377573.3061 -163.7020408 280.6631837 -521.894286
## Illiteracy 292.8680 -163.7020 0.3715306 -0.4815122 1.581776
## Life Exp -407.8425 280.6632 -0.4815122 1.8020204 -3.869480
## Murder 5663.5237 -521.8943 1.5817755 -3.8694804 13.627465
## HS Grad -3551.5096 3076.7690 -3.2354694 6.3126849 -14.549616
## HS Grad
## Population -3551.509551
## Income 3076.768980
## Illiteracy -3.235469
## Life Exp 6.312685
## Murder -14.549616
## HS Grad 65.237894cor(states)## Population Income Illiteracy Life Exp Murder
## Population 1.00000000 0.2082276 0.1076224 -0.06805195 0.3436428
## Income 0.20822756 1.0000000 -0.4370752 0.34025534 -0.2300776
## Illiteracy 0.10762237 -0.4370752 1.0000000 -0.58847793 0.7029752
## Life Exp -0.06805195 0.3402553 -0.5884779 1.00000000 -0.7808458
## Murder 0.34364275 -0.2300776 0.7029752 -0.78084575 1.0000000
## HS Grad -0.09848975 0.6199323 -0.6571886 0.58221620 -0.4879710
## HS Grad
## Population -0.09848975
## Income 0.61993232
## Illiteracy -0.65718861
## Life Exp 0.58221620
## Murder -0.48797102
## HS Grad 1.00000000 cor(states, method="spearman")## Population Income Illiteracy Life Exp Murder
## Population 1.0000000 0.1246098 0.3130496 -0.1040171 0.3457401
## Income 0.1246098 1.0000000 -0.3145948 0.3241050 -0.2174623
## Illiteracy 0.3130496 -0.3145948 1.0000000 -0.5553735 0.6723592
## Life Exp -0.1040171 0.3241050 -0.5553735 1.0000000 -0.7802406
## Murder 0.3457401 -0.2174623 0.6723592 -0.7802406 1.0000000
## HS Grad -0.3833649 0.5104809 -0.6545396 0.5239410 -0.4367330
## HS Grad
## Population -0.3833649
## Income 0.5104809
## Illiteracy -0.6545396
## Life Exp 0.5239410
## Murder -0.4367330
## HS Grad 1.0000000首个语句计算了方差和协方差，第二个语句则计算了Pearson积差相关系数，而第三个语句计算了Spearman等级相关系数。举例来说，我们可以看到收入和高中毕业率之间存在很强的正相关，而文盲率和预期寿命之间存在很强的负相关。 请注意，在默认情况下得到的结果是一个方阵（所有变量之间两两计算相关）。你同样可以
- 计算非方形的相关矩阵。观察以下示例：
x <- states[,c("Population", "Income", "Illiteracy", "HS Grad")]
y <- states[,c("Life Exp", "Murder")]
cor(x,y)## Life Exp Murder
## Population -0.06805195 0.3436428
## Income 0.34025534 -0.2300776
## Illiteracy -0.58847793 0.7029752
## HS Grad 0.58221620 -0.4879710当你对某一组变量与另外一组变量之间的关系感兴趣时，cor()函数的这种用法是非常实的。注意，上述结果并未指明相关系数是否显著不为0（即，根据样本数据是否有足够的证据得出总体相关系数不为0的结论）。由于这个原因，你需要对相关系数进行显著性检验。
3.1.2 偏相关
偏相关是指在控制一个或多个定量变量时，另外两个定量变量之间的相互关系。你可以使用ggm包中的pcor()函数计算偏相关系数。ggm包没有被默认安装，在第一次使用之前需要先进行安装。函数调用格式为：pcor(u,S)。其中的u是一个数值向量，前两个数值表示要计算相关系数的变量下标，其余的数值为条件变量（即要排除影响的变量）的下标。S为变量的协方差阵。这个示例有助于阐明用法：
library(ggm) ## Warning: package 'ggm' was built under R version 3.5.3## Loading required package: igraph## Warning: package 'igraph' was built under R version 3.5.3##
## Attaching package: 'igraph'## The following objects are masked from 'package:stats':
##
## decompose, spectrum## The following object is masked from 'package:base':
##
## union##
## Attaching package: 'ggm'## The following object is masked from 'package:igraph':
##
## pacolnames(states)## [1] "Population" "Income" "Illiteracy" "Life Exp" "Murder"
## [6] "HS Grad"pcor(c(1,5,2,3,6), cov(states))## [1] 0.3462724本例中，在控制了收入、文盲率和高中毕业率的影响时，人口和谋杀率之间的相关系数为0.346。偏相关系数常用于社会科学的研究中。
3.1.3 其他类型的相关
polycor包中的hetcor()函数可以计算一种混合的相关矩阵，其中包括数值型变量的Pearson积差相关系数、数值型变量和有序变量之间的多系列相关系数、有序变量之间的多分格相关系数以及二分变量之间的四分相关系数。多系列、多分格和四分相关系数都假设有序变量或二分变量由潜在的正态分布导出。请参考此程序包所附文档以了解更多。
3.2 相关性的显著性检验
3.2.1简单相关系数的显著性检验
在计算好相关系数以后，如何对它们进行统计显著性检验呢？常用的原假设为变量间不相关（即总体的相关系数为0）。你可以使用cor.test()函数对单个的Pearson、Spearman和Kendall相关系数进行检验。简化后的使用格式为： cor.test(x, y, alternative = , method = ) 其中的x和y为要检验相关性的变量，alternative则用来指定进行双侧检验或单侧检验（取值为“two.side”、“less”或“greater”），而method用以指定要计算的相关类型（“pearson”、“kendall” 或 “spearman” ）。 当 研 究 的 假 设 为 总 体 的 相 关 系 数 小 于 0 时，请使用alternative=“less” 。在研究的假设为总体的相关系数大于 0 时，应使用alternative=“greater”。在默认情况下，假设为alternative=“two.side”（总体相关系数不等于0）。参考代码清单:
#检验某种相关系数的显著性
cor.test(states[,3], states[,5])##
## Pearson's product-moment correlation
##
## data: states[, 3] and states[, 5]
## t = 6.8479, df = 48, p-value = 1.258e-08
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.5279280 0.8207295
## sample estimates:
## cor
## 0.7029752这段代码检验了预期寿命和谋杀率的Pearson相关系数为0的原假设。假设总体的相关度为0，则预计在一千万次中只会有少于一次的机会见到0.703这样大的样本相关度（即p=1.258e–08）。由于这种情况几乎不可能发生，所以你可以拒绝原假设，从而支持了要研究的猜想，即预期寿命和谋杀率之间的总体相关度不为0。遗憾的是，cor.test()每次只能检验一种相关关系。但幸运的是，psych包中提供的corr.test()函数可以一次做更多事情。corr.test()函数可以为Pearson、Spearman或Kendall相关计算相关矩阵和显著性水平。代码清单中给出了一个示例。
library(psych)## Warning: package 'psych' was built under R version 3.5.3 corr.test(states, use="complete")## Call:corr.test(x = states, use = "complete")
## Correlation matrix
## Population Income Illiteracy Life Exp Murder HS Grad
## Population 1.00 0.21 0.11 -0.07 0.34 -0.10
## Income 0.21 1.00 -0.44 0.34 -0.23 0.62
## Illiteracy 0.11 -0.44 1.00 -0.59 0.70 -0.66
## Life Exp -0.07 0.34 -0.59 1.00 -0.78 0.58
## Murder 0.34 -0.23 0.70 -0.78 1.00 -0.49
## HS Grad -0.10 0.62 -0.66 0.58 -0.49 1.00
## Sample Size
## [1] 50
## Probability values (Entries above the diagonal are adjusted for multiple tests.)
## Population Income Illiteracy Life Exp Murder HS Grad
## Population 0.00 0.59 1.00 1.0 0.10 1
## Income 0.15 0.00 0.01 0.1 0.54 0
## Illiteracy 0.46 0.00 0.00 0.0 0.00 0
## Life Exp 0.64 0.02 0.00 0.0 0.00 0
## Murder 0.01 0.11 0.00 0.0 0.00 0
## HS Grad 0.50 0.00 0.00 0.0 0.00 0
##
## To see confidence intervals of the correlations, print with the short=FALSE option参数use=的取值可为“pairwise”或“complete”（分别表示对缺失值执行成对删除或行删除）。参数method=的取值可为“pearson”（默认值）、“spearman”或“kendall”。这里可以看到，人口数量和高中毕业率的相关系数（–0.10）并不显著地不为0（p=0.5）。
3.2.2 其他显著性检验
在多元正态性的假设下，psych包中的pcor.test()函数可以用来检验在控制一个或多个额外变量时两个变量之间的条件独立性。使用格式为： pcor.test(r, q, n) 其中的r是由pcor()函数计算得到的偏相关系数，q为要控制的变量数（以数值表示位置），n为样本大小。在结束这个话题之前应当指出的是，psych包中的r.test()函数提供了多种实用的显著性 检验方法。此函数可用来检验： ???-某种相关系数的显著性； - 两个独立相关系数的差异是否显著； - 两个基于一个共享变量得到的非独立相关系数的差异是否显著； - 两个基于完全不同的变量得到的非独立相关系数的差异是否显著。 参阅help(r.test)以了解详情。
4.t 检验
在研究中最常见的行为就是对两个组进行比较。接受某种新药治疗的患者是否较使用某种现有药物的患者表现出了更大程度的改善？某种制造工艺是否较另外一种工艺制造出的不合格品更少？两种教学方法中哪一种更有效？如果你的结果变量是类别型的，那么可以直接使用7.3节中阐述的方法。这里我们将关注结果变量为连续型的组间比较，并假设其呈正态分布。 >为了阐明方法，我们将使用MASS包中的UScrime数据集。它包含了1960年美国47个州的刑罚制度对犯罪率影响的信息。我们感兴趣的结果变量为Prob（监禁的概率）、U1（1424岁年龄段城市男性失业率）和U2（3539岁年龄段城市男性失业率）。类别型变量So（指示该州是否位于南方的指示变量）将作为分组变量使用。数据的尺度已被原始作者缩放过。（注意，我原本打算将本节命名为“旧南方的罪与罚”，但是最后理智还是战胜了情感。）
4.1 独立样本的 t 检验
如果你在美国的南方犯罪，是否更有可能被判监禁？我们比较的对象是南方和非南方各州，因变量为监禁的概率。一个针对两组的独立样本t检验可以用于检验两个总体的均值相等的假设。这里假设两组数据是独立的，并且是从正态总体中抽得。检验的调用格式为： t.test(y ~ x, data) 其中的y是一个数值型变量，x是一个二分变量。调用格式或为： t.test(y1, y2) 其中的y1和y2为数值型向量（即各组的结果变量）。可选参数data的取值为一个包含了这些变量的矩阵或数据框。与其他多数统计软件不同的是，这里的t检验默认假定方差不相等，并使用Welsh的修正自由度。你可以添加一个参数var.equal=TRUE以假定方差相等，并使用合并方差估计。默认的备择假设是双侧的（即均值不相等，但大小的方向不确定）。你可以添加一个参数alternative=“less”或alternative=“greater”来进行有方向的检验。 在下列代码中，我们使用了一个假设方差不等的双侧检验，比较了南方（group 1）和非南方（group 0）各州的监禁概率：
library(MASS)
t.test(Prob ~ So, data=UScrime)##
## Welch Two Sample t-test
##
## data: Prob by So
## t = -3.8954, df = 24.925, p-value = 0.0006506
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.03852569 -0.01187439
## sample estimates:
## mean in group 0 mean in group 1
## 0.03851265 0.06371269你可以拒绝南方各州和非南方各州拥有相同监禁概率的假设（p<0.001）。
注意 由于结果变量是一个比例值，你可以在执行t检验之前尝试对其进行正态化变换。在本例中，所有对结果变量合适的变换（Y/1-Y、log(Y/1-Y)、arcsin(Y)、arcsin(sqrt(Y)）都会将检验引向相同的结论。数据变换详述于第8章。
4.2非独立样本的 t 检验
再举个例子，你可能会问：较年轻（1424岁）男性的失业率是否比年长（3539岁）男性的失业率更高？在这种情况下，这两组数据并不独立。你不能说亚拉巴马州的年轻男性和年长男性的失业率之间没有关系。在两组的观测之间相关时，你获得的是一个非独立组设计（dependent groups design）。前???后测设计（pre-postdesign）或重复测量设计（repeated measures design）同样也会产生非独立的组。非独立样本的t检验假定组间的差异呈正态分布。对于本例，检验的调用格式为：t.test(y1, y2, paired=TRUE) 其中的y1和y2为两个非独立组的数值向量。结果如下：
library(MASS)
sapply(UScrime[c("U1","U2")], function(x)(c(mean=mean(x),sd=sd(x))))## U1 U2
## mean 95.46809 33.97872
## sd 18.02878 8.44545with(UScrime, t.test(U1, U2, paired=TRUE))##
## Paired t-test
##
## data: U1 and U2
## t = 32.407, df = 46, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 57.67003 65.30870
## sample estimates:
## mean of the differences
## 61.48936差异的均值（61.5）足够大，可以保证拒绝年长和年轻男性的平均失业率相同的假设。年轻男性的失业率更高。事实上，若总体均值相等，获取一个差异如此大的样本的概率小于0.000 000 000 000 000 22（即2.2e–16）。
4.3多于两组的情况
如果想在多于两个的组之间进行比较，应该怎么做？如果能够假设数据是从正态总体中独立抽样而得的，那么你可以使用方差分析（ANOVA）。ANOVA是一套覆盖了许多实验设计和准实验设计的综合方法。就这一点来说，它的内容值得单列一章。
5.组间差异的非参数检验
如果数据无法满足t检验或ANOVA的参数假设，可以转而使用非参数方法。举例来说，若结果变量在本质上就严重偏倚或呈现有序关系，那么你可能会希望使用本节中的方法。 ## 5.1 两组的比较 若两组数据独立，可以使用Wilcoxon秩和检验（更广为人知的名字是Mann-Whitney U检验）来评估观测是否是从相同的概率分布中抽得的（即，在一个总体中获得更高得分的概率是否比另一个总体要大）。调用格式为：wilcox.test(y ~ x, data) 其中的y是数值型变量，而x是一个二分变量。调用格式或为：wilcox.test(y1, y2),其中的y1和y2为各组的结果变量。可选参数data的取值为一个包含了这些变量的矩阵或数据框。默认进行一个双侧检验。你可以添加参数exact来进行精确检验，指定alternative=“less”或alternative=“greater”进行有方向的检验。 如果你使用Mann-Whitney U检验回答上一节中关于监禁率的问题，将得到这些结果：
with(UScrime, by(Prob, So, median))## So: 0
## [1] 0.038201
## --------------------------------------------------------
## So: 1
## [1] 0.055552 wilcox.test(Prob ~ So, data=UScrime)##
## Wilcoxon rank sum test
##
## data: Prob by So
## W = 81, p-value = 8.488e-05
## alternative hypothesis: true location shift is not equal to 0你可以再次拒绝南方各州和非南方各州监禁率相同的假设（p<0.001）。Wilcoxon符号秩检验是非独立样本t检验的一种非参数替代方法。它适用于两组成对数据和 无法保证正态性假设的情境。调用格式与Mann-Whitney U检验完全相同，不过还可以添加参数paired=TRUE。让我们用它解答上一节中的失业率问题：
sapply(UScrime[c("U1","U2")], median)## U1 U2
## 92 34 with(UScrime, wilcox.test(U1, U2, paired=TRUE))## Warning in wilcox.test.default(U1, U2, paired = TRUE): cannot compute exact
## p-value with ties##
## Wilcoxon signed rank test with continuity correction
##
## data: U1 and U2
## V = 1128, p-value = 2.464e-09
## alternative hypothesis: true location shift is not equal to 0你再次得到了与配对t检验相同的结论。在本例中，含参的t检验和与其作用相同的非参数检验得到了相同的结论。当t检验的假设合理时，参数检验的功效更强（更容易发现存在的差异）。而非参数检验在假设非常不合理时（如对于等级有序数据）更适用。 ## 5.2多于两组的比较 在要比较的组数多于两个时，必须转而寻求其他方法。考虑state.x77数据集。它包含了美国各州的人口、收入、文盲率、预期寿命、谋杀率和高中毕业率数据。如果你想比较美国四个地区（东北部、南部、中北部和西部）的文盲率，应该怎么做呢？这称为单向设计（one-way design），我们可以使用参数或非参数的方法来解决这个问题。如果无法满足ANOVA设计的假设，那么可以使用非参数方法来评估组间的差异。如果各组 独立，则Kruskal-Wallis检验将是一种实用的方法。如果各组不独立（如重复测量设计或随机区组设计），那么Friedman检验会更合适。Kruskal-Wallis检验的调用格式为：kruskal.test(y ~ A, data) ,其中的y是一个数值型结果变量，A是一个拥有两个或更多水平的分组变量（grouping variable）。（若有两个水平，则它与Mann-Whitney U检验等价。）而Friedman检验的调用格式为：friedman.test(y ~ A | B, data),其中的y是数值型结果变量，A是一个分组变量，而B是一个用以认定匹配观测的区组变量（blocking variable）。在以上两例中，data皆为可选参数，它指定了包含这些变量的矩阵或数据框。 让我们利用Kruskal-Wallis检验回答文盲率的问题。首先，你必须将地区的名称添加到数据集中。这些信息包含在随R基础安装分发的state.region数据集中。
states <- data.frame(state.region, state.x77)
kruskal.test(Illiteracy ~ state.region, data=states)##
## Kruskal-Wallis rank sum test
##
## data: Illiteracy by state.region
## Kruskal-Wallis chi-squared = 22.672, df = 3, p-value = 4.726e-05显著性检验的结果意味着美国四个地区的文盲率各不相同（p<0.001）。 虽然你可以拒绝不存在差异的原假设，但这个检验并没有告诉你哪些地区显著地与其他地区不同。要回答这个问题，你可以使用Wilcoxon检验每次比较两组数据。一种更为优雅的方法是在控制犯第一类错误的概率（发现一个事实上并不存在的差异的概率）的前提下，执行可以同步进行的多组比较，这样可以直接完成所有组之间的成对比较。我写的函数wmc()可以实现这一目的，它每次用Wilcoxon检验比较两组，并通过p.adj()函数调整概率值。 说实话，我将本章标题中基本的定义拓展了不止一点点，但由于在这里讲非常合适，所以希 望你能够容忍我的做法。你可以从上下载到一个包含wmc()111函数的文本文件。代码清单7-17通过这个函数比较了美国四个区域的文盲率。
source("http://www.statmethods.net/RiA/wmc.txt")
states <- data.frame(state.region, state.x77)
wmc(Illiteracy ~ state.region, data=states, method="holm")## Descriptive Statistics
##
## West North Central Northeast South
## n 13.00000 12.00000 9.00000 16.00000
## median 0.60000 0.70000 1.10000 1.75000
## mad 0.14826 0.14826 0.29652 0.59304
##
## Multiple Comparisons (Wilcoxon Rank Sum Tests)
## Probability Adjustment = holm
##
## Group.1 Group.2 W p
## 1 West North Central 88.0 8.665618e-01
## 2 West Northeast 46.5 8.665618e-01
## 3 West South 39.0 1.788186e-02 *
## 4 North Central Northeast 20.5 5.359707e-02 .
## 5 North Central South 2.0 8.051509e-05 ***
## 6 Northeast South 18.0 1.187644e-02 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1source()函数下载并执行了定义wmc()函数的R脚本。函数的形式是wmc(y ~ A, data,method)，其中y是数值输出变量，A是分组变量，data是包含这些变量的数据框，method指定限制I类误差的方法。代码清单7.17使用的是基于Holm（1979）提出的调整方法，它可以很大程度地控制总体I类误差率（在一组成对比较中犯一次或多次I类错误的概率）。参阅help(p.adjust) 以查看其他可供选择的方法描述。wmc()函数首先给出了样本量、样本中位数、每组的绝对中位差。其中，西部地区（West）的文盲率最低，南部地区（South）文盲率最高。然后，函数生成了六组统计比较（南部与中北部（North Central）、西部与东北部（Northeast）、西部与南部、中北部与东北部、中北部与南部、东北部与南部）。可以从双侧p值（p）看到，南部与其他三个区域有明显差别，但当显著性水平p＜0.05时，其他三个区域间并没有统计显著的差别。
非参数多组比较很有用，但在R中的实现并不轻松。第21章中，你将有机会将wmc()函数拓展为一个完整的可做误差检验、信息图表的R包。