二変数の分割表(クロス表)
調査結果の報告においては,高度な分析をするだけではなく,収集した情報の基本的な集計結果を示す事も重要になる。
しかしRでは,度数分布表(単純集計表などとも呼ばれる)や分割表(連関表,クロス表)について,必要最低限の結果を出力するコマンドが基本で,色々な情報を纏めて提示する為には自分で少し手を加える必要がある。
以下では,学部学生の演習レヴェルであると便利だと思われる作表方法を紹介する。
模擬データの作成
以下は,二つの変数(性別を想定した3カテゴリと学歴を想定した5カテゴリ,NA混じり)を確率的に発生させて,それをデータ・フレイムにしている。
性別と教育年数の架空データ・フレイム作成コード
# 性別によって学歴分布を変化させる為に若干技巧的。
n <- 800
sex <- sample(c(1, 2, 3, NA),
size = n,
replace = T,
prob = c(.43, .48, .06, .03))
school <- ((!is.na(sex) & sex == 1) | is.na(sex)) *
sample(c(1:5, NA),
size = n,
replace = T,
prob = c(.13, .35, .10, .32, .05, .05)) +
(!is.na(sex) & sex >= 2) *
sample(c(1:5, NA),
size = n,
replace = T,
prob = c(.08, .25, .20, .42, .03, .02))
d01 <- data.frame(sex, school)
それぞれの変数にアクセスする為には,d01$sex,d01$schoolとする事になる。
二つのカテゴリカル変数の分割表(クロス表)
社会学でクロス表と呼ぶ事の多い分割表(連関表; contingency table)は二つのカテゴリカル変数の集計表であるが,これを作成するRの関数も table( ) である。table( ) に二つの変数を引数で与えると二変数の分割表になる。一変数の場合と同じく,useNA = オプションを指定しないとNAを除外した表になる。
table( ) で作成した分割表は表の本体だけであるが,集計結果を示す場合には,行周辺度数・列周辺度数なども併せて表示する。
table( ) で作成された分割表が t01 と云うオブジェクトに格納されているとすると,周辺度数を付加した表は,addmargins( ) 関数にt01を引数として与えれば作成される。
大抵の場合,度数ではなく相対度数の表も作成したい。関数は一変数の時と同じ prop.table( ) 関数である。
何もオプションを付けないと全体%の表になり,margin = 1 と云うオプションを付けると行%の表になる。margin = 2 では列%の表である。
オブジェクト名は適当であり,自由に付けてよい。
table(d01$sex, d01$school)
1 2 3 4 5
1 46 131 41 87 15
2 24 82 76 151 15
3 5 9 7 25 0
table(d01$sex, d01$school, useNA = "ifany")
1 2 3 4 5 <NA>
1 46 131 41 87 15 27
2 24 82 76 151 15 25
3 5 9 7 25 0 3
<NA> 7 10 3 7 2 2
prop.table(table(d01$sex, d01$school), margin = 1)
1 2 3 4 5
1 0.14375000 0.40937500 0.12812500 0.27187500 0.04687500
2 0.06896552 0.23563218 0.21839080 0.43390805 0.04310345
3 0.10869565 0.19565217 0.15217391 0.54347826 0.00000000
prop.table(table(d01$sex, d01$school), margin = 2)
1 2 3 4 5
1 0.61333333 0.59009009 0.33064516 0.33079848 0.50000000
2 0.32000000 0.36936937 0.61290323 0.57414449 0.50000000
3 0.06666667 0.04054054 0.05645161 0.09505703 0.00000000
prop.table(table(d01$sex, d01$school))
1 2 3 4 5
1 0.064425770 0.183473389 0.057422969 0.121848739 0.021008403
2 0.033613445 0.114845938 0.106442577 0.211484594 0.021008403
3 0.007002801 0.012605042 0.009803922 0.035014006 0.000000000
addmargins(table(d01$sex, d01$school))
1 2 3 4 5 Sum
1 46 131 41 87 15 320
2 24 82 76 151 15 348
3 5 9 7 25 0 46
Sum 75 222 124 263 30 714
with( ) と round( ) を使い,適宜オブジェクト保存して見易くしよう。
t01T <- with(d01, table(sex, school, useNA = "ifany"))
t01T
school
sex 1 2 3 4 5 <NA>
1 46 131 41 87 15 27
2 24 82 76 151 15 25
3 5 9 7 25 0 3
<NA> 7 10 3 7 2 2
t01 <- with(d01, table(sex, school))
t01
school
sex 1 2 3 4 5
1 46 131 41 87 15
2 24 82 76 151 15
3 5 9 7 25 0
pr01 <- prop.table(t01, margin = 1)
round(pr01*100, 1)
school
sex 1 2 3 4 5
1 14.4 40.9 12.8 27.2 4.7
2 6.9 23.6 21.8 43.4 4.3
3 10.9 19.6 15.2 54.3 0.0
pc01 <- prop.table(t01, margin = 2)
round(pc01*100, 1)
school
sex 1 2 3 4 5
1 61.3 59.0 33.1 33.1 50.0
2 32.0 36.9 61.3 57.4 50.0
3 6.7 4.1 5.6 9.5 0.0
ps01 <- prop.table(t01)
round(ps01*100, 1)
school
sex 1 2 3 4 5
1 6.4 18.3 5.7 12.2 2.1
2 3.4 11.5 10.6 21.1 2.1
3 0.7 1.3 1.0 3.5 0.0
t01m <- addmargins(t01)
t01m
school
sex 1 2 3 4 5 Sum
1 46 131 41 87 15 320
2 24 82 76 151 15 348
3 5 9 7 25 0 46
Sum 75 222 124 263 30 714
factor(因子,要因)型変数を用いてラベルを付ける
SPSSでは変数に変数ラベル,値に値ラベルを定義する事が出来,出力を見易くする事が出来るが,Rはこの点については不親切である。
もっとも特に値ラベルについては,実体としての数値が何であるかが分からなくなって却って初心者が間違える場合もあるので,うまく付けないと(或いは適切に表示オプションを設定しないと)便利ではなくなる。
ここでは,数値型変数から実質的にそれと同じfactor型(因子型,要因型)変数を作り,それで分割表や(下記の)モザイクプロットを表示する方法を示す。
数値型変数は小文字の変数名を持っているので,大文字の変数名で対応する因子型変数を作成する方針とする。
d01$SEX <- factor(d01$sex,
levels = 1:3,
labels = c("男性", "女性", "その他"))
新変数を作成したら必ず新旧変数の分割表で齟齬が無いかどうかを確認する。
with(d01, table(sex, SEX, useNA = "ifany"))
SEX
sex 男性 女性 その他 <NA>
1 347 0 0 0
2 0 373 0 0
3 0 0 49 0
<NA> 0 0 0 31
d01$SCHOOL <- factor(d01$school,
levels = 1:5,
labels = c("中学", "高校", "短大", "四大", "院"))
with(d01, table(school, SCHOOL, useNA = "ifany"))
SCHOOL
school 中学 高校 短大 四大 院 <NA>
1 82 0 0 0 0 0
2 0 232 0 0 0 0
3 0 0 127 0 0 0
4 0 0 0 270 0 0
5 0 0 0 0 32 0
<NA> 0 0 0 0 0 57
因子型変数で分割表を作成する。
t_SEX_SCHOOL <- with(d01, table(SEX, SCHOOL, useNA = "ifany"))
t_SEX_SCHOOL
SCHOOL
SEX 中学 高校 短大 四大 院 <NA>
男性 46 131 41 87 15 27
女性 24 82 76 151 15 25
その他 5 9 7 25 0 3
<NA> 7 10 3 7 2 2
分割表の行と列にラベルを付けなくても,因子型変数のラベルが用いられるので数値だけで表示されるより分かり易い。
同じ質問項目に数値型変数と因子型変数の二つを用意しておくと,場合に応じて使い分けられて便利である。
ここでは分割表にラベルを付ける単純な方法を紹介する。
with(d01, addmargins(table(sex)))
sex
1 2 3 Sum
347 373 49 769
with(d01, addmargins(table(sex, school)))
school
sex 1 2 3 4 5 Sum
1 46 131 41 87 15 320
2 24 82 76 151 15 348
3 5 9 7 25 0 46
Sum 75 222 124 263 30 714
SPSSの様な変数ラベルや値ラベルがついていないのでこれでは分りにくい。
性別の度数分布表をオブジェクトに格納して,変数値や変数にラベルを付けてみよう。
t10 <- with(d01, table(sex))
rownames(t10) <- c("1 男性", "2 女性", "3 その他")
t10
sex
1 男性 2 女性 3 その他
347 373 49
names(dimnames(t10)) <- "性別 sex"
t10
性別 sex
1 男性 2 女性 3 その他
347 373 49
上と同様にして,分割表の行と列にラベルを付けよう。
t11 <- with(d01, table(sex, school))
rownames(t11) <- c("1 男性", "2 女性", "3 その他")
colnames(t11) <- c("1 中学", "2 高校", "3 短大", "4 四大", "5 院")
names(dimnames(t11)) <- c("性別 sex", "最終学歴 school")
t11
最終学歴 school
性別 sex 1 中学 2 高校 3 短大 4 四大 5 院
1 男性 46 131 41 87 15
2 女性 24 82 76 151 15
3 その他 5 9 7 25 0
変数そのものにラベルを貼り付けるのとは違って集計表の行や列に名前を付けているだけなので,正直このやり方ではいちいち非常に面倒であり,SPSS的な変数ラベル,値ラベルを使いたいと思う事も多いだろう。
しかしそこは逆に,常に値とラベルの対応に注意させられる事で,ラベルに惑わされて実態としての値(数字)を取り違え,処理や解釈を誤ると云うリスクがなくなると前向きに考えておこう(実際初心者がSPSSを使う場合にこのリスクは小さくない)。
Rで変数ラベル,値ラベルを使用可能にする工夫も幾つか開発されているようが,ここでは極力追加アプリケイションや追加パッケイジをインストールしなくて済む方法を紹介しているので,割愛する。むしろ,factor型変数を作って分割表やグラフにする方法を推奨する。
度数と比率を共に含んだ分割表を出力する
実習などでは,レポートに貼り付ける分割表(クロス表)を作成するのが手間である,
Rのコンソール画面の出力は余りにそっけなくてそのままでは使えず,Excelなどで色々と手を加える事が多いだろう。
ここでは,なるべく加工する必要の少ない表をcsvオブジェクトでファイル出力する方法を工夫しよう。
ここではNAを含まない分割表で例示するが,若干の修正をすればNAを含む分割表にも対応するだろう。
度数分布表と比率の表
ftab01 <- with(d01, table(SEX, SCHOOL, useNA = "no")) # NAを含まない分割表
mtab01 <- addmargins(ftab01) # 行と列の周辺度数付き
prtab01 <- prop.table(ftab01, margin = 1) # 行比率の分割表
# pctab01 <- prop.table(ftab01, margin = 2) # 列比率の分割表
# 行と列の周辺度数の割合を計算する
mctab01 <- prop.table(mtab01[1:nrow(ftab01), ], margin = 2)
mrtab01 <- prop.table(mtab01[ , 1:ncol(ftab01)], margin = 1)
.ptab01 <- cbind(prtab01, mctab01[ , ncol(mctab01)])
ptab01 <- rbind(.ptab01, c(mrtab01[nrow(mrtab01), ], 1))
mtab01; ptab01
SCHOOL
SEX 中学 高校 短大 四大 院 Sum
男性 46 131 41 87 15 320
女性 24 82 76 151 15 348
その他 5 9 7 25 0 46
Sum 75 222 124 263 30 714
中学 高校 短大 四大 院
男性 0.14375000 0.4093750 0.1281250 0.2718750 0.04687500 0.44817927
女性 0.06896552 0.2356322 0.2183908 0.4339080 0.04310345 0.48739496
その他 0.10869565 0.1956522 0.1521739 0.5434783 0.00000000 0.06442577
0.10504202 0.3109244 0.1736695 0.3683473 0.04201681 1.00000000
度数分布表と比率表を組み合わせる工夫
gtab01 <- .gtab01 <- rbind(mtab01, round(ptab01*100, 1))
gtab01[seq(1, nrow(.gtab01), 2), ] <- .gtab01[1:nrow(mtab01), ]
gtab01[seq(2, nrow(.gtab01), 2), ] <- .gtab01[(nrow(mtab01) + 1):nrow(.gtab01), ]
row.names(gtab01)[seq(1, nrow(.gtab01), 2)] <- row.names(mtab01)
row.names(gtab01)[seq(2, nrow(.gtab01), 2)] <- ""
gtab01
中学 高校 短大 四大 院 Sum
男性 46.0 131.0 41.0 87.0 15.0 320.0
14.4 40.9 12.8 27.2 4.7 44.8
女性 24.0 82.0 76.0 151.0 15.0 348.0
6.9 23.6 21.8 43.4 4.3 48.7
その他 5.0 9.0 7.0 25.0 0.0 46.0
10.9 19.6 15.2 54.3 0.0 6.4
Sum 75.0 222.0 124.0 263.0 30.0 714.0
10.5 31.1 17.4 36.8 4.2 100.0
- 最終列の行周辺度数の下の%は,行周辺度数が全体に占める%
- 最終行,列周辺度数の下の%は,列周辺度数が全体に占める%
- 度数と%を一つの表に合併した為,度数にも少数点がついている。これはMS Excel上などで調整するしかない。
- %の小数点は,小数第2位を四捨五入して小数第1位までとしている。
- 学生がMS Excelで編集する事を想定して,文字コードは敢えて shift-JIS で保存した。インターネットやMacで標準の UTF-8 で保存したい場合には fileEncoding を変えるか削除すればよい。
- 行名・列名のSumは必要なら「小計」等に書き換えて貰う。
write.csv(gtab01, fileEncoding = "sjis",
file = paste0(names(dimnames(ftab01))[1], "×",
names(dimnames(ftab01))[2], ".csv"))
分割表をそのまま図示する
分割表をグラフ表示するには,伝統的な帯グラフよりはモザイク・プロットの方が良い。
そうでなければ横に並べた棒グラフだろう。
with(d01,
mosaicplot(SEX ~ SCHOOL,
las = 1,
main = "性別ごとの学歴",
col = terrain.colors(length(unique(SCHOOL))-1))
)
with(d01,
barplot(table(SCHOOL, SEX),
beside = T,
col = terrain.colors(length(unique(SCHOOL))-1),
legend = T)
)
with(d01,
mosaicplot(SEX ~ SCHOOL, shade = T)
)
因子型変数を用い,オプションを色々と設定した。
最初は複数のオプションを一度に設定するのではなく,一つ一つ追加してそれがどの様な設定を行うのかを確認して理解してから使うこと。
参考:超初心者向けのRガイドgー標準関数でのモザイクプロット
埋め込む統計量の準備コード
tab01 <- with(d01, table(SEX, SCHOOL))
chi01 <- chisq.test(tab01)
V <- sqrt(chi01$statistic/(min(dim(tab01) - 1) * sum(tab01)))
names(V) <- "Cramer's V"
with(d01,
mosaicplot(SEX ~ SCHOOL,
col = terrain.colors(5),
las = 1,
off = 3,
border = "#00000060",
cex = 1.0,
main = "性別と学歴のモザイクプロット",
sub = sprintf("Cramer's V = %.3f", V),
xlab = "性別",
ylab = "学歴")
)