【R言語】母比率の検定・比率の差の検定

R言語で母比率の検定を実行する方法を紹介します。

R言語ではt検定などと同様に、簡単に比率に関する検定を実行することが可能です。

今回は、母比率の検定に関する関数の説明やその実行例をまとめました。

この記事で紹介するプログラミングコードは以下のzipファイル中のRスクリプトにも保存しています。

R言語　母比率の検定　Rスクリプト

1 ファイル 1.32 KB

ダウンロード

また、母比率の検定の検定の手順を以下の記事で解説しています。

: 【統計学】母比率の検定・母比率の差の検定
通常の1標本の母比率の検定や2標本の母比率の差の検定を含む様々な母比率の検定を解説する。ここでは、Z検定を用いた大標本の下で母比率の検定について紹介し、検定統計量や棄却域の導出を行う。大標本におけ ...
続きを見る

R言語ではじめるプログラミングとデータ分析

Rでデータ解析を始めるならコレ

¥3,300 （2022/06/01 17:20時点 | Amazon調べ）

Amazon

楽天市場

Yahooショッピング

ポチップ

改訂2版 RユーザのためのRStudio[実践]入門〜tidyverseによるモダンな分析フローの世界

RStudioで解析したい方へ

¥3,278 （2022/06/02 22:52時点 | Amazon調べ）

Amazon

楽天市場

Yahooショッピング

ポチップ

母比率の検定

母比率の検定及び母比率の差の検定の概要とR言語の関数を以下にまとめました。

概要

母比率の検定と母比率の差の検定の検定統計量と棄却域は次のようになります。

母比率の検定

$x_1, \ldots, x_n$は$Bernoulli(p)$からの独立同一な標本であり、標本数$n$は十分に大きいと仮定する。このとき、次の「母比率$p$は特定の値$p_0$であるか」の仮説を検定する。

\begin{align}&H_0: p= p_0 \\&H_1: p\neq p_0\end{align}

検定統計量として次を用いる。

\begin{align} Z= \cfrac{\hat{p}- p_0}{\sqrt{ p_0 (1 - p_0) / n}} \sim N(0, 1). \end{align}

また、有意水準$\alpha$の検定の棄却域は次で与えられる。

\begin{align}(-\infty,- Z_{\alpha/2}) \cup (Z_{\alpha/2}, \infty),\end{align}

ここに、$\hat{p}$は次で与えられる母比率の推定値の確率変数であり、$Z_{\alpha/2}$は標準正規分布の上側$\alpha/2$の確率点である。

\begin{align}\hat{p} = \cfrac{1}{n}\sum_{i=1}^n X_i.\end{align}

母比率の差の検定

$x_{11} ,\ldots, x_{1n_1}$はパラメータ$p_1$のベルヌーイ分布$Bernolli(p_1)$からの独立同一な標本、$x_{21} ,\ldots, x_{2n_2}$はパラメータ$p_2$のベルヌーイ分布$Bernolli(p_2)$からの独立同一な標本とする。このとき次の仮説検定を考える。

\begin{align}&H_0:\ p_1= p_2\\ &H_1:\ p_1\neq p_2\end{align}

検定統計量として次を用いる。

\begin{align}Z = \cfrac{\hat{p}_1 -\hat{p}_2}{\sqrt{\hat{p}(1-\hat{p}) (1/n_1 + 1/n_2)}} \sim N(0,1),\end{align}

ここに

\begin{align}\hat{p} &= \cfrac{n_1\hat{p_1} + n_2\hat{p}_2}{n_1 + n_2},\\ \hat{p}_i &= \cfrac{1}{n_i}\sum_{j=1}^{n_i} x_j, i = 1,2.\end{align}

また、有意水準$\alpha$の棄却域は次で与えられる。

\begin{align}(-\infty, -Z(\alpha/2)) \cup (Z(\alpha/2), \infty),\end{align}

ここに、$Z(\alpha/2)$は標準正規分布の上側$\alpha/2$点である。

関数prop.test

母比率の検定を実行する関数であるprop.testについて紹介します。関数prop.testの詳細は以下の通りです。

