はじめに
Rでつかえる統計検定(t.test, var.test, cor.test, chisq.test)の使用方法のまとめを行った。自身で作ったデータセットで検証しながら各コマンドの出力について解説する。
この記事で用いるデータセット
量的変数のデータ
簡単に計算ができるデータセットが欲しかったので、こちらの記事を参考にRで実装してみる。
### 身長とBMIから体重を作成してデータフレーム化する
m_height <- rnorm(n=1000,mean=171.7,sd=6.6)
f_height <- rnorm(n=1000,mean=158.3,sd=5.7)
m_bmi <- rnorm(n=1000,mean=23.12,sd=4.24)
f_bmi <- rnorm(n=1000,mean=20.82,sd=3.42)
m_weight <- (m_height * 0.01)^2 * m_bmi
f_weight <- (f_height * 0.01)^2 * f_bmi
male <- data.frame(HEIGHT=m_height, WEIGHT=m_weight, SEX="male") #男性のみのデータ
female <- data.frame(HEIGHT=f_height, WEIGHT=f_weight, SEX="female") #女性のみのデータ
dat <- rbind(male, female) #結合したデータ
### できたデータをプロットしてみる
plot(0, 0, type = "n", xlim = c(min(m_weight, f_weight), max(m_weight, f_weight)), ylim = c(min(m_height, f_height), max(m_height, f_height)), xlab="weight [kg]", ylab="height [cm]" )
points(m_weight, m_height, col = "blue", pch = 1)
points(f_weight, f_height, col = "red", pch = 3)
legend("bottomright", legend = c("male", "female"), pch = c(1, 3), col = c("blue", "red"))
身長、体重で相関があるデータセットを作ることができた。
質的変数のデータ
質的変数も検定できるよう性別と血液型からなるデータセットを生成したい。参考として赤ちゃんの血液型の発現割合をこちらのサイトから調べて、人口比率をこちらのサイトから調べた。それぞれの割合になるよう乱数を判定してデータセットを生成した。
myfunc.mf <- function(x) {
a <- c()
for(i in 1:x){
tmp <- runif(1)
if (tmp < 0.513986) a <- c(a, "male")
else a <- c(a, "female")
}
return(a)
}
myfunc.bt <- function(x) {
a <- c()
for(i in 1:x){
tmp <- runif(1)
if (tmp < 0.39) a <- c(a, "A")
else if (0.39 <= tmp & tmp < 0.61) a <- c(a, "B")
else if (0.61 <= tmp & tmp < 0.90) a <- c(a, "O")
else a <- c(a, "AB")
}
return(a)
}
mf <- myfunc.mf(2000)
bt <- myfunc.bt(2000)
dat2 <- data.frame(SEX=mf, BLOOD=bt)
table(dat2)
BLOOD
SEX A AB B O
female 381 117 197 284
male 372 111 258 280
二つの変数の間に連関が無い質的変数のデータセットが完成した。
統計検定
t.testを用いた検定。
平均値の検定（One Sample t-test）
正規母集団からの無作為標本である場合、標本平均 $\bar{X}$ 、標本の不偏分散の正の平方根 $\hat{\sigma}$ 、標本の個数 $n$ を用いて
t = \frac{\bar{X}-\mu}{\hat{\sigma}/\sqrt{n}}
で計算される統計量 $t$ は帰無仮説のもとで自由度 $df=n-1$ の $t$ 分布に従う。この性質を利用して統計検定を行うRコマンドがt.testである。
男性の身長データについて$\mu = 170.0$ の検定を行ってみる。
t.test(male["HEIGHT"], mu=170.0)
One Sample t-test
data: male["HEIGHT"]
t = 8.2108, df = 999, p-value = 6.743e-16
alternative hypothesis: true mean is not equal to 170
95 percent confidence interval:
171.2686 172.0655
sample estimates:
mean of x
171.6671
データ、統計量・自由度・p値、対立仮説、95%の信頼区間、標本から予測される平均値が出力される。もちろん母集団の平均は男性 171.7 cmでデータを作成しているので $p=0.05$ の水準で帰無仮説（平均値は170である）は棄却される。
2つの平均値の比較(Two Sample t-test)
二つの標本の母集団の平均が一致するかどうかを調べたい時もt.testを使うことができる。検定統計量は、
t = \frac{\bar{X_1}-\bar{X_2}}{\sqrt{\frac{(n_1-1)\hat{\sigma_1}^2+(n_2-1)\hat{\sigma_2}^2}{n_1 + n_2 -2}\left(\frac{1}{n_1}+\frac{1}{n_2} \right)}}
となり、帰無仮説の下で自由度 $df=n_1+n_2-2$ の $t$ 分布に従う。ただしこの検定が使える条件は２つの母集団の分散が等質であることが要求される。分散の等質性の検定に詳細を記述する。
t.testを使って、男性の身長と女性の身長の母集団の平均値が等しいのか検定をしていく。
t.test(male["HEIGHT"], female["HEIGHT"], var.equal = TRUE)
Two Sample t-test
data: male["HEIGHT"] and female["HEIGHT"]
t = 48.889, df = 1998, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
12.57924 13.63064
sample estimates:
mean of x mean of y
171.6671 158.5621
データ、統計量・自由度・p値、対立仮説、95%の信頼区間、標本から予測される平均値が出力される。もちろん母集団の平均は男女で変えて行っているので $p=0.05$ の水準で帰無仮説（平均値が男女で一致する）は棄却される。
Welch の検定
分散の等質性の検定で分散の等質性が得られなかった場合に、二つの平均値を比較するには Welch の検定を行う。これも t.test で実現できる。先ほどのコマンドのオプションをvar.equal = FLASEとすればよいだけである。
t.test(male["HEIGHT"], female["HEIGHT"], var.equal = FALSE)
Welch Two Sample t-test
data: male["HEIGHT"] and female["HEIGHT"]
t = 48.889, df = 1955.5, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
12.57924 13.63065
sample estimates:
mean of x mean of y
171.6671 158.5621
データ、統計量・自由度・p値、対立仮説、95%の信頼区間、標本から予測される平均値が出力される。もちろん母集団の平均は男女で変えて行っているので $p=0.05$ の水準で帰無仮説（平均値が男女で一致する）は棄却される。
var.testを用いた検定
分散の等質性の検定
通常は2群の平均値の比較を行う前に分散が等質であるかの検定を行う。そのためのRコマンドがvar.testである。先ほど $t$ 検定で用いた男女の身長データに対して分散が等質であるかを調べてみる。var.testもデータフレームの入力を受け付けないので、列番号で指定する。
var.test(male[,1], female[,1])
F test to compare two variances
data: male[, 1] and female[, 1]
F = 1.3458, num df = 999, denom df = 999, p-value = 2.822e-06
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
1.188720 1.523576
sample estimates:
ratio of variances
1.345773
データ、統計量・自由度・p値、対立仮説、95%の信頼区間、標本から予測される分散の比が出力される。もちろん母集団の分散も男女で変えて行っているので $p=0.05$ の水準で帰無仮説（分散の比が１になる）は棄却される。
cor.testを用いた検定
無相関検定（Pearson's product-moment correlation）
相関係数 $r$ を用いて計算される統計量
t = \frac{r\sqrt{n-2}}{\sqrt{1-r^2}}
は、帰無仮説のもとで自由度 $df=n-2$ の $t$ 分布に従う。この性質を利用して相関の有無を検定するRコマンドがcor.testである。
cor.testを用いて、データの身長と体重に相関があるのか調べてみる。データはdat["WEIGHT"] などで入れるとデータフレーム扱いになってしまうので、列番号で指定した。
cor.test(dat[,1], dat[,2])
Pearson's product-moment correlation
data: dat[, 1] and dat[, 2]
t = 36.58, df = 1998, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.6063250 0.6588662
sample estimates:
cor
0.6333248
データ、統計量・自由度・p値、対立仮説、95%の信頼区間、標本から予測される相関係数が出力される。身長に基づきBMIから体重を作成しているので二つの数字に相関はある。検定では、 $p=0.05$ の水準で帰無仮説（身長と体重は無相関である）は棄却されている。
その他の無相関検定
cor.testのオプションで、method="spearman", method="kendall"を選べば、それぞれ「スピアマンの順位相関係数の無相関検定」、「ケンドールの順位相関係数の無相関検定」を行うことができる。因みにデフォルトで設定されているのはタイトルの通り、method="pearson"で「ピアソンの積率相関係数の無相関検定」である。
chisq.testを用いた検定
独立性の検定
質的変数についてクロス集計表を取った場合、その観測度数 $O_i$ 、期待度数 $E_i$ から計算される統計量である
\chi^2 = \sum_{i}^k \frac{\left(O_i - E_i\right)^2}{E_i}
は帰無仮説のもとで $\chi^2$ 分布に従う。ただし自由度はクロス集計表の「（行の数）-1」×「（列の数）-1」である。この性質を利用して独立性の検定を行うのがchisq.testである。
chisq.testを用いて性別と血液型に連関があるのか調べてみる。
tab2 <- table(dat2)
chisq.test(tab2)
Pearson's Chi-squared test
data: tab2
X-squared = 7.5932, df = 3, p-value = 0.05521
データ、統計量・自由度・p値、が出力される。検定では、 $p=0.05$ の水準で帰無仮説（性別と血液型は連関が無い）は棄却されないので、それぞれ独立なデータであると言える。