前回に引き続き、統計の備忘録。
例えば、カテゴリカル変数のXとYというデータの相関、オッズ比(超幾何分布による最尤推定値)、 95% 信頼限界などは、カイ二乗検定 odds.ratio()
やフィッシャーの正確率検定 fisher.test()
などでオッズ比(超幾何分布による最尤推定値)や 95% 信頼限界を求めることができるが、
下記のようなクロス集計表があった場合、
X1 | X2 | Total | |
Y1 | a | b | a+b |
Y2 | c | d | c+d |
Total | a+c | b+d | a+b+c+d |
オッズ比と信頼区間は、
odds.ratio()
でも求める事ができる。
また、これくらいの計算なら、エクセルでも可能 ▼
オッズ比とP値はエクセルで計算可。 X1 X2 Total Y1 a b a+b Y2 c d c+d Total a+c b+d a+b+c+d 上記のようなデータがあった場合、期待値は、 期待値 X1 X2 Tota …
今回はまず、odds.ratio()
を使って、オッズ比と95%信頼区間を求める方法。
オッズ比と95%信頼区間を求める
基本的には、青木先生のサイトで完結。
オッズ比とその信頼限界を求める。 R では,fisher.test 関数で,オッズ比(超幾何分布による最尤推定値)と 95% 信頼限界を求めることができる。 また,vcd ライブラリーにも oddsratio 関数が用意されている(各セルの数値に必ず 0.5 が加えられる。log=FALSE 引数をつけ,summary 関数を適用すること)。
まず、下記プログラムをRコンソールにペースト。
# オッズ比を計算する
odds.ratio <- function(a, b, c, d, # 4 つのセルの観察値
correct=FALSE) # 補正をするとき TRUE にする
{
cl <- function(x)
{
or*exp(c(1,-1)*qnorm(x)*sqrt(1/a+1/b+1/c+1/d)) # OR の計算
}
if (correct || a*b*c*d==0) { # a,b,c,d のいずれかが 0 のときには必ず補正する
a <- a+0.5
b <- b+0.5
c <- c+0.5
d <- d+0.5
}
or <- a*d/(b*c)
conf <- rbind(cl90=cl(0.05), cl95=cl(0.025), cl99=cl(0.005), cl999=cl(0.0005))
conf <- data.frame(conf)
colnames(conf) <- paste(c("下側","上側"), "信頼限界値" , sep="" )
rownames(conf) <- paste(c(90, 95, 99, 99.9), "%信頼区間" , sep="" )
list(or=or, conf=conf)
}
そして、上表の a, b, c, d の値を、下記コード a, b, c, d に当てはめて使用する。
補正なしの場合
odds.ratio(a, b, c, d)
# クロス集計表の a, b, c, d を当てはめる。例えば a=3, b=5, c=6, d=8 の場合、
odds.ratio(3, 5, 6, 8)
output
$or # オッズ比
[1] 0.8
$conf
下側信頼限界値 上側信頼限界値
90%信頼区間 0.17957613 3.563948
95%信頼区間 0.13488007 4.744956
99%信頼区間 0.07709193 8.301777
99.9%信頼区間 0.04027949 15.888980
上の例だと、
- オッズ比は0.8
- 95%信頼区間は、0.13-4.74
95%信頼区間の上側信頼限界値と下側信頼限界値が1をはさんでいるので、今回は有意差なし。
補正ありの場合
odds.ratio(a, b, c, d)
# クロス集計表の a, b, c, d を当てはめる。補正ありで、0.5を加える
odds.ratio(3, 5, 6, 8), correct=TRUE
output
$or # オッズ比
[1] 0.8321678
$conf
下側信頼限界値 上側信頼限界値
90%信頼区間 0.20235212 3.422269
95%信頼区間 0.15433409 4.487040
99%信頼区間 0.09089284 7.618899
99.9%信頼区間 0.04916982 14.083909
上の例だと、
- オッズ比は0.8321678
- 95%信頼区間は、0.15-4.48
補正後も、95%信頼区間の上側信頼限界値と下側信頼限界値が1をはさんでいるので、有意差なし。
Fisherの正確検定で求める方法
1.上記のクロス集計表をマトリックス表示する。
# マトリックスを作る
x < matrix(c(3,5,6,8), byrow=TRUE, nc=2)
x
output
x < matrix(c(3,5,6,8), byrow=TRUE, nc=2)
x
[,1] [,2]
[1,] 3 5
[2,] 6 8
2.R標準パッケージの fisher.test()
を使用する。
# マトリックスXを使って、フィッシャーの正確率検定
fisher.test(x)
output
> fisher.test(x)
Fisher's Exact Test for Count Data
data: x
p-value = 1
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
0.08839317 6.37801297
sample estimates:
odds ratio
0.8081216
結果は、
- オッズ比: 0.8
- P値: 1
- 95%信頼区間: 0.09-6.39
こちらも有意差なし。
References
オッズ比とその信頼限界を求める。 R では,fisher.test 関数で,オッズ比(超幾何分布による最尤推定値)と 95% 信頼限界を求めることができる。 また,vcd ライブラリーにも oddsratio 関数が用意されている(各セルの数値に必ず 0.5 が加えられる。log=FALSE 引数をつけ,summary 関数を適用すること)。
今回はR言語でフィッシャーの正確率検定を実行する方法をみていきます。前回紹介したカイ2乗検定の引き続き、カテゴリカルデータの検定方法であるフィッシャーの正確率検定の実行例を紹介します。
Rによるフィッシャーの直接確率検定の実行方法について.