卡方檢驗(yàn)詳解

“卡方”到底是什么

為什么要叫“卡方”并齐？因?yàn)樵恰癱hi-squared”，一半是音譯，一半是意譯况褪。其中撕贞，chi 是希臘字母 $\chi$ 的讀音，其實(shí)讀音更像是“開”窝剖，而不是“卡”麻掸。square表示平方，因此在英語中赐纱，卡方分布寫作 $\chi^2$ distribution脊奋。

卡方分布

在理解卡方檢驗(yàn)之前，應(yīng)當(dāng)理解卡方分布疙描〕舷叮卡方分布是一種連續(xù)概率分布。
如果一個(gè)隨機(jī)變量 $Z$ 服從標(biāo)準(zhǔn)正態(tài)分布起胰，即 $Z\sim N(0,1)$ 久又，那么 $Z^2$ 就服從自由度為1的卡方分布。記作 $Z^2 \sim \chi^2_1$ 或者 $Z^2 \sim \chi^2(1)$
而如果 $Z_1,Z_2,\cdots , Z_k$ 都服從標(biāo)準(zhǔn)正態(tài)分布效五，那么它們的平方和服從自由度為 $k$ 的卡方分布地消，記作：
$Z_1^2+Z_2^2+\cdots +Z_k^2 \sim \chi^2_k$
或者寫作 $Z_1^2+Z_2^2+\cdots +Z_k^2 \sim \chi^2(k)$ 。
對于非負(fù)自變量 $x$ 的自由度為 $k$ 的卡方分布的概率密度函數(shù) (簡稱"pdf")：
$f(x)=\dfrac{x^{k/2-1}e^{-x/2}}{2^{k/2}\Gamma(k/2)}$
（1）為什么 $x$ 非負(fù)畏妖？因?yàn)楦鶕?jù)定義脉执，卡方分布的自變量是一個(gè)平方和。
（2）這里的 $\Gamma$ 是一個(gè)函數(shù)戒劫。關(guān)于這個(gè)函數(shù)具體是什么半夷，以及上門的概率密度函數(shù)如何推導(dǎo)，這里不展開迅细，只需要知道有這么個(gè)函數(shù)即可。實(shí)在是好奇的，可以參考這里统阿。
（3）卡方分布的均值為 $\mu=k$ 彩倚，而標(biāo)準(zhǔn)差為 $\sigma=2k$ 。
（4）自由度越大盯质，該函數(shù)圖像越對稱。
（5）為什么這里 $Z$ 需要正態(tài)分布，我的理解是，如果零假設(shè)為真，那么觀測值和期望值之間的波動(dòng)程度喉祭，應(yīng)該是正態(tài)分布的蔽氨，或者說“噪聲”應(yīng)該是正態(tài)分布的鸯匹。

卡方檢驗(yàn)

卡方檢驗(yàn)有兩個(gè)用途：
擬合優(yōu)度檢驗(yàn) chi-squared test goodness of fit
獨(dú)立性檢驗(yàn) chi-squared test of independence

擬合優(yōu)度檢驗(yàn)

某新聞?wù)f某個(gè)籃球明星的原地兩連投的單次命中率是0.8染厅，根據(jù)歷次比賽的數(shù)據(jù)匯總得到下面的表格：

一次兩連投中的命中次數(shù)	0	1	2
觀察次數(shù)	5	82	251

意思是說涩馆，在比賽中，有5次兩連投是一次都沒中，有82次是在兩連投中命中1次。現(xiàn)在，我們來用卡方檢驗(yàn)驗(yàn)證新聞?wù)f的0.8的命中率是否正確。零假設(shè)如下：

$H_0$ ：兩連投的成功次數(shù)符合二項(xiàng)分布庆寺，且概率為 $p=0.8$

（1）先根據(jù)零假設(shè)計(jì)算“期望”的命中次數(shù)分布：
由于總的觀察次數(shù)為 $5+82+251=338$ 琅轧，于是在 $p=0.8$ 成立的前提下，可以計(jì)算每種兩連投結(jié)果的期望次數(shù)：
0次命中: $338 \times C_2^00.8^0 (1-0.8)^2=13.52$
1次命中: $338 \times C_2^10.8^1 (1-0.8)^1=108.16$
2次命中: $338 \times C_2^20.8^2 (1-0.8)^0=216.32$

一次兩連投中的命中次數(shù)	0	1	2
實(shí)際的觀察次數(shù) (observed)	5	82	251
期望的觀察次數(shù) (expected)	13.52	108.16	216.32

顯然憋沿，期望的觀察次數(shù)和實(shí)際的觀察次數(shù)是有偏差的，那么問題在于這個(gè)偏差是否大到具有統(tǒng)計(jì)顯著性甥绿，進(jìn)而可以否定零假設(shè)则披。

