分割表

statsmodels は、独立性、対称性、均一性を評価する方法、層別母集団からの表のコレクションを操作する方法など、分割表を分析するためのさまざまなアプローチをサポートしています。

ここで説明するメソッドは主に双方向テーブル用です。多元表は対数線形モデルを使用して分析できます。 statsmodels には現在、対数線形モデリング用の専用 API がありませんが、ポアソン回帰を statsmodels.genmod.GLM この目的に使用できます。

分割表は、各観測値が複数の変数ごとに 1 つのカテゴリに属する​​データ セットを記述する多元表です。たとえば、2つの変数があり、1つは \(r\) レベルで、もう1つは \(c\) レベルである場合、 \(r\times c\) 分割表があります。テーブルは、テーブルの特定のセルに入る観測値の数で記述できます。たとえば、 \(T_{ij}\) は、最初の変数に対してレベル \(i\) を持ち、2番目の変数に対してレベル \(j\) を持つ観測値の数です。各変数は、順序付きまたは順序なしのいずれかである有限数のレベル(またはカテゴリ)を持たなければならないことに注意してください。異なるコンテキストでは、分割表の軸を定義する変数は、 カテゴリ変数 または 要素変数 と呼ばれることがあります。それらは 公称 (それらのレベルが順序付けされていない場合)または 序数 (それらのレベルが順序付けされている場合)のいずれかです。

分割表の基礎となる母集団は、 分布テーブル \(P_{i, j}\) によって記述されます。 \(P\) の要素は確率であり、 \(P\) のすべての要素の合計は1です。分割表を分析するメソッドは、 \(T\) のデータを使用して \(P\) のプロパティについて学習します。

statsmodels.stats.Table は、分割表を扱うための最も基本的なクラスです。分割表セルカウントを含む矩形状の配列のようなオブジェクトから直接 Table オブジェクトを作ることができます:

In [1]: import numpy as np

In [2]: import pandas as pd

In [3]: import statsmodels.api as sm

In [4]: df = sm.datasets.get_rdataset("Arthritis", "vcd").data

In [5]: df.fillna({"Improved":"None"}, inplace=True)

In [6]: tab = pd.crosstab(df['Treatment'], df['Improved'])

In [7]: tab = tab.loc[:, ["None", "Some", "Marked"]]

In [8]: table = sm.stats.Table(tab)

あるいは、生データを渡して、Table クラスにセル数の配列を構築させることもできます:

In [9]: data = df[["Treatment", "Improved"]]

In [10]: table = sm.stats.Table.from_data(data)

独立性

独立性 は、行と列の要素が独立して発生するプロパティです。 関連性 は独立性の欠如です。結合分布が独立している場合は、行と列の周辺分布の外積として記述できます:

\[P_{ij} = \sum_k P_{ij} \cdot \sum_k P_{kj} \quad \text{for all} \quad i, j\]

観測されたデータに最適な独立分布を取得し、最も強く独立性に違反している特定のセルを識別する残差を表示できます:

In [11]: print(table.table_orig)
Improved   Marked  None  Some
Treatment                    
Placebo         7    29     7
Treated        21    13     7

In [12]: print(table.fittedvalues)
Improved      Marked  None      Some
Treatment                           
Placebo    14.333333  21.5  7.166667
Treated    13.666667  20.5  6.833333

In [13]: print(table.resid_pearson)
Improved     Marked      None      Some
Treatment                              
Placebo   -1.936992  1.617492 -0.062257
Treated    1.983673 -1.656473  0.063758

この例では、行と列が独立している集団からの標本と比較して、プラセボ/改善なしおよび処置/顕著な改善のセルでは観察が多すぎ、プラセボ/顕著な改善および処置/改善なしのセルでは観察が少なすぎます。これは、処置の明らかな利点を反映しています。

表の行と列が順序付けられていない(つまり、名目上の要素である)場合、独立性を正式に評価するための最も一般的なアプローチは、ピアソンの \(\chi^2\) 統計を使用することです。 \(\chi^2\) 統計へのセルごとの寄与を見て、依存の証拠がどこから来ているのかを知ることはしばしば有用です。

In [14]: rslt = table.test_nominal_association()

In [15]: print(rslt.pvalue)
0.0014626434089526352

