Rで行う基本的な統計検定とその選択フローチャート1
今回は生物学で用いられている基本的な統計検定をRで行う際の手順と実際のコードを特集します。自身のデータセットの特徴と実験目的から適切な統計検定法を選択し、それをRで実行できるようになるための指針となれば幸いです。
第一回:平均値の差の検定(2群)
第一回:平均値の差の検定(2群)
第二回:一元配置の分散分析
第三回:一元配置の多重比較
第四回:二元配置の分散分析とその後の多重比較
第五回:比率の差の検定
第六回:サンプルサイズの設計法
第七回:回帰分析
第一回 平均値の差の検定(2群)
1.1 適切な検定法の選択
自分の取得した標本が対応のあるデータなのか、正規分布に従うのか、2標本の分散が等しいか、を検証して適切な統計検定法を選択できるような手順と必要となる検定のソースコード(R)を紹介します。
1.1.1 ヒストグラムとQQプロット
1.1.2 Shapiro-Wilk test
1.1.3 Kolmogorov-Smirnov test
1.1.4 F-test
1.2 パラメトリック検定法のRソースコード
パラメトリック検定のソースコードを紹介します。
1.2.1 Paired T-test(対応のあるデータの検定)
1.2.2 Student's T-test
1.2.3 Welch's T-test
1.3 ノンパラメトリック検定法のRソースコード
ノンパラメトリック検定のソースコードを紹介します。
1.3.1 paired wilcoxon's signed-rank test
1.3.2 wilcoxon's signed-rank test
1.1 適切な検定法の選択
まずは最も頻出の2標本の平均値の差の検定についてです。以下に検定法選択のためのフローチャートをのせておきます。選択の上で重要なのが、測定値について(i)母集団が正規分布に従うと仮定できるか:正規性、と(ii)母集団の分散が標本間で等しいか:等分散性、の2つになります。
1.1.1 正規性の検証(可視化)
データが正規分布しているか否かはまずヒストグラムを描いて分布が釣鐘状になっているか確認します。
あるいはQQプロットを描く方法もあります。QQプロットは(作成中)。
データの分布がどうなっているかを可視化しておくことは正規性の確認だけでなく、データの特徴を掴む上で重要ですのでデータを取得した際には必ず行うようにしましょう。
1.1.2 正規性の検証(Shapiro-Wilk test)
上にあげた手法は直感的で分かりやすい一方で、どの程度ならば正規分布といっていいかの判断が難しい場合があります。その場合は、正規性の確認を統計検定によって判断することもできます。
Shapiro-Wilk testは一般に大サンプルでない場合(2000以下?)に用いられる正規性検定です。
Rで実際に検定を行なってみましょう。帰無仮説は(H0)は標本分布が正規分布に従うことです。
score=c(57, 67, 23, 50, 52, 51, 48, 49, 54, 44, 31, 53, 40, 60)
shapiro.test(score)
結果
Shapiro-Wilk normality test
data: score
W = 0.93902, p-value = 0.4059
以上から帰無仮説は棄却されず、データは正規分布に従うと仮定してよいと考えられます。
参照:https://data-science.gr.jp/implementation/ist_r_shapiro_wilk_test.html
1.1.3 正規性の検証(Kolmogorov-Smirnov test)
一方でデータが大サンプル(2000以上?)である場合は、Kolmogorov-Smirnov testを行うことが推奨されている。
Rで実際に検定を行なってみましょう。帰無仮説は(H0)は標本分布が正規分布に従うことです。
今回は大標本ではなく、上記のShapiro-Wilk testと同じデータセットで検定してみましょう。
score=c(57, 67, 23, 50, 52, 51, 48, 49, 54, 44, 31, 53, 40, 60)
ks.test(x=score,y="pnorm",mean=mean(score),sd=sd(score))
結果
One-sample Kolmogorov-Smirnov test
data: score
D = 0.19668, p-value = 0.5839
alternative hypothesis: two-sided
以上から帰無仮説は棄却されず、データは正規分布に従うと仮定してよいと考えられます。
参照:https://data-science.gr.jp/implementation/ist_r_kolmogorov_smirnov_test.html
Shapiro-Wilk testと比べるとp値が異なります。
Shapiro-Wilk testとKolmogorov-Smirnov testはサンプル数によって使い分けることが推奨されています。サンプル数がどの程度であればどちらの検定を使うべきかについては今回は割愛します。
詳しくは-> https://blog.statsbeginner.net/entry/2014/08/13/115744 や https://qiita.com/purple_jp/items/f22b0e3fec8f2d548d75などを参照してください。
生物学における実験では大標本はあまり得ることができないので基本的にはShapiro-Wilk testを行うものと考えてください。
1.1.4 等分散性の検証(F-test)
データが正規分布に従いそうと分かれば、次は等分散性を検証する。もしデータに正規性がない場合はノンパラメトリック手法を使う。等分散性の検証はWelchのt検定かstudentのt検定かのどちらを使うべきかの選択に関わる。2群の分散が等しいかどうかを検定するにはF-testを用いる。多群の場合はBartlet testを用いる(詳細は多重比較で説明)。
Rで実際に検定を行なってみましょう。帰無仮説は(H0)は2標本の分散が等しいことです。
va=c(301, 311, 325, 291, 388, 412, 325, 361, 287)
vb=c(197, 180, 247, 260, 247, 199, 179, 134, 163, 200)
var.test(va,vb)
結果
F test to compare two variances
data: va and vb
F = 1.2001, num df = 8, denom df = 9, p-value = 0.7859
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
0.2925646 5.2290589
sample estimates:
ratio of variances
1.200087
以上から帰無仮説は棄却されず、2標本の分散は等しいと仮定してよいと考えられます。
参照:https://data-science.gr.jp/implementation/ist_r_f_test.html
1.2 検定手法をRで実行してみる(パラメトリック編)
母集団が正規分布に従うと仮定できる場合、パラメトリックな検定手法を使って、平均値が2群間で差があるか検定することになります。
1.2.1 Paired T-test(対応のあるデータの検定)
対応のあるデータに対しては、paired t-testを用います。
さっそくRで実行してみましょう。帰無仮説は(H0)は対応のある2標本の平均値に差がないことです。
va=c(301, 311, 325, 291, 388, 412, 325, 361, 287)
vb=c(197, 180, 247, 260, 247, 199, 179, 134, 163)
t.test(va,vb,paired=TRUE)
結果
Paired t-test
data: va and vb
t = 6.524, df = 8, p-value = 0.0001834
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
85.84529 179.71026
sample estimates:
mean of the differences
132.7778
以上から帰無仮説は棄却され、2標本の平均値には差があるといえます。
1.2.2 Student's T-test
正規性に従うが、対応がなく、母分散が等しいデータに対してはstudent's t-testを用います。
Rで実行してみましょう。帰無仮説は(H0)は2標本の平均値に差がないことです。
va=c(301, 311, 325, 291, 388, 412, 325, 361, 287)
vb=c(197, 180, 247, 260, 247, 199, 179, 134, 163, 200)
t.test(va,vb,var.equal=T,paired=F)
結果
Two Sample t-test
data: va and vb
t = 6.8649, df = 17, p-value = 2.743e-06
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
92.01697 173.67192
sample estimates:
mean of x mean of y
333.4444 200.6000
以上から帰無仮説は棄却され、2標本の平均値には差があるといえます。
さらにvaがvbより大きいことを検定するには、片側検定を使います。片側検定を行いたい場合は、引数にalternative = "less"またはalternative = "greater"を指定しておきます。
片側検定の帰無仮説は先ほどと同じく、2標本の平均値に差がないことですが、 対立仮説(H1)は2標本の平均値に差がある、ではなく、vaの平均値はvbよりも高いとなります。
t.test(va,vb,var.equal=T,paired=F,alternative = "greater")
片側検定の結果
Two Sample t-test
data: va and vb
t = 6.8649, df = 17, p-value = 1.372e-06
alternative hypothesis: true difference in means is greater than 0
95 percent confidence interval:
99.18096 Inf
sample estimates:
mean of x mean of y
333.4444 200.6000
以上から帰無仮説は棄却され、vaの平均値はvbよりも高いといえます。
1.2.3 Welch's T-test
正規性に従うが、対応がなく、母分散が等しくないデータに対してはstudent's t-testを用います。
Rで実行してみましょう。帰無仮説は(H0)は2標本の平均値に差がないことです。
va=c(301, 311, 325, 291, 388, 412, 325, 361, 287)
vb=c(197, 180, 247, 98, 227, 199, 79, 134, 163, 200)
wilcox.test(va, vb, paired=FALSE)
結果
以上から帰無仮説は棄却され、2標本の平均値には差があるといえます。
student's t-testと同様に、片側検定を行うことができます。
1.3 検定手法をRで実行してみる(ノンパラメトリック編)
1.3.1 paired wilcoxon's signed-rank test
Rで実行してみましょう。帰無仮説は(H0)は2標本の平均値に差がないことです。
va=c(341, 311, 325, 391, 388, 392, 325, 361, 387)
vb=c(210, 69, 247, 98, 247, 169, 79, 114, 183)
wilcox.test(va,vb, paired=TRUE)
結果
Wilcoxon signed rank exact test
data: va and vb
V = 45, p-value = 0.003906
alternative hypothesis: true location shift is not equal to 0
以上から帰無仮説は棄却され、2標本の平均値には差があるといえます。
1.3.2 wilcoxon's signed-rank test
Rで実行してみましょう。帰無仮説は(H0)は2標本の平均値に差がないことです。
va=c(341, 311, 325, 391, 388, 392, 324, 361, 387)
> vb=c(197, 88, 247, 98, 227, 199, 79, 134, 163, 200, 289, 89, 319)
> wilcox.test(va, vb, paired=FALSE)
結果
Wilcoxon rank sum exact test
data: va and vb
W = 116, p-value = 8.041e-06
alternative hypothesis: true location shift is not equal to 0
以上から帰無仮説は棄却され、2標本の平均値には差があるといえます。