「自然科学の統計学」第6章演習問題1の二項分布に対する検定の検定力に関する問題を解いて、Pythonで可視化してみた

元教師

こんにちは！データサイエンティストの青木和也（https://twitter.com/kaizen_oni）です！

今回の記事では、統計学の青本「自然科学の統計学」の第6章-演習問題1「二項分布に対する検定の検定力」を解いた上で、検出力関数のグラフをPythonで可視化をしていきたいと思います。

検出力、検出力関数あたりの話は馴染みのない方には理解しにくいかもしれませんので、適宜解説を入れています。

検出力まわりについて聞いたことがない方も「なるほどね！」と理解いただきながら読み進めていただけると幸いです！

問題文
解答
まとめ

問題文

<二項分布に対する検定の検出力>

コインにゆがみがないかどうかを検定するのに、そのコインを4回投げて4回とも表が出るか４回とも裏が出たらゆがみがあるとする検定方式を考える。

この検定方式の(1)検出力関数を求め(2)グラフにせよ
東京大学教養学部統計学教室『自然科学の統計学』(東京大学出版社/2001) 第6章演習問題 P199

帰無仮説・対立仮説とは

帰無仮説も対立仮説も、どちらも何らかしらの仮説です。

例えば、コインに歪みがある/ない・男女間で給与水準に差がある/ないなどが挙げられます。

そして、対立仮説の方が本当に言いたい仮説（今回であればコインに歪みがある）であることが多いです。

一方で、帰無仮説の方は「帰無（＝無に帰す）」とあるように、本当に言いたい仮説の話をするために否定したい（=無に帰したい）仮説であることが多いです。

今回の場合でいえば、帰無仮説と対立仮説は以下のように整理できます。

帰無仮説：コインに歪みはない

対立仮説：コインに歪みはある

検定力関数とは

検定力関数とは、ある検定方式Aを用いて、$\theta$という値に対して帰無仮説を棄却する（=対立仮説を採択する）確率$\beta_A (\theta)$のことを指します。

$$ \beta_A (\theta) = P_{\theta}(検定方法Aを用いて帰無仮説を棄却)$$

$\theta$には帰無仮説を棄却する確率を求めるために使うさまざまな値が入ります。

以下はその1例です。

コイン投げ（二項分布）：コインの表が出る確率$p$
ある高校の身長（正規分布）：身長の平均値$\mu$
ある高校の男女の身長比較（正規分布）：男女の身長の平均値$\mu_{男},~\mu_{女}$

解答

それでは今回の問題に対する解答を順を追って説明していきましょう。

検定力関数を求める

検定力関数を求めるために以下の3ステップで考えを進めていきます。

帰無仮説・対立仮説を整理する
帰無仮説を棄却する時の状況を整理する
検定力関数を求める
グラフを書く

順を追って解説していきます。

帰無仮説・対立仮説を整理する

今回の問題では「コインに歪みがあるかどうか」を検定したいとのことでした。

なので、「コインに歪みがある」というのが本当に言いたい仮説だと考えます。

すると、帰無仮説・対立仮説は以下のように整理できます。

帰無仮説：コインに歪みはない

対立仮説：コインに歪みはある

帰無仮説を棄却する時の状況を整理する

帰無仮説・対立仮説を整理したところで、検定力関数、つまり「帰無仮説を棄却する関数」を求めるために、「帰無仮説を棄却するときってどんな時？」というのを考えてみましょう。

ここで、問題文には以下のように書かれています。

コインにゆがみがないかどうかを検定するのに、そのコインを4回投げて4回とも表が出るか４回とも裏が出たらゆがみがあるとする検定方式を考える。
東京大学教養学部統計学教室『自然科学の統計学』(東京大学出版社/2001) 第6章演習問題 P199

以上の一文から、帰無仮説を棄却するのは以下の2パターンだと考えることができます。

コインを４回投げた時、４回とも表が出た時
コインを４回投げた時、４回とも裏が出た時

検定力関数を求める

帰無仮説を棄却する状況について理解を深めることができたので、実際に検定力関数を求めてみましょう。

検定力関数を求めるための前提条件を2つほど整理しておきます。

コイン投げの試行は二項分布と考えることができる
コインの表が出る確率を$p$とおく

帰無仮説を棄却する状況は以下の2パターンでした。

コインを４回投げた時、４回とも表が出た時
コインを４回投げた時、４回とも裏が出た時

