2×2 分割表のフィッシャーの正確確率検定     Last modified: Jun 28, 2004

目的

2×2 分割表において,フィッシャーの正確確率検定を行う。
R には fisher.test があるが,この関数では fisher.test が計算する Fisher の方法による P 値と同時に,Pearson の方法による P 値も計算する(オッズ比を基準にするものも計算する)。

使用法

my.fisher(x)

引数

x      2×2 分割表を表す,2×2 行列
 
ソース

インストールは,以下の 1 行をコピーし,R コンソールにペーストする
source("http://aoki2.si.gunma-u.ac.jp/R/src/my_fisher.R", encoding="euc-jp")

# 2×2 分割表のフィッシャーの正確確率検定
my.fisher <- function(x)                                     # 2×2 分割表
{
        odds.ratio <- function(a, b, c, d)                   # オッズ比を求める関数
        {
                if (a*b*c*d == 0) {                             # セルのどれかが 0 の場合
                        a <- a+0.5                           # それぞれに 0.5 を加える
                        b <- b+0.5
                        c <- c+0.5
                        d <- d+0.5
                }
                res <- a*d/(b*c)                             # オッズ比
                return(max(res, 1/res))                         # どちらか大きい方を返す
        }

        stats <- function(i) {
                e <- (a <- i[1]) + (b <- i[2])
                f <- (c <- i[3]) + (d <- i[4])
                n <- (g <- a+c) + (h <- b+d)
                return(c(n*(a*d-b*c)^2/(e*f*g*h),               # カイ二乗値
                         odds.ratio(a, b, c, d),                # オッズ比
                         exp(lchoose(e, a)+lchoose(f, c)-lchoose(n, g))))       # 生起確率
        }

        ct <- colSums(x)                                     # 列和
        rt <- rowSums(x)                                     # 行和
        n <- sum(x)                                          # 総合計
        mx <- min(rt[1], ct[1])                                      # a が取り得る最大値
        mi <- max(0, rt[1]+ct[1]-n)                          # a が取り得る最小値
        A <- mi:mx                                           # a のベクトル
        B <- rt[1]-A                                         # b のベクトル
        C <- ct[1]-A                                         # c のベクトル
        D <- ct[2]-B                                         # d のベクトル
        Cell <- cbind(A, B, C, D)                            # 行方向の 4 つの数値が一つの分割表
        result<-apply(Cell, 1, stats)                                # 行方向の 4 つの数値に,stats 関数を適用
        Chi.sq <- result[1,]                                 # カイ二乗値ベクトル
        OddsRatio <- result[2,]                                      # オッズ比ベクトル
        Probability <- result[3,]                            # 生起確率ベクトル
        Cum.Prob1 <- cumsum(Probability)                     # 生起確率の累積和
        Cum.Prob2 <- rev(cumsum(rev(Probability)))           # 逆方向からの生起確率の累積和
        Pearson <- Chi.sq >= stats(c(x))[1]-1e-15            # カイ二乗値から判定した P 値に算入すべきセル
        Fisher <- Probability <= stats(c(x))[3]+1e-15             # フィッシャーの定義による P 値に算入すべきセル
        OR <- OddsRatio >= stats(c(x))[2]-1e-15                      # オッズ比から判定した P 値に算入すべきセル
        p.Pearson <- sum(Probability[Pearson])                       # カイ二乗値から判定した P 値
        p.Fisher <- sum(Probability[Fisher])                 # フィッシャーの定義による P 値
        p.OR <- sum(Probability[OR])                         # オッズ比から判定した P 値
        return(list(result=cbind(Cell, Chi.sq, Probability, OddsRatio, Pearson, Fisher, OR, Cum.Prob1, Cum.Prob2), p.Pearson=p.Pearson, p.Fisher=p.Fisher, p.OR=p.OR))
}


使用例

my.fisher(matrix(c(10, 13, 16, 61), nc=2, byrow=TRUE))

出力結果例

> my.fisher(matrix(c(10, 13, 16, 61), nc=2, byrow=TRUE))

# 一行がそれぞれ,可能な2×2分割表
# A, B, C, D は各セルの値
# Chi.sq は,カイ二乗値(連続性の修正なし)
# Probability は,分割表の生起確率
# Pearson は,カイ二乗値が観察された分割表のカイ二乗値以上の場合に 1 を取る(Pearson の基準による)
# Fisher は,生起確率が観察された分割表の生起確率以下の場合に 1 を取る(Fisher の基準による)
# Cum.Prob1, Cum.Prob2 は両端からの累積確率(片側検定の時に利用できる)

