Rで独立性のカイ二乗検定 そのまま使える自作関数

[公開日]2016/12/21[更新日]2017/10/25 [カテゴリー]R言語入門 Written by  IMIN

今回は、R言語を用いて、カイ二乗検定を行う方法をご紹介します。カイ二乗検定はカイ二乗分布を用いた、検定手法です。カイ二乗分布について詳しくは、カイ二乗分布のわかりやすいまとめにて、まとめました。

カイ二乗検定について、理論を知りたい方はこちら
独立性のカイ二乗検定 例題を用いてわかりやすく解説

独立性の検定(カイ二乗検定)を行なって、その結果をcsvファイルに出力する自作関数も作成しました。当記事の下の方につけたので、そちらも合わせてお役立てください。

データ数が少なくて、カイ二乗分布への近似が不正確になるときは、フィッシャーの正確確率検定を使います。それはこちらに記載しています→Rでフィッシャーの正確確率検定 そのまま使える自作関数例

※当サイトのプログラムをコピーアンドペーストすると、エラーが出ることがあるようです(特にmacの場合)。ペースト時にコマンド+alt+Vとすると、エラーなく実行出来ることがあります。(コマンド+Vではなく。)

 
 

カイ二乗検定を行う関数"chisq.test()"

R言語では、カイ二乗検定を行うのに、あらかじめ"chisq.test()"という関数が用意されています。以下のように、分割表にしたデータフレームや行列を、引数に入れて使います。

プログラム

a = matrix(c(10,47,18,15),2,2) #aという行列を作成
chisq.test(a) #カイ二乗検定を実施

Pearson's Chi-squared test with Yates' continuity correction
data: a
X-squared = 11.681, df = 1, p-value = 0.0006315

上記のプログラムではaという行列を引数にとって、カイ二乗検定を行なっています。この表示されている結果の見方は、

X-squared:カイ二乗統計量
df:自由度
p-value:p値

となります。p値があらかじめ設定していた、有意水準よりも小さければ、帰無仮説を棄却し、対立仮説である「二つの変数は独立ではない」という仮説を採択します。

Rでのカイ二乗検定の結果の見方

chisq.testによるカイ二乗検定の結果の見方を紹介します。まずは、操作を簡単にするために、検定の結果を何かしらの変数に代入しておきます。ここではresultに代入します。

プログラム

result = chisq.test(a)

期待度数の見方

期待度数を見るには、カイ二乗検定の結果のうしろに$expectedをつけるか、[7]とします。

プログラム

result$expected
result[7]

[,1] [,2]
[1,] 17.73333 10.26667
[2,] 39.26667 22.73333

ピアソン残差の見方

ピアソン残差を見るには、カイ二乗検定の結果のうしろに$residualsをつけるか、[8]とします。

プログラム

result$residuals
result[8]

[,1] [,2]
[1,] -1.836418 2.413526
[2,] 1.234112 -1.621941

調整残差の見方

調整残差を見るには、カイ二乗検定の結果のうしろに$stdresをつけるか、[9]とします。

プログラム

result$stdres
result[9]

[,1] [,2]
[1,] -3.653939 3.653939
[2,] 3.653939 -3.653939

独立性のカイ二乗検定の自作関数

カイ二乗検定を行って、csvファイルへの結果の出力までを行う自作関数を作成しました。

chi_test(解析用データ(2列),出力ファイル名)

として使用します。

プログラム

chi_test <- function(dat,filename
{
ans = chisq.test(table(dat))
x2 = ans$statistic
dof = ans$parameter
p_value = ans$p.value
ex = ans$expected
res = ans$residuals
std = ans$stdres
name= paste(colnames(dat)[1],"vs",colnames(dat)[2])
result = cbind(table(dat),x2,dof,p_value,ex,res,std)
result[2:nrow(result),c((ncol(table(dat))+1):(ncol(table(dat))+3))] = ""
write.table(matrix(c(name,colnames(result)),nrow=1),filename,append=T,quote=F,sep=",",row.names=F,col.names=F,fileEncoding="CP932")
write.table(result,filename,append=T,quote=F,sep=",",row.names=T,col.names=F,fileEncoding="CP932")
write.table("",filename,append=T,quote=F,sep=",",row.names=F,col.names=F,fileEncoding="CP932")
}

出力結果は、

解析に用いた分割表、カイ二乗統計量、自由度検定のp値、期待度数、ピアソン残差、調整残差

の順番になっています。

独立性の検定の自作関数(p値,カイ二乗統計量、自由度を前に出力)

出力の順番を変えたプログラムも掲載して欲しいとの要求があったので、以下に掲載します。やってることは上記の自作関数と同じですが、出力の順番を変えました。

カイ二乗統計量、自由度、p値、解析に用いた分割表、期待度数、ピアソン残差、調整残差

の順に出力します。こちらの関数のメリットとしては、分割表の大きさが異なる解析結果を多数同じcsvファイルに出力したときにp値の列を揃う点です。

プログラム

chi_test2 <- function(dat,filename) 

{

ans = chisq.test(table(dat))

x2 = ans$statistic

dof = ans$parameter

p_value = ans$p.value

ex = ans$expected

res = ans$residuals

std = ans$stdres

name = paste(colnames(dat)[1],"(行) vs",colnames(dat)[2],"(列)")

result = cbind(x2,dof,p_value,rownames(table(dat)),table(dat),ex,res,std)

result[2:nrow(result),c(1:3)] = ""

rownames(result)[1:nrow(result)] = ""

write.table(matrix(c(name,colnames(result)),nrow=1),filename,append=T,quote=F,sep=",",row.names=F,col.names=F,fileEncoding="CP932")

write.table(result,filename,append=T,quote=F,sep=",",row.names=T,col.names=F,fileEncoding="CP932")

write.table("",filename,append=T,quote=F,sep=",",row.names=F,col.names=F,fileEncoding="CP932")

}

使用例と注意事項

プログラム

dat = df[,c(1,2)]
chi_test(dat,"独立性.csv")

以上のプログラムで元データの1列目と2列目の変数に関して、独立性の検定を行い、結果の出力をします。

独立性の検定を何種類も行う場合は、for文を使うのが便利です。

プログラム

for(i in 3:5)
{
for(j in 6:12)
{
dat = df[,c(i,j)]
chi_test(dat,"独立性.csv")
}
}

注意事項

上で示した自作関数は、分割表を作る前のデータを引数に取ります。データの中の各カテゴリーの度数を計算して、分割表を関数内で作っています。分割表をすでに作成済みの方は、関数内の

table(dat)

と書かれた部分を全て

dat

と書き換えることでプログラムを使えるはずです。

  • スポンサーリンク

  • コメントを残す

    メールアドレスが公開されることはありません。 * が付いている欄は必須項目です

    CAPTCHA