首先, 让我们获取一些数据。 MASS软件包包含1993年在美国销售的93辆汽车的数据。这些数据存储在Cars93对象中, 并且每辆汽车包含27个特征, 其中一些是分类的。因此, 让我们加载MASS软件包, 并查看cars93中包含的车辆类型：
library(MASS)
Cars93$Type
## [1] Small Midsize Compact Midsize Midsize Midsize Large Large
## [9] Midsize Large Midsize Compact Compact Sporty Midsize Van
## [17] Van Large Sporty Large Compact Large Small Small
## [25] Compact Van Midsize Sporty Small Large Small Small
## [33] Compact Sporty Sporty Van Midsize Large Small Sporty
## [41] Sporty Small Compact Small Small Sporty Midsize Midsize
## [49] Midsize Midsize Midsize Large Small Small Compact Van
## [57] Sporty Compact Midsize Sporty Midsize Small Midsize Small
## [65] Compact Van Midsize Compact Midsize Van Large Sporty
## [73] Small Compact Sporty Midsize Large Compact Small Small
## [81] Small Compact Small Small Sporty Midsize Van Small
## [89] Van Compact Sporty Compact Midsize
## Levels: Compact Large Midsize Small Sporty Van
那里有6种类型的汽车。 table函数告诉我们每种类型有多少种：
table(Cars93$Type)
##
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
prop.table将其转换为分数：
prop.table(table(Cars93$Type))
##
## Compact Large Midsize Small Sporty Van
## 0.17204301 0.11827957 0.23655914 0.22580645 0.15053763 0.09677419
与汽车的起源相同：
table(Cars93$Origin)
##
## USA non-USA
## 48 45
prop.table(table(Cars93$Origin))
##
## USA non-USA
## 0.516129 0.483871
太好了, 我们发现我们的数据集包含相似数量的美国和非美国汽车, 并且最普遍的类型是中型和小型。但是, 也许美国和非美国在类型上有所不同？
让我们来看看汽车的种类和来源。我们可以再次使用table, 但是现在有两个参数。第一个将成为行变量, 第二个将成为列变量：
table(Cars93$Type, Cars93$Origin)
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
现在, 我们看到了每个人都知道的：美国人喜欢大型车辆！上表显示了两个类别变量(类型和来源)的联合分布。这样的表称为列联表。
rowSums和colSums函数确实是不言自明的
(tab1<-table(Cars93$Type, Cars93$Origin))
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
rowSums(tab1)
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
colSums(tab1)
## USA non-USA
## 48 45
与表格嵌套的prop.table给出频率
prop.table(table(Cars93$Type, Cars93$Origin))
##
## USA non-USA
## Compact 0.07526882 0.09677419
## Large 0.11827957 0.00000000
## Midsize 0.10752688 0.12903226
## Small 0.07526882 0.15053763
## Sporty 0.08602151 0.06451613
## Van 0.05376344 0.04301075
转换为百分数只是乘以100：
prop.table(table(Cars93$Type, Cars93$Origin))*100
##
## USA non-USA
## Compact 7.526882 9.677419
## Large 11.827957 0.000000
## Midsize 10.752688 12.903226
## Small 7.526882 15.053763
## Sporty 8.602151 6.451613
## Van 5.376344 4.301075
注意, 这是一个联合概率分布, 从中我们可以看到, 例如, 约有7.5％的汽车是小型汽车, 并且是美国原产的。
通常, 我们对一个变量在另一个变量创建的组中的分布感兴趣。在这里, 美国和(分别)非美国汽车之间的汽车类型分布似乎很有趣。为此, 我们使用prop.table函数的margin参数。它告诉分组变量在行(margin = 1)或列(margin = 2)中的哪个位置：
prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100
##
## USA non-USA
## Compact 14.583333 20.000000
## Large 22.916667 0.000000
## Midsize 20.833333 26.666667
## Small 14.583333 31.111111
## Sporty 16.666667 13.333333
## Van 10.416667 8.888889
现在, 我们可以轻松地看到, 在美国以外的地区, 小型汽车的使用频率是数据集中美国部分地区的两倍。
还要注意, 百分比在列中总计为100, 而在联合分配表(不带边距参数的百分比)中, 100是整个表的总和。
(tab2<-prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100)
##
## USA non-USA
## Compact 14.583333 20.000000
## Large 22.916667 0.000000
## Midsize 20.833333 26.666667
## Small 14.583333 31.111111
## Sporty 16.666667 13.333333
## Van 10.416667 8.888889
colSums(tab2)
## USA non-USA
## 100 100
列联表中最常见的问题是行和列变量是否独立。回答这个问题的最基本方法是运行卡方检验。本教程将对此进行详细介绍。让我们检查Type和Origin是否独立：
chisq.test(Cars93$Type, Cars93$Origin)
## Warning in chisq.test(Cars93$Type, Cars93$Origin): Chi-squared
## approximation may be incorrect
##
## Pearson's Chi-squared test
##
## data: Cars93$Type and Cars93$Origin
## X-squared = 14.08, df = 5, p-value = 0.01511
显然, 它们不是, 但是我们也得到了卡方近似值可能是错误的警告。这是因为卡方统计量仅近似遵循卡方分布。我们拥有的观察越多, 近似值就越好。每当预期计数之一小于5时, chisq.test函数就会引发上述警告(有关”预期计数”的信息, 请参见上面链接的教程)。
Fisher精确检验是卡方检验的替代方法, 主要在卡方逼近不令人满意时使用。让我们运行它：
fisher.test(Cars93$Type, Cars93$Origin)
##
## Fisher's Exact Test for Count Data
##
## data: Cars93$Type and Cars93$Origin
## p-value = 0.007248
## alternative hypothesis: two.sided
结果与卡方的结果非常相似, 但并非总是如此。 Fisher精确测试的一个主要缺点是, 对于大表(或大样本), 它在计算上变得无效。
另一种选择是所谓的G检验。它的统计量也近似为卡方分布, 但对于小样本, 此近似值比卡方检验所使用的近似值更近。对于G测试, 我们可以使用DescTools包中的GTest函数。结果再次与之前的两个测试非常相似：类型和来源不是独立的。
library(DescTools)
GTest(Cars93$Type, Cars93$Origin)
##
## Log likelihood ratio (G-test) test of independence without
## correction
##
## data: Cars93$Type and Cars93$Origin
## G = 18.362, X-squared df = 5, p-value = 0.002526
在2×2列联表中, 可以通过Yates的连续性校正来增强卡方检验。它简单地从每个|减去0.5。观察-预期|卡方统计中的术语。如果你现在迷路了, 请再次参考本教程。而且, R会在必要时自动应用Yates校正。让我们看看美国和非美国汽车的手动变速箱版本的可用性：
(tab3<-table(Cars93$Man.trans.avail, Cars93$Origin))
##
## USA non-USA
## No 26 6
## Yes 22 39
chisq.test(tab3)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tab3
## X-squared = 15.397, df = 1, p-value = 8.712e-05
耶茨的校正是自动应用的, R告诉了我们。
让我们一次查看Man.trans.avail, Origin和Type。表根据提供的第三变量对数据集进行拆分：
table(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
## , , = Compact
##
##
## USA non-USA
## No 2 0
## Yes 5 9
##
## , , = Large
##
##
## USA non-USA
## No 11 0
## Yes 0 0
##
## , , = Midsize
##
##
## USA non-USA
## No 9 4
## Yes 1 8
##
## , , = Small
##
##
## USA non-USA
## No 0 0
## Yes 7 14
##
## , , = Sporty
##
##
## USA non-USA
## No 0 0
## Yes 8 6
##
## , , = Van
##
##
## USA non-USA
## No 4 2
## Yes 1 2
ftable提供了更紧凑的视图：
ftable(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
## Compact Large Midsize Small Sporty Van
##
## No USA 2 11 9 0 0 4
## non-USA 0 0 4 0 0 2
## Yes USA 5 0 1 7 8 1
## non-USA 9 0 8 14 6 2
从上面的表结果中, 应该清楚Origin和Man.trans.avail之间的关系随类型而不同。例如, 对于小型和运动型汽车, 没有关联：每辆美国以及每辆非美国的小型或运动型汽车都带有手动版本。另一方面, 大多数中型美国汽车没有手动版本, 而大多数此类非美国汽车却没有。为了说明层中可能存在的不同关系, 可以使用Cochran-Mantel-Haenszel检验：
mantelhaen.test(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
##
## Mantel-Haenszel chi-squared test with continuity correction
##
## data: Cars93$Man.trans.avail and Cars93$Origin and Cars93$Type
## Mantel-Haenszel X-squared = 8.0153, df = 1, p-value = 0.004638
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 2.226531 76.891307
## sample estimates:
## common odds ratio
## 13.08438
传递给mantelhaen.test的第三个参数标识层。将以上结果与没有分层的结果进行比较(上面的Yates校正示例)。协会仍然存在, 但是证据却很薄弱。
一旦我们发现了变量之间的一些关联, 就该测量它们的强度了。有多种类型的可能的措施。其中有许多描述。现在, 我将重点介绍两种使用最广泛的方法。
V基于卡方统计量：
$$ V = \ sqrt {\ frac {\ chi ^ 2 / N} {\ min(C-1, R-1)}}, $$其中：
- N是列联表(所有单元格的总和)的总计,
- C是列数
- R是行数。
V∈[0; 1]。 V越大, 变量之间的关系越强。 V = 0可以解释为独立性(因为且仅当χ2= 0时, V = 0)。 V的主要缺点是缺乏精确的解释。 V = 0.6是强, 中还是弱关联？
DescTools的CramerV函数可以为我们计算出它：
CramerV(Cars93$Type, Cars93$Origin)
## [1] 0.3890967
再一次：这是强, 中还是弱的联系？
古德曼(Goodman)和克鲁斯卡尔(Kruskal)拉姆达(lambda)是基于变化成比例减少的度量的一个示例。这些措施旨在模仿R2-线性回归的确定系数：
- 他们从[0, 1]取值
- 它们表示由独立变量解释的一部分变化。
让我们使用列变量作为独立变量。然后, 古德曼和克鲁斯卡尔的lambda公式为：$$ \ lambda = \ frac {L-\ sum_j L_j} {L}, $$其中：
- Lj是第j列中的非模态频率之和,
- L是”总计”列中非模态频率的总和
最好用一个例子来解释。让我们看一下Type和Origin, 后者是独立的：
table(Cars93$Type, Cars93$Origin)
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
US列的模式为11, 因此L1 = 7 + 10 + 7 + 8 + 5 = 37；非US列的模式为14, 因此L2 = 9 + 0 + 12 + 6 + 4 = 31。是
rowSums(table(Cars93$Type, Cars93$Origin))
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
模式为22, 因此L = 16 + 11 + 21 + 14 + 9 = 71和$$ \ lambda = \ frac {71-(37 + 31)} {71} = 0.042 $$
我们可以改用DescTools包中的Lambda函数：
Lambda(Cars93$Type, Cars93$Origin, direction='row')
## [1] 0.04225352
direction参数确定因变量在行中或列中的位置。
来源仅解释类型变化的4％。从公式中注意到, lambda将变异定义为属于最大群与不属于最大群之间的二分法。
值得一提的是, 每列的模态类别相同时, lambda为零。看看这张桌子：
(lambdaTab<-cbind(c(0, 100), c(49, 51)))
## [, 1] [, 2]
## [1, ] 0 49
## [2, ] 100 51
chisq.test(lambdaTab)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: lambdaTab
## X-squared = 62.279, df = 1, p-value = 2.981e-15
模态类别对于每一列都是相同的(第二行)。因此, 尽管存在显着可见的关联, 但lambda必须为零：
Lambda(lambdaTab, direction='row')
## [1] 0
如果你想了解有关R中统计信息的更多信息, 请参加srcmini的R中统计模型(第1部分)课程。