$result
       A  B  C  D       Chi.sq  Probability   OddsRatio Pearson Fisher OR    Cum.Prob1    Cum.Prob2
 [1,]  0 23 26 51 1.049491e+01 3.317552e-04   24.184466       1      1  1 0.0003317552 1.000000e+00
 [2,]  1 22 25 52 7.278386e+00 3.815185e-03   10.576923       1      1  1 0.0041469406 9.996682e-01
 [3,]  2 21 24 53 4.648818e+00 1.979577e-02    4.754717       0      1  1 0.0239427135 9.958531e-01
 [4,]  3 20 23 54 2.606207e+00 6.158685e-02    2.839506       0      0  0 0.0855295626 9.760573e-01
 [5,]  4 19 22 55 1.150553e+00 1.287725e-01    1.900000       0      0  0 0.2143020653 9.144704e-01
 [6,]  5 18 21 56 2.818568e-01 1.922390e-01    1.350000       0      0  0 0.4065410157 7.856979e-01
 [7,]  6 17 20 57 1.173914e-04 2.124746e-01    1.005882       0      0  0 0.6190156452 5.934590e-01
 [8,]  7 16 19 58 3.053351e-01 1.779344e-01    1.335526       0      0  0 0.7969500639 3.809844e-01
 [9,]  8 15 18 59 1.197510e+00 1.146018e-01    1.748148       0      0  0 0.9115518929 2.030499e-01
[10,]  9 14 17 60 2.676642e+00 5.730091e-02    2.268908       0      0  0 0.9688528074 8.844811e-02
[11,] 10 13 16 61 4.742731e+00 2.235675e-02    2.932692       1      1  1 0.9912095577 3.114719e-02
[12,] 11 12 15 62 7.395777e+00 6.818481e-03    3.788889       1      1  1 0.9980280387 8.790442e-03
[13,] 12 11 14 63 1.063578e+01 1.623448e-03    4.909091       1      1  1 0.9996514866 1.971961e-03
[14,] 13 10 13 64 1.446274e+01 3.004940e-04    6.400000       1      1  1 0.9999519805 3.485134e-04
[15,] 14  9 12 65 1.887666e+01 4.292771e-05    8.425926       1      1  1 0.9999949082 4.801948e-05
[16,] 15  8 11 66 2.387753e+01 4.683023e-06   11.250000       1      1  1 0.9999995913 5.091772e-06
[17,] 16  7 10 67 2.946536e+01 3.844272e-07   15.314286       1      1  1 0.9999999757 4.087489e-07
[18,] 17  6  9 68 3.564015e+01 2.327847e-08   21.407407       1      1  1 0.9999999990 2.432162e-08
[19,] 18  5  8 69 4.240190e+01 1.012107e-09   31.050000       1      1  1 1.0000000000 1.043154e-09
[20,] 19  4  7 70 4.975060e+01 3.043931e-11   47.500000       1      1  1 1.0000000000 3.104672e-11
[21,] 20  3  6 71 5.768626e+01 6.002118e-13   78.888889       1      1  1 1.0000000000 6.074018e-13
[22,] 21  2  5 72 6.620888e+01 7.145379e-15  151.200000       1      1  1 1.0000000000 7.189975e-15
[23,] 22  1  4 73 7.531845e+01 4.449177e-17  401.500000       1      1  1 1.0000000000 4.459634e-17
[24,] 23  0  3 74 8.501499e+01 1.045635e-19 1000.428571       1      1  1 1.0000000000 1.045635e-19

$p.Pearson	# Pearson の基準による P 値(上記 Pearson の列が 1 になっている生起確率の和)
[1] 0.03529413

$p.Fisher	# Fisher の基準による P 値(上記 Fisher の列が 1 になっている生起確率の和)
[1] 0.05508991

$p.OR	# オッズ比を基準としたときの P 値(上記 OR の列が 1 になっている生起確率の和)
[1] 0.05508991

# 上記三つの P 値はほとんどの場合に同じになるが,上述の例のように違う結果になることもままある。
# R の fisher.test は,当然ながら(?) Fisher の基準によるものである。
> fisher.test(matrix(c(10, 13, 16, 61), nc=2, byrow=TRUE))

	Fisher's Exact Test for Count Data

data:  matrix(c(10, 13, 16, 61), nc = 2, byrow = TRUE) 
p-value = 0.05509
  以下略

・ 解説ページ


・ 直前のページへ戻る  ・ E-mail to Shigenobu AOKI

Made with Macintosh