prop.test(x, n, p = NULL, alternative = c("two.sided", "less", "greater"), conf.level = 0.95, correct = TRUE)

関数prop.testの引数
x	成功回数を表すnumeric型のベクトル。または成功と失敗の回数から成る2次行列。
n	試行回数を表すnumeric型のベクトル。xが行列である場合無視される。
p	成功する確率。ベルヌーイ分布のパラメータ。
alternative	対立仮説を指定するために用いる。"two.sided"で両側検定、"greater"で右片側検定、"less"で左片側検定を実行できる。
conf.level	信頼区間の信頼度。
correct	TRUEで連続補正を適用する。FALSEで適用しない。

実行例

関数prop.testを用いて母比率の検定と2群の母比率の差の検定の実行例を紹介します。

母比率の検定

まず1標本検定である母比率の検定についてみていきます。

今、次のような100回のベルヌーイ試行によるデータが得られたとします。

> samples <- sample(c(0, 1), 100, replace = TRUE) #100回のベルヌーイ試行
> head(samples, 20)
 [1] 1 1 1 1 0 1 0 1 0 0 1 1 0 0 0 1 0 1 0 0

> samples <- sample(c(0, 1), 100, replace = TRUE) #100回のベルヌーイ試行

> head(samples, 20)

[1] 1 1 1 1 0 1 0 1 0 0 1 1 0 0 0 1 0 1 0 0

samplesは0と1の二値のデータから成ります。関数tableを用いると度数分布表を作成することができ、0と1の数を把握することが可能です。

> sampleTable <- table(samples)
> sampleTable
samples
 0  1 
52 48

> sampleTable <- table(samples)

> sampleTable

samples

0 1

52 48

この0と1が発生する確率はそれぞれ等しいか（p=0.5であるか）という仮説検定を実行したいと思います。

この仮説検定を実行するには、以下のようにprop.testの第一引数に0が発生した数、第二引数に試行回数を指定します。

2行目のように度数分布表を渡す方法でも実行することができます。

testResult <- prop.test(sampleTable["0"], length(samples), correct = FALSE)
prop.test(sampleTable, n = length(samples), correct = FALSE) #上記と同じ結果

1 2	testResult <- prop.test(sampleTable["0"], length(samples), correct = FALSE) prop.test(sampleTable, n = length(samples), correct = FALSE) #上記と同じ結果

testResultを参照すると、母比率の検定の検定統計量やp値などの結果がコンソールに出力されます。

> testResult

	1-sample proportions test without continuity correction

data:  sampleTable["0"] out of length(samples), null probability 0.5
X-squared = 0.16, df = 1, p-value = 0.6892
alternative hypothesis: true p is not equal to 0.5
95 percent confidence interval:
 0.4231658 0.6153545
sample estimates:
   p 
0.52 

>

> testResult

1-sample proportions test without continuity correction

data: sampleTable["0"] out of length(samples), null probability 0.5

X-squared = 0.16, df = 1, p-value = 0.6892

alternative hypothesis: true p is not equal to 0.5

95 percent confidence interval:

0.4231658 0.6153545

sample estimates:

0.52

ポイント

「p-value = 0.6892 > 0.05」より有意水準0.05の下で、このベルヌーイ試行の成功確率は0.5であることがいえました。（成功と失敗の確率が異なるとはいえない）

検定結果を取得したい場合、以下のようにtestResultの後ろに「$」を付けます。

検定統計量や、p値、信頼区間といった検定において重要な指標を取得することができます。

stat <- testResult$statistic     #検定統計量
pValue <- testResult$p.value     #p値
estimator <- testResult$estimate #母比率の推定値
confInt <- testResult$conf.int   #母比率の信頼区間

stat <- testResult$statistic #検定統計量

pValue <- testResult$p.value #p値

estimator <- testResult$estimate #母比率の推定値

confInt <- testResult$conf.int #母比率の信頼区間

これらの検定結果をcsvファイルへと保存する例を最後に紹介します。

関数data.frameを使いこれらの結果をデータフレームにまとめた後に、write.csvを用いることで以下の画像のように結果を保存することができます。

