R语言基础 | 卡方检验与列联表、独立性检验

---

专注系列化、高质量的R语言教程

推文索引 | 联系小编 | 付费合集

---

上篇介绍了三大检验的t检验和F检验，本篇接着介绍卡方检验。相比于前两者，我们其实更早地接触到卡方检验，它在高中数学教材中就已经出现，但用的却相对较少。

本篇目录如下：

• 1 22列联表下的卡方检验

• 2 一般情况的卡方检验

• 3 chisq.test函数

1 22列联表下的卡方检验

卡方检验适用于计数事件的独立性检验。在高中阶段见到的就是22列联表（contingency table）。

以北师大版高中数学教材为例，为了研究吸烟与患肺癌是否有联系，有如下调查数据：

人数	患肺癌	未换肺癌
吸烟	56	1932
不吸烟	23	4567

如果吸烟者和不吸烟者的患肺癌概率具有显著差异，那么就可以认为吸烟和患肺癌之间存在联系，即二者不独立。

更一般的表示是：

频数	A1	A2
B1	a	b
B2	c	d

我们可以通过卡方统计量来检验变量A和B是否有关联。卡方统计量的计算公式为：

的数值越大，越有把握判定变量A和B存在关联（即不独立）。具体地，

• 当时，没有充分证据表明A和B有关联，即认为A和B相互独立；

• 当时，有90%的把握认为A和B有关联；

• 当时，有95%的把握认为A和B有关联；

• 当时，有99%的把握认为A和B有关联。

以上就是我们在高中阶段学习到的卡方检验的大致内容。

根据公式计算卡方统计量：

a = 56; b = 1932
c = 23; d = 4567 

(a+b+c+d)*(a*d-b*c)^2/((a+b)*(c+d)*(a+c)*(b+d))
## [1] 62.69832

= 62.69832 > 6.635，因此我们可以有99%的把握认为吸烟和患肺癌存在关联。

在stats工具包中，用于卡方检验的函数是chisq.test()：

M <- as.table(rbind(c(56, 1932), c(23, 4567))) 
dimnames(M) <- list(cancer = c("Y", "N"),
                    smoke = c("Y", "N")) 
M 
##       smoke
## cancer    Y    N
##      Y   56 1932
##      N   23 4567

chisq.test(M, correct = F)
##  Pearson's Chi-squared test
## 
## data:  M
## X-squared = 62.698, df = 1, p-value = 2.409e-15

2 一般情况的卡方检验

将22列联表的边际频数补充完整，并记：

频数	A1	A2	总计
B1	a	b	a+b
B2	c	d	c+d
总计	a+c	b+d	n

使用频率代替概率，则有

，
，
，
。

若变量A和B没有关联，即相互独立，则应有

上式可看作事件的期望频率，期望频数为，实际频率为，实际频数为。显然，期望频率和实际频率的差距越小，越可以认为变量A和B相互独立。

对于更一般的列联表（和为行、列数），也可以使用类似的方法计算期望频率和期望频数。使用表示实际频数，表示实际频率，表示期望频数，表示期望频率，表示所有事件的频数之和，则卡方统计量为：

若变量A和B相互独立，则服从卡方分布，自由度为。

3 chisq.test函数

chisq.test()函数的语法结构如下：

chisq.test(x, y = NULL, correct = TRUE,
           p = rep(1/length(x), length(x)), rescale.p = FALSE,
           simulate.p.value = FALSE, B = 2000)

例1.1

参数x为列联表数据，一般应为矩阵，元素为非负整数。

## 同第1节的示例
m <- matrix(c(56, 1932, 23, 4567),
            nrow = 2, byrow = T)
m

chisq.test(m, correct = F)

##  Pearson's Chi-squared test
## 
## data:  m
## X-squared = 62.698, df = 1, p-value = 2.409e-15

例1.2

若参数x为单行或单列矩阵，或原子向量，chisq.test()函数会进行拟合优度检验（goodness-of-fit test），即检验x中的频数分布是否符合给定的概率p，若p参数未给定，则默认为等概率。

chisq.test(x = c(56, 1932, 23, 4567), 
           p = c(0.1, 0.35, 0.05, 0.5))
##  Chi-squared test for given probabilities
## 
## data:  c(56, 1932, 23, 4567)
## X-squared = 1391.2, df = 3, p-value < 2.2e-16

• x = c(56, 1932, 23, 4567)表示某事件被分为四组，各组频数分别为56、1932、23、4567；

• p = c(0.1, 0.35, 0.05, 0.5)表示各组的期望概率；

• p-value < 2.2e-16说明实际频率分布显著异于期望频率分布。

例2.1

参数p的元素之和应为1；如不为1，可设置参数rescale.p = TRUE，否则会报错。

chisq.test(c(56, 1932, 23, 4567), 
           p = c(10, 35, 5, 50),
           rescale.p = T)  
##  Chi-squared test for given probabilities
## 
## data:  c(56, 1932, 23, 4567)
## X-squared = 1391.2, df = 3, p-value < 2.2e-16

例2.2

如果是未进行统计的原始数据，可以使用table()函数进行统计：

set.seed(1111)
x = rpois(50,5)
x
##  [1]  5  4  8  3  6 10  8  2  5  3  0  6  6  6  8  3  8  5  8  7  3  4  7 11  2
## [26]  4  1  4  6  3  7  8  5  4  3  1  5  2  4  1  8  5  6  5  3  9  4  8  9  5

table(x) 
## x
##  0  1  2  3  4  5  6  7  8  9 10 11 
##  1  3  3  7  7  8  6  3  8  2  1  1

chisq.test(table(x))
##  Chi-squared test for given probabilities
## 
## data:  table(x)
## X-squared = 21.04, df = 11, p-value = 0.03296

• 由于未设置参数p，默认概率相同，即检验x是否服从均匀分布；

• p-value = 0.03296 < 0.05，因此不能认为x是等概率分布的。

例3

若simulate.p.value = FALSE（默认情况），则p-value是根据理论上的卡方分布的概率密度函数进行计算的；若simulate.p.value = TRUE，则使用蒙特卡洛模拟计算p-value，参数B为模拟次数。

m <- matrix(c(56, 1932, 23, 4567),
            nrow = 2, byrow = T)
chisq.test(m, correct = F)
##  Pearson's Chi-squared test
## 
## data:  m
## X-squared = 62.698, df = 1, p-value = 2.409e-15

chisq.test(m, simulate.p.value = T,
           B = 999)
##  Pearson's Chi-squared test with simulated p-value (based on 999
##  replicates)
## 
## data:  m
## X-squared = 62.698, df = NA, p-value = 0.001