In [16]: print(table.chi2_contribs)
Improved     Marked      None      Some
Treatment                              
Placebo    3.751938  2.616279  0.003876
Treated    3.934959  2.743902  0.004065

順序付けられた行と列の係数を持つ表の場合、順序を尊重する対立仮説に対してより多くの検出力を得るために、 線形×線形 関連検定を使用できます。線形×線形関連検定の検定統計量は次のとおりです

\[\sum_k r_i c_j T_{ij}\]

ここで \(r_i\)\(c_j\) は行と列のスコアです。多くの場合、これらのスコアは0,1,... の順に設定されます。これが 'コクラン・アーミテージのトレンド検定' です。

In [17]: rslt = table.test_ordinal_association()

In [18]: print(rslt.pvalue)
0.023644578093923983

一連の \(2\times 2\) テーブルを作成し、それらのオッズ比を計算することで、 \(r\times x\) テーブル内の関連を評価することができます。これには2つの方法があります。 ローカルオッズ比 は、隣接する行と列のカテゴリから \(2\times 2\) テーブルを作成します。

In [19]: print(table.local_oddsratios)
Improved     Marked      None  Some
Treatment                          
Placebo    0.149425  2.230769   NaN
Treated         NaN       NaN   NaN

In [20]: taloc = sm.stats.Table2x2(np.asarray([[7, 29], [21, 13]]))

In [21]: print(taloc.oddsratio)
0.14942528735632185

In [22]: taloc = sm.stats.Table2x2(np.asarray([[29, 7], [13, 7]]))

In [23]: print(taloc.oddsratio)
2.230769230769231

累積オッズ比 は、可能な各ポイントで行と列の要素を二分することによって、 \(2\times 2\) テーブルを構成します。

In [24]: print(table.cumulative_oddsratios)
Improved     Marked      None  Some
Treatment                          
Placebo    0.185185  1.058824   NaN
Treated         NaN       NaN   NaN

In [25]: tab1 = np.asarray([[7, 29 + 7], [21, 13 + 7]])

In [26]: tacum = sm.stats.Table2x2(tab1)

In [27]: print(tacum.oddsratio)
0.18518518518518517

In [28]: tab1 = np.asarray([[7 + 29, 7], [21 + 13, 7]])

In [29]: tacum = sm.stats.Table2x2(tab1)

In [30]: print(tacum.oddsratio)
1.0588235294117647

モザイク プロットは、二元テーブルの依存関係を非公式に評価するためのグラフィカルなアプローチです。

In [31]: from statsmodels.graphics.mosaicplot import mosaic

In [32]: fig, _ = mosaic(data, index=["Treatment", "Improved"])

対称性と均一性

対称性 とは、すべての \(i\)\(j\) に対して \(P_{i, j}=P_{j, i}\) であるという性質です。 均質性 とは、行係数と列係数の周辺分布が同一であるという性質です。つまり

\[\sum_j P_{ij} = \sum_j P_{ji} \forall i\]

これらのプロパティを適用するためには、テーブル \(P\) (および \(T\) )は正方形でなければならず、行と列のカテゴリは同一で、同じ順序でなければならないことに注意してください。

説明のために、データセットをロードし、分割表を作成し、行と列のマージンを計算します。 Table クラスには、 \(r\times c\) 分割表を分析するためのメソッドが含まれています。以下にロードされたデータセットには、人の左目と右目の視力の評価が含まれています。まずデータをロードし、分割表を作成します。

In [33]: df = sm.datasets.get_rdataset("VisualAcuity", "vcd").data

In [34]: df = df.loc[df.gender == "female", :]

In [35]: tab = df.set_index(['left', 'right'])

In [36]: del tab["gender"]

In [37]: tab = tab.unstack()

In [38]: tab.columns = tab.columns.get_level_values(1)

In [39]: print(tab)
right     1     2     3    4
left                        
1      1520   234   117   36
2       266  1512   362   82
3       124   432  1772  179
4        66    78   205  492

次に、分割表から SquareTable オブジェクトを作成します。

In [40]: sqtab = sm.stats.SquareTable(tab)

In [41]: row, col = sqtab.marginal_probabilities

In [42]: print(row)
right
1    0.255049
2    0.297178
3    0.335295
4    0.112478
dtype: float64

In [43]: print(col)
right
1    0.264277
2    0.301725
3    0.328474
4    0.105524
dtype: float64