resultTable <- data.frame(stat = stat, p.value = pValue, estimate = estimator, 
                          c.i.lower = confInt[1], c.i.upper = confInt[2], row.names = NULL)
write.csv(resultTable, "母比率の検定.csv", row.names = FALSE)

resultTable <- data.frame(stat = stat, p.value = pValue, estimate = estimator,

c.i.lower = confInt[1], c.i.upper = confInt[2], row.names = NULL)

write.csv(resultTable, "母比率の検定.csv", row.names = FALSE)

母比率の検定結果

母比率の差の検定

最後に母比率の差の検定の実行方法について紹介します。

母比率の差の検定の検定統計量と棄却域は次のようになります。

母比率の差の検定は母比率の検定と同様に、関数prop.testを用いることでできます。

まず、データセットを用意します。パッケージepitoolsをインストールし、データセットTitanicをデータフレームへと変換します。

library(epitools)
dataset <- expand.table(Titanic)
dataset$Survived <- factor(dataset$Survived, levels = c("Yes", "No"))

library(epitools)

dataset <- expand.table(Titanic)

dataset$Survived <- factor(dataset$Survived, levels = c("Yes", "No"))

datasetは次のようなClass、Sex、Age、Survivedの4つの項目から成ります。

> head(dataset)
  Class  Sex   Age Survived
1   1st Male Child      Yes
2   1st Male Child      Yes
3   1st Male Child      Yes
4   1st Male Child      Yes
5   1st Male Child      Yes
6   1st Male Adult       No

> head(dataset)

Class Sex Age Survived

1 1st Male Child Yes

2 1st Male Child Yes

3 1st Male Child Yes

4 1st Male Child Yes

5 1st Male Child Yes

6 1st Male Adult No

SurvivedとSexに関する二元分割表は次のようにして作れます。

> twoWayTable <- table(dataset$Sex, dataset$Survived)
> twoWayTable
        
          Yes   No
  Male    367 1364
  Female  344  126

> twoWayTable <- table(dataset$Sex, dataset$Survived)

> twoWayTable

Yes No

Male 367 1364

Female 344 126

ポイント

この表から、Female（女性）よりMale（男性）の方が生存率が低いことが推測できます。

女性のNoの割合が全体の半分未満なのに対して、男性のNoの割合がほとんどを占めているのが分かります。

男性の生存率と女性の生存率が等しいかどうか仮説検定をしていきたいと思います。

母比率の差の検定を行うには、関数prop.testの第一引数に各グループの生存数、第二引数に各グループの総数を指定します。

2行目のように、行列を引数に渡す方法でも母比率の差の検定を実行できます。

testResult <- prop.test(twoWayTable["Yes", ], colSums(twoWayTable), correct = FALSE)
prop.test(twoWayTable, correct = FALSE) #上記と同じ結果

1 2	testResult <- prop.test(twoWayTable["Yes", ], colSums(twoWayTable), correct = FALSE) prop.test(twoWayTable, correct = FALSE) #上記と同じ結果

testResultを参照すると、母比率の差の検定の結果が出力されます。

> testResult

	2-sample test for equality of proportions without
	continuity correction

data:  twoWayTable["Yes", ] out of colSums(twoWayTable)
X-squared = 456.87, df = 1, p-value < 2.2e-16
alternative hypothesis: two.sided
95 percent confidence interval:
 -0.5643339 -0.4754635
sample estimates:
   prop 1    prop 2 
0.2120162 0.7319149

> testResult

2-sample test for equality of proportions without

continuity correction

data: twoWayTable["Yes", ] out of colSums(twoWayTable)

X-squared = 456.87, df = 1, p-value < 2.2e-16

alternative hypothesis: two.sided

95 percent confidence interval:

-0.5643339 -0.4754635

sample estimates:

prop 1 prop 2

0.2120162 0.7319149

ポイント

testResultを参照するとp値<2.2e-16より、男性の生存率と女性の生存率は等しいとは言えないことが分かりました。

性別間で生存率に差があることが仮説検定により示されました。