1の状況について考えてみると、

表が出る確率は$p$なので、確率$p$が4回の中から$4$回選ばれることは以下の式で整理されます。

$$ (４回表が出る確率) = {}_4\mathrm{C}_{4} p^4 = p^4$$

2の状況について考えてみると、

裏が出る確率は$1-p$なので、確率$1-p$が4回の中から$4$回選ばれることは以下の式で整理されます。

$$ (４回裏が出る確率) = {}_4\mathrm{C}_{4} (1-p)^4 = (1-p)^4$$

この2つの事象は排反事象、つまり同時に起こらないことであるので、足し算をすることができます。

よって、この検定方法（4回表か裏が出たら歪みあり）をAとすると、求める検定力関数は以下のような式で表されます。

$$ \beta_{A}(\mu) = p^4 + (1-p)^4$$

グラフを書く

求めた検定力関数のグラフを書いていきましょう。

この時、検定力関数の変数はコインの表が出る確率$p$です。

確率$p$が$0〜1$の値を動いた時の検定力関数（コインに歪みがある場合に本当に歪みを見抜ける確率）がどのような値を動くのかをみていきましょう。

グラフの可視化コード

グラフの可視化に際しては以下のようなコードで可視化をすることができます。

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import pandas as pd

# コインの表が出る確率を0~1の間で100個作成
p = np.linspace(0,1,100)
# それぞれの確率における検定力を求める
beta = p**4 + (1-p) ** 4

# Seabornで可視化するためにデータフレームに格納
df = pd.DataFrame(columns = ['Probability', 'Power Function'])
df['Probability'] = p
df['Power Function'] = beta

# 検定力の最小値を求めておく
min_num = min(beta)

# グラフの可視化、ついでに検定力関数の最低値も可視化
sns.lineplot(df, x = 'Probability', y = 'Power Function')
plt.axhline(min_num, color  ='gray', ls = 'dotted')
plt.text(-0.14, min_num-0.01, f'{min_num: .2f}')

上記の図から、以下のようなことがわかります。

表が出る確率が100%または裏が出る確率が100%の時は完璧に検出できる
表が出る確率と裏が出る確率が100%から低下していく時は検出力も低下していき、$p=0.5$の時に最小値を取る
コインに歪みがない時（$p=0.5$）のときも13％の確率で「コインに歪みがある」と誤判定してしまう

実装コードでは以下のステップで検定力関数を可視化しています。

コインを投げて表が出る確率$p$と検定力のリストを作成し、データフレームに格納する
検定力の最小値を求める
Seabornライブラリを使用して可視化を行う

順を追って解説をしていきます。

コインを投げて表が出る確率$p$と検定力のリストを作成し、データフレームに格納する

# コインの表が出る確率を0~1の間で100個作成
p = np.linspace(0,1,100)
# 検定力関数を定義
beta = p**4 + (1-p) ** 4

# Seabornで可視化するためにデータフレームに格納
df = pd.DataFrame(columns = ['Probability', 'Power Function'])
df['Probability'] = p
df['Power Function'] = beta

np.linspaceを使用して、0以上1以下の確率を等間隔に100個作成していきます。

上記で作成した確率を$\beta_{A}(p) = p^4 + (1-p)^4$にそれぞれ代入して、検定力も同様に100個作成します。

最後にSeabornで可視化を行うために、コイン投げの確率(Probability)と検定力(Power Function)をデータフレームに格納します。

検定力の最小値を求める

問題の主題ではありませんが、帰無仮説を誤って棄却してしまう（=コインが歪んでいないにもかかわらず、コインが歪んでいると判定してしまう）確率も求めておきたいので、検定力の最小値を求めておきます。

# 検定力の最小値を求めておく
min_num = min(beta)

なお、このような「帰無仮説が正しいにもかかわらず、帰無仮説を誤って棄却してしまう」ことを第一種の過誤と呼びます。

Seabornライブラリを使用して可視化を行う

Seabornライブラリのlineplotメソッドを利用して、コイン投げの確率と検定力関数のグラフを出力します。

ついでに検定力の最小値についても水平線で可視化をしておきます。

# グラフの可視化、ついでに検定力関数の最低値も可視化
sns.lineplot(df, x = 'Probability', y = 'Power Function')
plt.axhline(min_num, color  ='gray', ls = 'dotted')
plt.text(-0.14, min_num-0.01, f'{min_num: .2f}')