\(2 \times 2\)交叉表是類別資料分析最常見的形式,值得再詳細介紹分析方法。這裡同樣引用Effectiveness of the Pfizer-BioNTech and Oxford-AstraZeneca vaccines on covid-19 related symptoms, hospital admissions, and mortality in older adults in England: test negative case-control study,將Table 4稍作修改作為本部分的分析資料:
Pfizer-BionTech | Oxford-AstraZeneca | ||||
---|---|---|---|---|---|
未住院 | 住院 | 未住院 | 住院 | ||
未接種 | 7527 | 1365 | 未接種 | 7527 | 1365 |
已接種 | 3063 | 421 | 已接種 | 616 | 72 |
卡方檢定與Yates校正
除了table()搭配summary()可用於卡方分析外,chisq.test()是個更方便、更直接的指令,同時可以設定Yates校正。卡方檢定中當儲存格的期望值小於5且樣本大於50時,可以透過Yates校正修正。原則上當觀察值小於(大於)期望值時,則透過觀察值加(減)0.5進行校正。chisq.test()預設會自動進行Yates校正,也可以透過correct=FALSE參數關閉此預設值。
為了說明chisq.test(),我們以上述covid-19疫苗比較表為基礎,創造出三個\(2 \times 2\)交叉表,分別是Pfizer-BionTech的住院交叉表、Oxford-AstraZeneca的住院交叉表、以及混合Pfizer-BionTech與Oxford-AstraZeneca的住院交叉表。我們將用這三個\(2 \times 2\)交叉表來驗證Pfizer-BionTech與Oxford-AstraZeneca疫苗是否能降低重症住院的風險。
> pfizer<-matrix(c(7527,1365,3063,421), ncol=2, byrow=TRUE, dimnames=list(c("Unvaccinated","Vaccinated"), c("Unhospitalized", "Hospitalized")))
> pfizer
Unhospitalized Hospitalized
Unvaccinated 7527 1365
Vaccinated 3063 421
> az<-matrix(c(7527,1365,616,72), ncol=2, byrow=TRUE, dimnames=list(c("Unvaccinated","Vaccinated"), c("Unhospitalized", "Hospitalized")))
> az
Unhospitalized Hospitalized
Unvaccinated 7527 1365
Vaccinated 616 72
> pfizer_az<-matrix(c(3063,421,616,72), ncol=2, byrow=TRUE, dimnames=list(c("Pfizer-BionTech","Oxford-AstraZeneca"), c("Unhospitalized", "Hospitalized")))
> pfizer_az
Unhospitalized Hospitalized
Pfizer-BionTech 3063 421
Oxford-AstraZeneca 616 72
首先用chisq.test()檢定Pfizer-BionTech疫苗與住院重症的卡方分析,可以發現是否接種疫苗與重症住院達到統計顯著,兩者有統計上的差異。
> chisq.test(pfizer)
Pearson's Chi-squared test with Yates' continuity correction
data: pfizer
X-squared = 21.373, df = 1, p-value = 3.78e-06
Oxford-AstraZeneca疫苗的檢定結果,同樣顯示接種疫苗與重症住院兩者有統計上的差異。
> chisq.test(az)
Pearson's Chi-squared test with Yates' continuity correction
data: az
X-squared = 11.576, df = 1, p-value = 0.0006682
接著比較Pfizer-BionTech與Oxford-AstraZeneca兩款疫苗對預防重症住院的效果,卡方檢定未達統計顯著,顯示兩款疫苗對重症住院的預防效果是一樣的。
> chisq.test(pfizer_az) #有Yates校正
Pearson's Chi-squared test with Yates' continuity correction
data: pfizer_az
X-squared = 1.2935, df = 1, p-value = 0.2554
> chisq.test(pfizer_az, correct=FALSE) #無Yates校正
Pearson's Chi-squared test
data: pfizer_az
X-squared = 1.4447, df = 1, p-value = 0.2294
事實上我們可以用prop.table()比較兩款疫苗,可以發現無論是哪一款疫苗,只要接種過就能預防重症住院。其中接種Pfizer-BionTech只有12%的人住院,接種Oxford-AstraZeneca只有10%的人住院。
> prop.table(pfizer_az, 1)
Unhospitalized Hospitalized
Pfizer-BionTech 0.8791619 0.1208381
Oxford-AstraZeneca 0.8953488 0.1046512
卡方檢定與概似比檢定
除了卡方統計量外,概似比檢定也是常見的統計方法。與卡方檢定相同,概似比統計量越大,表示越有機會拒絕虛無假設。透過R vcd套件中的assocstats()可以很容易計算概似比。
> library(vcd)
> assocstats(pfizer_az)
X^2 df P(> X^2)
Likelihood Ratio 1.4855 1 0.22292
Pearson 1.4447 1 0.22939
Phi-Coefficient : 0.019
Contingency Coeff.: 0.019
Cramer's V : 0.019
概似比檢定與卡方檢定的結果相似,統計量僅1.48未達統計顯著,表示Pfizer-BionTech與Oxford-AstraZeneca的差異,沒有大到可以拒絕虛無假設。independence_table()可顯示接受虛無假設的期望值,與實際觀察值進行比較可以發現兩者差異不大,這也是未達統計顯著的原因。
> expected_value<-independence_table(pfizer_az) #期望值
> expected_value
Unhospitalized Hospitalized
Pfizer-BionTech 3072.3001 411.6999
Oxford-AstraZeneca 606.6999 81.3001
> pfizer_az #觀察值
Unhospitalized Hospitalized
Pfizer-BionTech 3063 421
Oxford-AstraZeneca 616 72
> expected_value-pfizer_az #期望值與觀察值的差距
Unhospitalized Hospitalized
Pfizer-BionTech 9.300096 -9.300096
Oxford-AstraZeneca -9.300096 9.300096
費雪精確性檢定(Fisher's exact test)
費雪(Ronald Fisher)在1935年的教科書《實驗設計》(The Design of Experiments)第二章中,描述了一個經典的案例。這個案例起源於一名女同事宣稱,她可以分辨奶茶到底是茶先倒入,還是牛奶先倒入。一個聽起來荒謬的說法,費雪設計了一個「茶實驗」(lady tasting tea)來印證真實性。費雪的實驗設計是這樣子:準備8杯奶茶,其中4杯先倒入茶,另外4杯先倒入牛奶,然後隨機排列讓女同事試喝這8杯奶茶後挑選4杯來猜測,最後實驗結果無論是牛奶或茶都猜中3杯,各有1杯猜錯。把實驗結果整理成交叉表後長得像下面這樣:
牛奶 | 茶 | 合計 | |
牛奶 | 3 | 1 | 1 |
茶 | 1 | 3 | 4 |
合計 | 4 | 4 |
因為這個實驗一共有8杯奶茶要選4杯,所以一共有\(C_{4}^{8}=70\)個組合。猜對的機率會是下面這張表。費雪的女同事猜對3杯,看起來好像真的可以分辨到底是不是牛奶先加入,但其實如果根據機率就會發現,其實4杯猜中1杯、4杯猜中2杯、4杯猜中3杯的機率都還蠻高的。
猜對次數 | 0 | 1 | 2 | 3 | 4 |
猜對機率 | \(\dfrac{C_{4}^{4}}{C_{4}^{8}}=1.4\%\) | \(\dfrac{C_{1}^{4}C_{3}^{4}}{C_{4}^{8}}=22.9\%\) | \(\dfrac{C_{2}^{4}C_{2}^{4}}{C_{4}^{8}}=51.4\%\) | \(\dfrac{C_{3}^{4}C_{1}^{4}}{C_{4}^{8}}=22.9\%\) | \(\dfrac{C_{4}^{4}}{C_{4}^{8}}=1.4\%\) |
這項著名的茶實驗模擬樣本,可以在範例檔案中fishertea.csv下載。其中茶先倒入=1,牛奶先倒入=2,可以看到guess各分別猜錯一次,我們可以用fisher.test()來執行精確性檢定,獲得不拒絕虛無假設的結果:
> fisher<-read.csv("c:/Users/USER/downloads/fishertea.csv", header=T, sep=",")
> fisher
real guess
1 1 1
2 1 1
3 1 1
4 1 2
5 2 2
6 2 2
7 2 1
8 2 2
> attach(fisher)
> tea_test<-table(real, guess)
> fisher.test(tea_test)
Fisher's Exact Test for Count Data
data: tea_test
p-value = 0.4857
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
0.2117329 621.9337505
sample estimates:
odds ratio
6.408309