（2）我們來構(gòu)造卡方檢驗(yàn)統(tǒng)計(jì)量（chi-squared test statistic）：
$\chi^2=\sum_{cell}\dfrac{(observed-expected)^2}{expected}$
這個(gè)值是把表里每個(gè)格子的實(shí)際值和期望值進(jìn)行對比承璃。為什么要用平方程癌？目的在于規(guī)避正負(fù)號的影響舷嗡。為什么要除以期望值？目的在于消除數(shù)量絕對值的影響锐峭。例如你預(yù)算3塊錢的水中鼠，商家加價(jià)50元，那么這個(gè)波動(dòng)是你無法忍受的沿癞，而你預(yù)算20萬的車援雇，商家加價(jià)50元，則變得可以忍受抛寝。也就是說熊杨，除以期望值目的在于聚焦于變化率曙旭，而不是變化量盗舰。
之后，把這些“變化率”加總得到 $\chi^2$ 桂躏。而計(jì)算自由度有一個(gè)公式：
$DF = (R-1)(C-1)$
其中 R 表示行數(shù)钻趋，C 表示列數(shù)。對于本例：
$DF=(2-1)(3-1)=2$
從另一個(gè)角度解釋為什么 $DF=2$ ：前面的定義是如果是 $k$ 個(gè)符合標(biāo)準(zhǔn)正態(tài)分布的 $Z_i^2$ 相加剂习，則自由度是 $k$ 蛮位，但是這里自有兩個(gè)格子可以自由變化较沪，第三個(gè)格子可以用總觀察數(shù)減出來，例如 $338-13.52-108.16=216.32$ 失仁。
因此尸曼，真正自由的只有2個(gè)格子，所以自由度是2萄焦。

好了控轿，將格子的數(shù)據(jù)代入，求出檢驗(yàn)統(tǒng)計(jì)量：
$\chi^2=\dfrac{(5-13.52))^2}{13.52}+\dfrac{(82-108.16)^2}{108.16}+\dfrac{(251-216.32)^2}{216.32}=17.26$

（3）根據(jù)自由度為2的卡方分布拂封，找到檢驗(yàn)統(tǒng)計(jì)量對應(yīng)的位置：

DF=2的卡方分布

不難理解茬射，隨著統(tǒng)計(jì)量增大，表示預(yù)期的分布和實(shí)際的分布的差異也就越來越大冒签。
另外在抛，由于通常意義上，p值是越小越能推翻零假設(shè)萧恕，那么顯然我們需要用右側(cè)的面積來表示p值刚梭，這里用Python計(jì)算來代替查表：

from scipy.stats import chisquare

observed = [5, 82, 251] 
expected = [13.52, 108.16, 216.32]

output = chisquare(observed, expected)
print('statistic: {:.2f}, pvalue: {:.4f}'.format(output.statistic, output.pvalue))

輸出：statistic: 17.26, pvalue: 0.0002
由于p值很小（假設(shè)我們的顯著性水平的0.05）票唆，那么我們可以推翻零假設(shè)望浩。

進(jìn)一步的，我們來探索下惰说，該運(yùn)動(dòng)員的兩連投的成功次數(shù)分?jǐn)?shù)是否真的符合二項(xiàng)分布磨德。零假設(shè)：

$H_0$ ：兩連投的成功次數(shù)符合二項(xiàng)分布。

既然符合二項(xiàng)分布吆视，那么我們需要先估算一下最合理的 $p$ 概率典挑，那當(dāng)然是用總命中數(shù)除以總投籃數(shù)來計(jì)算了：
$\hat{p}=\dfrac{0\times5+1\times82+2\times252}{2\times(5+82+252)}=0.86$
然后，用該概率值重復(fù)之前的計(jì)算啦吧，也就是先計(jì)算出一個(gè)期望的表格：

一次兩連投中的命中次數(shù)	0	1	2
實(shí)際的觀察次數(shù) (observed)	5	82	251
期望的觀察次數(shù) (expected)	6.26	79.48	252.26

注意您觉，這里的 $DF\ne 2$ ，這是因?yàn)槭谧遥覀兠繌臄?shù)據(jù)估計(jì)一個(gè)參數(shù)琳水，那么我們就損失一個(gè)自由度。這里用了一個(gè)平均命中的概率般堆，因此自由度只有 $DF=(R-1)(C-1)-1=1$ 在孝。

這時(shí)候，在使用 Python 進(jìn)行計(jì)算時(shí)淮摔，注意調(diào)整默認(rèn)的自由度：

# 如果根據(jù)數(shù)據(jù)推斷均值私沮，則損失一個(gè)額外的自由度，因此 ddof = 1
# 自由度 = number of cells - 1 - ddof
observed = [5, 82, 251] 
expected = [6.26, 79.48, 252.26]
output = chisquare(observed, expected, ddof=1)
print(output.statistic, output.pvalue)