summary メソッドは、対称性と均質性の検定手順の結果を出力します。

In [44]: print(sqtab.summary())
            Statistic P-value DF
--------------------------------
Symmetry       19.107   0.004  6
Homogeneity    11.957   0.008  3
--------------------------------

data という名前のデータフレームに個々の処置記録があれば、SquareTable.from_data クラスメソッドを使って生のデータを渡すことで同じ分析を行うこともできます。

sqtab = sm.stats.SquareTable.from_data(data[['left', 'right']])
print(sqtab.summary())

単一の 2x2 テーブル

sm.stats.Table2x2 クラスには、個々の2x2テーブルを操作するためのいくつかのメソッドが用意されています。 summary メソッドは、テーブルの行と列の間の関連のいくつかの測度を表示します。

In [45]: table = np.asarray([[35, 21], [25, 58]])

In [46]: t22 = sm.stats.Table2x2(table)

In [47]: print(t22.summary())
               Estimate   SE   LCB   UCB  p-value
-------------------------------------------------
Odds ratio        3.867       1.890 7.912   0.000
Log odds ratio    1.352 0.365 0.636 2.068   0.000
Risk ratio        2.075       1.411 3.051   0.000
Log risk ratio    0.730 0.197 0.345 1.115   0.000
-------------------------------------------------

リスク比は対称ではないため、転置された表を分析すると異なる結果が得られることに注意してください。

In [48]: table = np.asarray([[35, 21], [25, 58]])

In [49]: t22 = sm.stats.Table2x2(table.T)

In [50]: print(t22.summary())
               Estimate   SE   LCB   UCB  p-value
-------------------------------------------------
Odds ratio        3.867       1.890 7.912   0.000
Log odds ratio    1.352 0.365 0.636 2.068   0.000
Risk ratio        2.194       1.436 3.354   0.000
Log risk ratio    0.786 0.216 0.362 1.210   0.000
-------------------------------------------------

階層化された 2x2 テーブル

層別化は、同じ行と列の因子によって定義された分割表の集合がある場合に行われます。下の例では、中国のいくつかの地域のそれぞれにおける喫煙と肺がんの同時分布を反映した2x2の表の集合があります。たとえ周辺確率が階層間で異なっていても、すべての表が共通のオッズ比を持っている可能性があります。 'Breslow-Day' 法は、データが共通のオッズ比と一致しているかどうかを検定します。これは以下に 定数 OR検定 として表示されます。Mantel-Haenszel法は、この共通のオッズ比が1に等しいかどうかを検定します。これは以下に OR=1 の検定 として表示されます。共通のオッズ比とリスク比を推定し、それらの信頼区間を得ることもできます。 summary 法は、これらの結果をすべて表示します。個々の結果は、クラスのメソッドと属性から取得できます。

In [51]: data = sm.datasets.china_smoking.load_pandas()

In [52]: mat = np.asarray(data.data)

In [53]: tables = [np.reshape(x.tolist(), (2, 2)) for x in mat]

In [54]: st = sm.stats.StratifiedTable(tables)

In [55]: print(st.summary())
                   Estimate   LCB    UCB 
-----------------------------------------
Pooled odds           2.174   1.984 2.383
Pooled log odds       0.777   0.685 0.868
Pooled risk ratio     1.519              
                                         
                 Statistic P-value 
-----------------------------------
Test of OR=1       280.138   0.000 
Test constant OR     5.200   0.636 
                       
-----------------------
Number of tables    8  
Min n             213  
Max n            2900  
Avg n            1052  
Total n          8419  
-----------------------

モジュールリファレンス

Table(table[, shift_zeros])

双方向の分割表。

Table2x2(table[, shift_zeros])

2x2 分割表で実行できる分析。

SquareTable(table[, shift_zeros])

正方形分割表を分析する方法。

StratifiedTable(tables[, shift_zeros])

2x2 分割表のコレクションを分析します。

mcnemar(table[, exact, correction])

均一性のマクネマー検定。

cochrans_q(x[, return_object])

同一二項比率に対するコクランの Q 検定。

こちらも参照

Scipy には、現在統計モデルには含まれていないフィッシャーの直接確率検定など、分割表を分析するための関数がいくつかあります。


最終更新日: 2025年01月28日