はじめに
本稿は以下のBingチャットでの一節について確かめるためのメモです.
「p値の分布が一般的には一様分布となる」
というのが本当なのかどうか直感的に理解できませんでした.
同一のガンマ分布から実際本当に無作為抽出した２つの標本でこの検定を何度も行うと、p値の分布は一般的には一様分布になります。 ただし、標本数が少ない場合は、p値の分布が一様分布にならないことがあります。(Bing)
メモ
そもそもp値とはどうやって計算するのかをBingに聞いてみました.
p値とは、帰無仮説が正しいと仮定した時に、実際に観測されたデータ以上に対立仮説が正しいことを支持するデータが観測される確率のことです。p値が小さいほど、帰無仮説が正しいという仮説を棄却することができます。例えば、p値が0.05以下であれば、そのデータは有意差があると判断されます。(Bing)
そうなんですね.
それでは、帰無仮説が正しい時にある確率変数が正規分布に従う場合、p値はどのように計算されるのでしょうか.
正規分布に従う確率変数のp値を求めるには、標準正規分布表を用いる方法があります1。標準正規分布表は、平均が0、標準偏差が1の正規分布に従う確率変数Zの累積分布関数P(Z≦z)を表したものです。まず、標準正規分布表P(Z≧z)のz値を求めます。次に、求めたz値を使って、P(Z≧z)を求めます。最後に、P(Z≧z)を2倍して、p値を求めます。例えば、平均が0、標準偏差が1の正規分布に従う確率変数Zがあるとき、Z=1.96であれば、P(Z≧1.96)=0.025であり、p値は0.05です。(Bing)
説明としてはかなり違和感を感じるものの、雰囲気を掴めてきたかもしれません.
正規分布は左右対称なので二倍すれば良いということでしょうか.
はい、正規分布は左右対称ですので、p値を求めるためには、P(Z≧z)を2倍することができます。(Bing)
そうするとP(Z>z)が0.5のときはp値は1ということですね.
また、p値の取りうる範囲は0から1までということでしょうか.
P(Z>z)が0.5のとき、p値は1です。
はい、p値は0から1までの範囲で取りうる値です。(Bing)
統計検定2級の記憶を頼りにすると、標準正規分布からのが実現値zが0を境にp値の計算方法が変化しそうです. あるいは絶対値をとって計算する気がします.
Bingに質問することで、大雑把ではありますがp値に関する良い復習になりました.
こんなところで良いでしょう.
p値についてグラフで確かめるには標準正規分布を用いた方がよさそうです.
下図は標準正規分布の確率密度関数と累積分布関数です.
"""
はい、scipy.stats.norm.cdfを使用して、標準正規分布の累積分布関数を描画することができます。
以下のコードを使用して、標準正規分布の確率密度関数と累積分布関数を描画できます。
"""
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
x = np.linspace(-3, 3, 100)
plt.plot(x, norm.pdf(x))
plt.plot(x, norm.cdf(x))
plt.show()
次に、オレンジ色のプロットを累積分布関数F(x)ではなく、1-F(x)にしてみたいと思います.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
x = np.linspace(-3, 3, 100)
plt.plot(x, norm.pdf(x))
plt.plot(x, 1 - norm.cdf(x))
plt.show()
良い感じになりました.
z=0のときp値は1になるみたいなので2倍してみましょう.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
x = np.linspace(-3, 3, 100)
plt.plot(x, norm.pdf(x))
plt.plot(x, 2 * (1 - norm.cdf(x)))
plt.show()
zを絶対値に変換してから、オレンジを計算してみましょう.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
x = np.linspace(-3, 3, 100)
plt.plot(x, norm.pdf(x))
plt.plot(x, 2 * (1 - norm.cdf(abs(x))))
plt.show()
間違っているかもしれませんが上図で完成にしたいと思います.
横軸に確率変数の実現値をとると標準正規分布のp値はこのようなプロットになるのですね.
一般的に何と表現するのでしょうか. 私にはメガネフレームのように見えます.
縦軸と横軸を逆転させてみましょう.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
x = np.linspace(-3, 3, 100)
p_for_linspace_x = 2 * (1 - norm.cdf(abs(x)))
pdf_for_linspace_x = norm.pdf(x)
plt.plot(pdf_for_linspace_x, x)
plt.plot(p_for_linspace_x, x)
plt.show()
p=1以外では、同じpに対応するzがペアで存在するということでしょうか.
このあと色々考えてみたのですが行き詰まってしまいました.
心機一転、linspaceではなく、無作為抽出したいと思います.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
import pandas as pd
fig, (ax1, ax2) = plt.subplots(1, 2)
z = np.random.normal(loc=0, scale=1, size=100000)
ax1.hist(z)
ax1.set_xlabel('z from Z')
ax1.set_ylabel('Freqency of 100,000 times')
p = 2 * (1 - norm.cdf(abs(z)))
ax2.hist(p)
ax2.set_xlabel('P-Value')
ax2.set_ylabel('Freqency of 100,000 tiems')
plt.show()
確かにp値は一様に分布しているように見えます.
z=0付近の値の出現頻度が高いので、p=1付近とp=0付近ではp値の出現率が異なりそうなものですが、初学者にありがちな勘違いということでしょうか.
正規分布は左右に無限に広がっているので、左側のグラフの横軸を同じぐらいのp値となる区間に中央から左右2方向に向かって区切っていくと、左右に広がるにつれて1区間の長さが大きくなっていくようなイメージを持ちました.
まとめ
- 標準正規分布に従う確率変数を何度も無作為抽出してp値の出現頻度を確認したところ、確かに一様な分布になる.
- 心情としては納得できない.