本欄で説明したデータ解析の方法のほとんどは、量的なデータを対象にしたものである。今月号から数回の誌面を用いてカテゴリカルデータ (categorical data) について説明を行うことにする。
我々の周りには多くの物事の性質を数え上げる計数 (count) データがある。例えば、性別(男、女)、血液型(A,B,O,AB)やアンケート調査票の中の質問における選択肢などを集計したデータは、いずれも計数データである。このような計数データをカテゴリカルデータと呼び、性別や血液型などのカテゴリーをカテゴリカル変数 (categorical variable, 略して変数) と呼ぶ。
カテゴリカル変数は、データの尺度によって2種類に分けることができる。1つはカテゴリカル変数が順序関係を持たない性別や血液型などのような吊義 (nominal) 尺度データであり、もう1つはアンケート調査質問票などでよく用いられている満足度に関する選択肢 (非常に満足する、満足する、なんともいえない、満足しない、非常に満足しない) のようなカテゴリーの順序に意味を持っている順序 (ordinal) 尺度データである。
順序尺度は、カテゴリーの並べの順序に意味を持っているので、順序尺度データのために開発されたデータ解析方法は、吊義尺度には用いられないことに注意することが必要である。
量的データは、ある一定の区間を区切りカテゴリカルデータに変換して解析を行ったり、カテゴリカルデータのカテゴリーに順序のスコアを付与(あるいはダミー変数0,1を導入)して量的のデータとして解析したりする場合もある。
本稿では、Rによるカテゴリカルデータの解析のための操作方法や分割表の統計量や独立性の検定方法について説明を行う。カテゴリカルデータに関連する関数は幾つかのパッケージに用意されているが、本稿ではパッケージ stats、vcd の中の関数を用いる。
2.カテゴリカルデータの操作
まずカテゴリカルデータの記録形式と基本操作について説明する。パッケージ vcd の中には86個体5変数の薬物実験に関するカテゴリカルデータフレーム Arthritis がある。その変数を次に示す。
ID: |
患者のID |
Treatment: |
実験グループの因子 (Placebo, Treated) |
Sex: |
性別 (Female, Male) |
Age: |
患者の年齢 |
mproved: |
観察結果の因子 (None, Some, Marked) |
次のコマンドでデータを確認して見よう。関数 head を用いてデータ先頭の一部分を返すことができる。デフォルトでは先頭6行を返すように指定されている。
> library(vcd);data(Arthritis)
> head(Arthritis,n=4)
|
ID |
Treatment |
Sex |
Age |
Improved |
1 |
57 |
Treated |
Male |
27 |
Some |
2 |
46 |
Treated |
Male |
29 |
None |
3 |
77 |
Treated |
Male |
30 |
None |
4 |
17 |
Treated |
Male |
32 |
Marked |
このようなデータの記録は1行が1つの個体 (ここでは1人の被実験者) であり、1列が1つの変数である。通常アンケートの調査票などの集計にも基本的にはこのような形式を用いる。
このようなデータ形式はしばしば分割表 (contingency table) に集計して分析を行う。上記の形式から分割表に集計する関数としてパッケージ stats には関数 xtabs、パッケージ vcd には関数 structable がある。後者はフラットな分割表を返す。それぞれの書き式を次に示す。
xtabs(formula = ~., data …)
structable(formula, data,…)
引数 formula には、応答変数と説明変数や分割表の表側、表頭の変数を指定する。データ Arthritis の変数 Treatment を表側、変数 Improved を表頭とした分割表を作成する例を示す。
> (ArtTI <- xtabs(~ Treatment + Improved, data = Arthritis))
Improved
Treatment |
None |
Some |
Marked |
Placebo |
29 |
7 |
7 |
Treated |
13 |
7 |
21 |
このよう分割表を2元分割表 (two-way table) と呼ぶ。
引数 subset に条件を付けることにより、その条件に基づいた分割表を作成することができる。例えば、女性に限定した上記のような分割表を返すのには引数 subset = Sex = = "Female" を追加する。
データ Arthritis の中の変数 Age は、量的な変数として扱うことも可能であり、カテゴリカルデータとして扱うことも可能である。
仮に量的のデータとして扱うときには次のコマンドのように引数 formula を指定すると、分割表の各セルには年齢の合計が返される。
> (ArtTI.age<-xtabs(Age~ Treatment +Improved, data = Arthritis))
Improved
Treatment |
None |
Some |
Marked |
Placebo |
1439 |
402 |
403 |
Treated |
648 |
397 |
1193 |
分割表の各セルの度数は異なるので合計より、度数当たりの平均度数を求めた方がよい場合もある。そのコマンドを次に示す。
> ArtTI.age/ArtTI
Improved
Treatment |
None |
Some |
Marked |
Placebo |
49.62069 |
57.42857 |
57.57143 |
Treated |
49.84615 |
56.71429 |
56.80952 |
次のコマンドのように書くことで2つの変数 Sex と Treatment を表側 (あるいは表頭) に配置することができる。
> xtabs(~ I(Sex:Treatment) + Improved, data = Arthritis)
Improved
I(Sex:Treatment) |
None |
Some |
Marked |
Female:Placebo |
19 |
7 |
6 |
Female:Treated |
6 |
5 |
16 |
Male:Placebo |
10 |
0 |
1 |
Male:Treated |
7 |
2 |
5 |
引数 formula の説明変数を3つに指定すると次のような3元分割表が返される。3元分割表には3つの添え字で分割表を操作する。添え字の順序は引数 formula で指定した変数の順序と一致する。
> (ArtTIS<- xtabs(~ Treatment +Improved+Sex, data = Arthritis))
, , Sex = Female
Improved
Treatment |
None |
Some |
Marked |
Placebo |
19 |
7 |
6 |
Treated |
6 |
5 |
16 |
, , Sex = Male
Improved
Treatment |
None |
Some |
Marked |
Placebo |
10 |
0 |
1 |
Treated |
7 |
2 |
5 |
> dim(ArtTIS)
[1] 2 3 2
> ArtTIS[,,1] #女の分割表
<結果は省略>
> ArtTIS[,,2] #男の分割表
<結果は省略>
多元分割表をフラット分割表に変換する関数としてパッケージ stats にはftable、パッケージ vcd には structable がある。その例を次に示す。
> ftable(ArtTIS)
Sex Female Male
Treatment Improved
Placebo |
None |
19 |
10 |
|
Some |
7 |
0 |
|
Marked |
6 |
1 |
Treated |
None |
6 |
7 |
|
Some |
5 |
2 |
|
Marked |
16 |
5 |
表側、表頭の組み合わせなどは引数 row.vars , col.vars を次のように指定することができる。
> ftable(ArtTIS,row.vars ="Improved")
Treatment Placebo Treated |
|
Sex |
Female |
Male |
Female |
Male |
Improved |
|
|
|
None |
19 |
10 |
6 |
7 |
Some |
7 |
0 |
5 |
2 |
Marked |
6 |
1 |
16 |
5 |
関数 structable の使用例を次に示す。
> structable(ArtTIS)
|
Improved |
None |
Some |
Marked |
Treatment |
Sex |
|
|
|
Placebo |
Female |
19 |
7 |
6 |
|
Male |
10 |
0 |
1 |
Treated |
Female |
6 |
5 |
16 |
|
Male |
7 |
2 |
5 |
場合によっては、多元分割表の変数を併合して分析すると都合がよい場合がある。Rには配列を併合する関数 margin.table がある。3元分割表 ArtTIS を2元に併合する例を次に示す。
> margin.table(ArtTIS,1:2)
Improved
Treatment |
None |
Some |
Marked |
Placebo |
29 |
7 |
7 |
Treated |
13 |
7 |
21 |
Rでは関数 prop.table を用いて分割表の各セルの比率を求めることができる。分割表の度数の総合計を基準とする各セルの比率は次のように求める。
> ArtTI.p<-prop.table(ArtTI)
> round(ArtTI.p,3)
Improved
Treatment |
None |
Some |
Marked |
Placebo |
0.345 |
0.083 |
0.083 |
Treated |
0.155 |
0.083 |
0.250 |
列の合計を基準とするときには引数 margin=1,行の合計を基準とするときには引数 margin=2 を追加する。
> prop.table(ArtTI, margin=1)
Improved
Treatment |
None |
Some |
Marked |
Placebo |
0.6744186 |
0.1627907 |
0.1627907 |
Treated |
0.3170732 |
0.1707317 |
0.5121951 |
説明の便利のため、まず2元分割表の一般形式を表1に示す。
表1 2元分割表の一般形式と周辺度数
|
b1 |
b2 |
… |
bj |
… |
b1 | 計 |
a1 | n11 |
n12 |
… |
n1j |
… |
n1c |
n1+ |
a2 | n21 |
n22 |
… |
n2j |
… |
n2c |
n2+ |
… |
… |
… |
… |
… |
… |
… |
… |
ai |
ni1 |
ni2 |
… |
nij |
… |
nic |
ni+ |
… |
… |
… |
… |
… |
… |
… |
… |
ar |
nr1 |
nr2 |
… |
nrj |
… |
nrc |
nr+ |
計 |
n1+ |
n2+
| … |
ni+ |
… |
nr+ |
n++ |
2元分割表の分析に用いる統計量の中で最も広く知られているのはピアソンの χ2 統計量 (Pearson's chi-squared statistics) である。略してカイ2乗統計量と呼ぶことにする。カイ2乗統計量は次の式で定義されている。
この式で得られた統計量は自由度 ( r - 1 ) ( c - 1 ) の χ2 分布に従うことが知られている。式の中の nij は分割表のi 行j 列セルの度数であり、Eij はi 行j 列セルの期待度数 (expected frequencies) である。
期待度数 Eij は、i 行の和 ni+、j 列和 n+j、分割表の度数の総合計 n++ で次のように求める。分割表の ni+ , n+j , n++ を周辺度数とも呼ぶ。
パッケージ vcd には、2元分割表の周辺度数を求める関数 mar_table、期待度数を求める関数 independence_table がある。
第2節で作成した2元分割表 ArtTI を用いた例を次に示す。
> mar_table(ArtTI)
Improved
Treatment |
None |
Some |
Marked |
TOTAL |
Placebo |
29 |
7 |
7 |
43 |
Treated |
13 |
7 |
21 |
41 |
TOTAL |
42 |
14 |
28 |
84 |
> independence_table(ArtTI)
Improved
Treatment |
None |
Some |
Marked |
Placebo |
21.5 |
7.166667 |
14.33333 |
Treated |
20.5 |
6.833333 |
13.66667 |
分割表の各セルの値
をピアソン残差 (Person residuals) と呼ぶ。ピアソン残差は近似的に標準正規分布に従うことが知られている。
パッケージ stats には、適合検定や分割表の独立性の検定に必要となるカイ2乗統計量を求める関数 chisq.test がり、カイ2乗統計量だけではなく、期待度数、ピアソン残差も返す。2元分割表 ArtTI を用いた例を次に示す。
> (ArtTI.chi<-chisq.test(ArtTI))
Pearson's Chi-squared test
data: ArtTI
X-squared = 13.055, df = 2, p-value = 0.001463
> summary(ArtTI.chi)
|
Length |
Class |
Mode |
|
1 |
-none- |
numeric |
|
1 |
-none- |
numeric |
p.value |
1 |
-none- |
numeric |
|
1 |
-none- |
character |
data.name |
1 |
-none- |
character |
|
6 |
xtabs |
numeric |
|
6 |
-none- |
numeric |
|
6 |
xtabs |
numeric |
> ArtTI.chi$expected
Improved
Treatment |
None |
Some |
Marked |
Placebo |
21.5 |
7.166667 |
14.33333 |
Treated |
20.5 |
6.833333 |
13.66667 |
> round(ArtTI.chi$residuals,2)
Improved
Treatment |
None |
Some |
Marked |
Placebo |
1.62 |
-0.06 |
-1.94 |
Treated |
-1.66 |
0.06 |
1.98 |
関数 chisq.test のデフォルトは、2 × 2 分割表のときには,連続性の補正を行うように引数が "correct=TRUE" に指定されている。連続修正は次の式を用いている。
分割表のカイ2自乗統計量の近似値として尤度比統計量 (likelihood ratio statistics)
がある。分割表の独立検定にはカイ2乗統計量が多く用いられているが、モデルの比較や選択には尤度比統計量が多く用いられている。
パッケージ vcd の中の関数 assocstats は分割表のカイ2乗統計量、尤度比統計量、ファイ連関係数 φ、ピアソンの連関係数 C、グラメール (Gramer) の連関係数 V を返す。これらの連関係数は次のように定義されている。式の中の k は行と列の数の中の少ない方の数を取る。
分割表 ArtTI を用いた例を次に示す。
> assocstats (ArtTI)
|
X^2 |
df |
P(> X^2) |
Likelihood Ratio |
13.530 |
2 |
0.0011536 |
Pearson |
13.055 |
2 |
0.0014626 |
Phi-Coefficient |
: |
0.394 |
Contingency Coeff. |
: |
0.367 |
Cramer's V |
: |
0.394 |
表2のような2×2の分割表は、しばしば特殊なケースとして扱われている。
表2 2×2の分割表
|
b1 | b2 |
計 |
a1 |
n11 |
n12 |
n1+ |
a2 |
n21 |
n22 |
n2+ |
計 |
n+1 |
n+2 |
n++ |
2×2の分割表のピアソンカイ2乗値は

で計算されるが、イェーツ (Yates) はよりカイ2乗分布に近似するように次のように連続性の補正を施して用いることを提唱した。

セルの度数が非常に低く、期待値が5以下のセルが全セルの25%以上を含むときには上記の式によるカイ2乗独立性の検定は上適切である。そこでフィッシャー (Fisher) は超幾何分布を用いて、2×2の分割表の正確検定統計量の計算式を導出した。

この確率を用いた検定をフィッシャーの正確検定 (Fisher's exact test) と呼ぶ。Rにはフィッシャーの正確検定値を求める関数 fisher.test がある。
ここではフィッシャーが1935年に公表した有吊な紅茶の実験データを用いることにする。フィッシャーの紅茶の実験は「紅茶にミルクを注ぐ《と「ミルクに紅茶を注ぐ《という異なる2つの方法で入れた紅茶を、その味から紅茶の入れ方を推定するテストである。2つの入れ方によるそれぞれ4杯の紅茶を、味からその入れ方を識別することが可能であると吊乗る女性に識別テストを行ったデータを表3に示す。
表3 フィッシャーの紅茶実験データ
(表の中のミルク、紅茶は先に入れたもの)
実際 |
識別結果 |
ミルク |
紅茶 |
合計 |
ミルク |
3 |
1 |
4 |
紅茶 |
1 |
3 |
4 |
合計 |
4 |
4 |
8 |
本例では、フィッシャーの両側検定は
p(0) +
p(1) +
p(3) +
p (4) 、片側検定は
p(0) +
p (1) (左側検定)、あるいは
p(3) +
p (4) (右側検定)値を用いる。
表3のデータを用いたフィッシャーの正確検定値を求めるコマンドを次に示す。
> (TeaTasting<-matrix(c(3,1,1,3),nr = 2,dimnames = list(実際 = c("ミルク","紅茶"),識別 = c("ミルク", "紅茶"))))
識別
実際 |
ミルク |
紅茶 |
ミルク |
3 |
1 |
紅茶 |
1 |
3 |
> fisher.test(TeaTasting, alternative = "greater") #片側検定統計量を算出
Fisher's Exact Test for Count Data
data: TeaTasting
p-value = 0.2429
<後略>
この結果では、独立性の仮説が棄却されないため、残念ながらこの小規模のテストでは、紅茶の味からその入れ方を識別することが可能であるという女性の主張は統計的に認めることはできない。
表側 a と表頭 b が独立であるならば、
が1になるはずである。この比をオッズ比 (Odds Ratios) と呼ぶ。フィッシャーの正確検定はオッズ比
である仮説検定である。関数 fisher.test では、両側検定がデフォルトに指定されている。
2×2の分割表の場合はオッズ比の検定が多く用いられている。パッケージ vcd にはオッズ比に関連する関数 oddsratio, woolf_test などがある。誌面上の都合によりその使用例は省略する。
カテゴリカルデータによく用いられている幾つか検定に関する関数を表4に示す。パッケージ Hmisc の中の関数 summary.formula は、分割表の作成や比率と検定統計量などを求めることができる。
表4 検定に関する幾つかの関数
関数 |
機能 |
パッケージ |
assocstats |
分割表の統計量と検定:カイ2乗、尤度比、連関係数 |
vcd |
binom.test |
二項分布の検定 |
stats |
coindep_test |
条件付き独立検定 |
vcd |
fisher.test |
Fisherの直接検定 |
stats |
goodfit |
適合検定 |
vcd |
kruskal.test |
ランク和の検定 |
stats |
mantelhaen.test |
Cochran-Mantel -Haenszelの3元分割表のカイ2乗検定 |
stats |
mcnemArtest |
正方形の分割表のMcnemarのカイ2乗検定 |
stats |
oddsratio |
オッズ比と検定統計量 |
vcd |
prop.test |
比率の検定 |
stats |
summary.formula |
分割表の統計量など |
Hmisc |
woolf_test |
同種の2×2×k分割表のWoolfのオッズ比検定 |
vcd |