這里的 ddof 就是額外損失的自由度和橙，本意是“delta degree of freedom”
輸出：statistic: 0.34, pvalue: 0.56
可以看到p值很大仔燕，因此不足以推翻零假設(shè)造垛，也就是說該運(yùn)動(dòng)員的投籃命中次數(shù)可能真的是二項(xiàng)分布。

獨(dú)立性檢驗(yàn)

下面表格表示喝酒頻率和與警察發(fā)生麻煩的頻數(shù)晰搀。
以第一列為例五辽，表示從不喝酒的人中，4992人不發(fā)生麻煩外恕，71人會(huì)發(fā)生麻煩奔脐。
現(xiàn)在的問題是，能否從以下數(shù)據(jù)推斷說喝酒頻率和與警察發(fā)生麻煩這兩個(gè)事件相互獨(dú)立吁讨？

	Never	Occasional	Frequent	Row total
Trouble with police	71	154	398	623
No trouble with police	4992	2808	2737	10537
Column total	5063	2962	3135	11160

我們的零假設(shè)應(yīng)該如何設(shè)計(jì)髓迎？如果要說明兩者相互獨(dú)立，那么上表的分布應(yīng)該滿足乘法公式建丧。也就是說兩個(gè)獨(dú)立事件一起發(fā)生的概率等于分別發(fā)生的概率之積排龄。
于是我們有：
$p(trouble)=\dfrac{623}{11160}$ 發(fā)生麻煩的總?cè)藬?shù)除以總?cè)藬?shù)
$p(never)=\dfrac{5063}{11160}$ 不喝酒的總?cè)藬?shù)除以總?cè)藬?shù)
進(jìn)一步，根據(jù)總?cè)藬?shù)算出不喝酒而發(fā)生麻煩的人數(shù)的期望（下標(biāo)表示零假設(shè)）：
$E_0(trouble \cap never)=p(trouble)\cdot p(never)\cdot 11160=282.6$
用類似的算法翎朱，計(jì)算每一個(gè)格子在零假設(shè)成立的情況下的值橄维，寫在原表數(shù)據(jù)下的括號里：

	Never	Occasional	Frequent	Row total
Trouble with police	71 (282.6)	154 (165.4)	398 (175.0)	623
No trouble with police	4992 (4780.4)	2808 (2796.6)	2737 (2960.0)	10537
Column total	5063	2962	3135	11160

仔細(xì)觀察可以看出，其實(shí)每個(gè)格子就是對應(yīng)的：
$p = (row\_total \times column\_total) / total$
另外可以看到拴曲，零假設(shè)下的各個(gè)格子的行列之和與原來相同争舞。這不是偶然的，我們用字母代替計(jì)算一下就知道了：

	Never	Occasional	Frequent	Row total
Trouble with police	a	b	c	a+b+c
No trouble with police	d	e	f	d+e+f
Column total	a+d	b+e	c+f	a+b+c+d+e+f

于是第一列的兩個(gè)格子應(yīng)該是：
$\begin{aligned}E_0(trouble \cap never)&=\dfrac{(a+d)(a+b+c)}{a+b+c+d+e+f}\cdot\dfrac{(a+d)(d+e+f)}{a+b+c+d+e+f} \\ &=\dfrac{(a+d)(a+b+c+d+e+f)}{a+b+c+d+e+f}=a+d\end{aligned}$
對于其他格子澈灼、行的總和竞川，都一樣，這里不多說了叁熔。

好委乌，繼續(xù)分析趾盐。我們直接用上表計(jì)算卡方統(tǒng)計(jì)量和p值：

from scipy.stats import chisquare

observed = [71, 154, 398, 4992, 2808, 2737]
expected = [282.6, 165.4, 175.0, 4780.4, 2796.6, 2960.0]

output = chisquare(observed, expected)
print('statistic: {:.2f}, pvalue: {:.4f}'.format(output.statistic, output.pvalue))

這部分計(jì)算方法和擬合優(yōu)度是一樣的妥泉，就不贅述了。計(jì)算發(fā)現(xiàn)這個(gè)p值非常小网沾，接近0心软，因此我們可以推翻零假設(shè)壕吹。也就是說，喝酒的頻率和被警察找麻煩的并不是獨(dú)立的删铃，而是相關(guān)的耳贬。

費(fèi)希爾精確檢驗(yàn) Fisher's exact test