testResult中の検定統計量や、p値、信頼区間などの検定結果のは以下のようにして取得できます。

stat <- testResult$statistic     #検定統計量
pValue <- testResult$p.value     #p値
estimator <- testResult$estimate #p1とp2の推定値
confInt <- testResult$conf.int   #p1-p2の信頼区間

stat <- testResult$statistic #検定統計量

pValue <- testResult$p.value #p値

estimator <- testResult$estimate #p1とp2の推定値

confInt <- testResult$conf.int #p1-p2の信頼区間

関数data.frameを用いてこれらの結果をデータフレームに格納し、write.csvを用いることで次の画像のように検定結果をcsvファイルへ保存することが可能です。

resultTable <- data.frame(stat = stat, p.value = pValue, 
                          sample.prop.1 = estimator[1], sample.prop.2 =  estimator[2],
                          c.i.lower = confInt[1], c.i.upper = confInt[2], row.names = NULL)
write.csv(resultTable, "母比率の差の検定.csv", row.names = FALSE)

resultTable <- data.frame(stat = stat, p.value = pValue,

sample.prop.1 = estimator[1], sample.prop.2 = estimator[2],

c.i.lower = confInt[1], c.i.upper = confInt[2], row.names = NULL)

write.csv(resultTable, "母比率の差の検定.csv", row.names = FALSE)

母比率の差の検定結果

まとめ

R言語で母比率の検定を実行する関数とその実行例を

標準で母比率の検定の関数が実装されており、簡単に検定を実施することができます。

また、t検定と同様に母比率の差（2群の差）の検定も同じ関数を用いて実行することが可能です。

usagi-san

統計学とゲームとかをメインに解説していくよ。数式とかプログラミングコードにミスがあったり質問があったりする場合はコメントで受け付けます。すぐに対応します。

2024/04/28

【R言語】F検定テンプレートスクリプト

2024/04/28

【R言語】経験密度関数・経験分布関数のプロット　関数densityとecdfの使い方

2024/04/10

【R言語】ベータ関数とガンマ関数　関数beta, gammaの使い方

usagi-sanの記事をもっと見る

-R言語
-R言語, プログラミング, 統計学, 統計解析

comment コメントをキャンセル

: R言語
【R言語】factor型の使い方・まとめ
この記事では、R言語特有のfactor型（因子型）の操作方法について解説していきます。この記事では以下の操作方法を扱っています。 factorの作成からlevelsの変更などを紹介しています。この ...

: R言語
【R言語】指数関数・ネイピア数
R言語でネイピア数eおよび指数関数について解説します。この記事では関数expを使った指数関数の計算の仕方を見ていきます。自然対数については次の記事を参照ください。この記事で紹介するコードは以 ...

: R言語
R言語　自動でデータ整形してくれる関数【データ整形不要】
こんにちは、usagi-sanです。今回はデータ解析に役立つ関数を作ってみました。自分がアルバイトをする中で、データ整形（データクレンジング）に用いる時間が勤務時間の9割に及び、肝心な統計解析にミ ...

: R言語
【R言語】関数lmによる線形回帰　単回帰分析・重回帰分析　回帰分析で予測する
R言語の回帰分析の実行方法を紹介していきます。回帰分析の実行方法だけでなく、回帰分析により作ったモデルでテストデータの予測も行っていきます。実行方法では、回帰分析の結果の見方やその抽出方法を詳しく解 ...

: R言語
【R言語】中央値・四分位数・最大値・最小値の計算　記述統計
R言語で中央値や四分位数、最大値、最小値を求める関数やその実行例を解説します。データの中央値や四分位数、最大値、最小値の計算の例や、これらの要約統計量と箱ひげ図の関係も紹介します。データの要約の流れ ...

独立なカイ2乗分布の積の分布【統計学】

【統計学】母比率の検定・母比率の差の検定

【R言語】母比率の検定・母比率の差の検定

【統計学】母比率の検定・母比率の差の検定

母比率の検定

概要

関数prop.test

実行例

母比率の検定

母比率の差の検定

まとめ