關(guān)于獨(dú)立性檢驗(yàn)，有一個(gè)比卡方檢驗(yàn)更精準(zhǔn)的檢驗(yàn)泳姐，叫 fisher's exact test效拭。它通過直接計(jì)算否定零假設(shè)的概率，也就直接得到了一個(gè)準(zhǔn)確的p值胖秒。有一個(gè)經(jīng)典的女士品茶的統(tǒng)計(jì)學(xué)故事缎患。有一個(gè)女士號稱可以區(qū)分出一杯茶是先倒入了奶還是先倒入了茶。統(tǒng)計(jì)學(xué)家 Fisher 為了驗(yàn)證她的說法阎肝，做了一個(gè)實(shí)驗(yàn)挤渔。拿了8杯茶，4杯是先茶后奶风题，4杯是先奶后茶判导。

	女士判斷先奶后茶	女士判斷先茶后奶	總數(shù)
實(shí)際先奶后茶	4	0	4
實(shí)際先茶后奶	0	4	4
總數(shù)	4	4	8

實(shí)驗(yàn)結(jié)果是全部說對了。那么問題是沛硅，這是否具有統(tǒng)計(jì)顯著性呢眼刃？比如說一個(gè)人猜對了一次硬幣，他的預(yù)測能力靠譜嗎摇肌？
我們假設(shè)女士的判斷是完全隨機(jī)的擂红，這個(gè)是我們的零假設(shè)。那么8杯里面抽中4杯全對的概率是：
$p=\dfrac{C_4^4}{C_8^4}=\dfrac{1}{70}=0.014$
如果顯著性水平是0.01围小，那么我們不能推翻零假設(shè)昵骤，即不敢確定這位女士真有這個(gè)識(shí)別能力。如果顯著性水平定在0.05肯适，則我們可以認(rèn)為她確實(shí)有這個(gè)識(shí)別能力变秦。
如果讓實(shí)驗(yàn)結(jié)果有更大的說服力呢？一個(gè)簡單的辦法就是增加茶的數(shù)量框舔，比如我們設(shè)定為兩種茶各10杯蹦玫，要求10杯都判斷正確，那么p值為多少呢刘绣？
這個(gè)算起來比較麻煩钳垮，這里我寫一個(gè) python 腳本來計(jì)算：

def factorial(n:int):
    """階乘函數(shù)"""
    s = 1
    while n > 1:
        s *= n
        n -= 1
    return s

def C(m:int, n:int):
    """組合數(shù)公式。C(m, n)表示$C_m^n$"""
    return factorial(m) / (factorial(m-n)*factorial(n))

print('p-value = {:.4f}%'.format(C(10,10)/C(20,10)*100))

計(jì)算結(jié)果：
$p=\dfrac{C_{10}^{10}}{C_{20}^{10}}=\dfrac{1}{184756}=0.0005\%$
這個(gè)p值就小得很夸張了额港，基本可以斷定零假設(shè)不成立了饺窿。

那么，回到實(shí)驗(yàn)本身移斩，如果女士只選對了三杯肚医，那么在零假設(shè)的前提下，這個(gè)發(fā)生的概率是多少向瓷？
$p=\dfrac{C_4^3C_4^1}{C_8^4}=\dfrac{16}{70}=22.9\%$
這個(gè)概率比較大了肠套，原大于通常使用的顯著性水平 0.05，因此我們沒有辦法推翻零假設(shè)猖任。為什么要這樣 $C_4^3$ 乘 $C_4^1$ 呢你稚？這個(gè)是因?yàn)橐还灿羞@么多種取法。你把所有可能的取法羅列處理，就是16種刁赖，然后除以總的取法數(shù)搁痛，就是隨機(jī)取到這樣結(jié)果的概率。

比較 Fisher's exact test 和 chi-squared test宇弛，可以參考這篇文章鸡典。

條件	chi-squared test	Fisher's exact test
最小樣本量	大	小
準(zhǔn)確度	近似	精確
列聯(lián)表	任意維度	一般是2x2
解釋	pearson residuals	odds ratio

一般來說，兩者都適用的情況下枪芒，應(yīng)該優(yōu)先選擇 Fisher's exact test彻况，因?yàn)樗蔷_值。如果實(shí)驗(yàn)觀察的數(shù)量很芯俗佟（小于10）纽甘，應(yīng)該不使用 chi-squared test。

下面使用一個(gè)腳本來計(jì)算：

from scipy.stats import fisher_exact

# 列聯(lián)表
contigency_table = [[0, 4], [4, 0]]

# alternative 表示備擇假設(shè)
# alternative='less' 表示零假設(shè)是相等抽碌，備擇假設(shè)是
oddsratio, pvalue = fisher_exact(contigency_table, alternative='less')
print(pvalue)
# 輸出：0.014285714285714268

參考

An Introduction to the Chi-Square Distribution
Chi-square tests: Goodness of Fit for the Binomial Distribution
Chi-square Tests of Independence (Chi-square Tests for Two-Way Tables)
Python Programming: scipy.stats.chisquare
Python Programming: scipy.stats.fisher_exact