「入門はじめての統計解析」をPythonで写経 Vol.11 ～ ４章「はじめての統計的検定」④相関係数の検定：無相関の検定、母相関係数の検定、２つの母相関係数の差の検定
４章「はじめての統計的検定」
書籍の著者 石村貞夫 先生
この記事は、書籍「入門はじめての統計解析」４章「はじめての統計的検定」の Python写経活動 を取り扱います。
書籍の図・表・計算を淡々とPython化する写経シリーズです。
この記事は４章の統計的検定テーマのうち相関係数に関する３つの検定を実践します。
・無相関の検定
・母相関係数の検定
・２つの母相関係数の差の検定
ChatGPTの活用も継続してまいります！
では書籍を開いて統計解析の旅に出発です🚀
はじめに
このブログシリーズは、書籍「入門はじめての統計解析」（東京図書、「テキスト」と呼びます）の Python 写経を通じて得た「統計解析の楽しさ」をご紹介します。
書籍の紹介と引用表記はリンク先の記事に掲載しています。
４章 はじめての統計的検定
この記事は４章の以下のSectionを取り扱います。
4.9 相関係数の検定
記事に用いるデータは、テキストに掲載されたデータそのものを引用しています。
データ件数の少ないものはコード上でデータを登録し、データ件数の多いものはCSVファイル化してデータを読み込みしています。
４章で用いるライブラリをインポートします。
### インポート
# 数値計算
import math # python標準ライブラリ
import numpy as np
import pandas as pd
# 統計
import scipy.stats as stats
import pingouin as pg
from statsmodels.stats.proportion import proportions_ztest # 母比率のz検定
from statsmodels.stats.weightstats import ttest_ind # 2標本のt検定
# 描画
import matplotlib.pyplot as plt
import seaborn as sns
plt.rcParams['font.family'] = 'Meiryo'イントロダクション
前回に引き続いて「統計的検定」を頑張って取り組んでいきましょう。
ChatGPTによる「やさしいイントロ」をご堪能ください。
「この2つの数字、なんとなく一緒に動いてる気がする」
「売上と広告費って、やっぱり関係あるんじゃない？」
そんな“気がする”を、データで確かめてくれるのが相関係数です。
相関係数は、2つのデータの間にどれくらい“関係”があるかを、$${−1}$$ から $${1}$$ の数値で表してくれる指標。
でも実際のデータから計算された相関係数が、「たまたまそうなっただけ」なのか「本当に関係がある」のか――
それを見極めるために使うのが、今回のテーマである相関係数に関する検定です。
🍀 今回取り上げる３つの検定は…
無相関の検定：相関係数がゼロ（＝関係なし）と言えるかどうかを調べる
母相関係数の検定：特定の理論的な相関値と差があるかを調べる
２つの母相関係数の差の検定：グループAとBで、相関の強さに差があるかを比べる
数字と数字の“つながり”に目を向けることで、見えてくるものがあります。
今回は、「関係があると言えるのか？」を見きわめる検定の世界に、やさしく足を踏み入れてみましょう。
Section 4.9 相関係数の検定
■ 相関係数の振り返り
相関係数は２つのデータの関係の強さを表す統計量です。
相関係数を直感的にイメージできる「散布図」で確認しましょう。
## 設定と準備
# 標本サイズ
N = 10
# 乱数生成器の初期化
rng = np.random.default_rng(seed=1)
## データの準備
# x の作成：一様分布乱数
x = rng.uniform(size=N)
# y の作成：y = x + ε, ε ~ Normal(0, 1)
y = x + rng.normal(size=N)
# xとyの相関係数の算出
xy_corrcoef = np.corrcoef(x, y)[0, 1]
## 可視化
# 描画領域の設定
fig, ax = plt.subplots()
# x,yの散布図の描画
ax.plot(x, y, 'o', ms=8, alpha=0.8)
# 修飾
ax.set(xlabel='$x$', ylabel='$y$', title=f'相関係数 = {xy_corrcoef:.3f}')
ax.grid(lw=0.5);【実行結果】
データ $${x,y}$$ の間には右上がりの傾向が見えます。
相関係数 $${0.787}$$ は正の強めの相関を示しています。
無相関の検定
■ 無相関の検定：検定統計量の公式 p.190
無相関の検定は、２つのデータの間に相関関係が無い（無相関）かどうかの検定です。
❓️母集団分布は…❓️
無相関の検定における母集団分布の仮定について、テキストは明示していないようです。
ChatGPTの回答は「２つの母集団分布に２変量正規分布を仮定する」です。
理論的な裏付けを取れていないので…
ChatGPTとのやりとりを後段の「無相関の検定における母集団分布の仮定：ChatGPT問答」に掲載いたします。
【仮説】
仮説を次のように設定します。
帰無仮説 $${H_0}$$：２つの変数 $${x}$$ と $${y}$$ は無相関である
対立仮説 $${H_1}$$：２つの変数 $${x}$$ と $${y}$$ は相関がある
片側検定は無いそうです。
◆ ◆ ◆
無相関の検定における検定統計量の公式をテキストよりお借りします。
標本サイズ $${N}$$ の標本相関係数が $${r}$$ のとき、無相関の検定統計量
$$
T(r, N) = \cfrac{r\sqrt{N-2}}{\sqrt{1 - r^2}}
$$
は自由度 $${N-2}$$ の $${t}$$ 分布に従います。
公式に則って「無相関の検定関数」を定義します。
確率計算には scipy.stats を利用します。
### 無相関の検定 p.190
# 標本相関係数算出関数
def calc_corr(x_list1, x_list2):
# 標本サイズ ※２つの標本サイズは同じ
N = len(x_list1)
# ２つの標本平均
x_bar, y_bar = sum(x_list1) / N, sum(x_list2) / N
# 分子、分母１，分母２の計算
numerator = sum((x - x_bar) * (y - y_bar) for x, y in zip(x_list1, x_list2))
denominator1 = (sum([(x - x_bar)**2 for x in x_list1]))**(1/2)
denominator2 = (sum([(y - y_bar)**2 for y in x_list2]))**(1/2)
# 標本相関係数r
r = numerator / (denominator1 * denominator2)
# 戻り値: 標本相関係数r, 標本サイズN
return r, N
# 無相関のt検定統計量の算出関数
def t_stat_non_corr(x_list1, x_list2):
# 標本相関係数の算出
r, N = calc_corr(x_list1, x_list2)
# 戻り値: t検定統計量
return r * (N - 2)**(1/2) / (1 - r**2)**(1/2), r
# 無相関のt検定関数 ※確率計算はscipy.stats利用
def non_corr_ttest(x_list1, x_list2, alpha=0.05, alternative='two-sided'):
# 標本サイズ ※２つの標本サイズは同じ
N = len(x_list1)
# t検定統計量と標本相関係数の算出
t_val, r = t_stat_non_corr(x_list1, x_list2)
# 自由度(N-2)のt分布の設定
t_dist = stats.t(df=N - 2)
# 棄却限界値c_valとp値p_valの算出
match alternative:
case 'two-sided':
c_val = t_dist.ppf(q=1 - alpha/2)
c_val = -c_val, c_val
p_val = t_dist.sf(x=abs(t_val)) * 2
case 'less':
c_val = t_dist.ppf(q=alpha)
p_val = t_dist.cdf(x=t_val)
case 'greater':
c_val = t_dist.ppf(q=1 - alpha)
p_val = t_dist.sf(x=t_val)
return {'t_value': t_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val,
'r': r}次のテストデータで無相関の検定をやってみましょう。
$$
\begin{array}{cc}
変数1 & 変数2 \\
\hline
1 & 1.5 \\
2 & 1.8 \\
3 & 2.0 \\
\end{array}
$$
検定の条件です。
帰無仮説 $${H_0}$$ 「２つのデータは無相関である」
有意水準 5%（$${\alpha=0.05}$$）
両側検定
# テスト
x_list1, x_list2 = [1, 2, 3], [1.5, 1.8, 2.0]
non_corr_ttest(x_list1, x_list2)【実行結果】
「有意水準 5% で有意とは言えず、帰無仮説を棄却できない」です。
相関係数は $${0.993}$$ ですが、検定の結果、「２つのデータは無相関」を受容します。
出力内容を簡単に説明します。
t_value：検定統計量 $${T}$$ の実現値
c_value：棄却限界値
帰無仮説を棄却できる検定統計量 $${T}$$ のしきい値
$${T}$$ の実現値が棄却限界値の外のとき、帰無仮説を棄却できるalpha：有意水準 $${\alpha}$$
p_value：$${p}$$ 値
帰無仮説のもとで検定統計量が実現値になる確率
$${p}$$ 値が有意水準以下（または未満）のとき、帰無仮説を棄却できるr：標本相関係数
scipy.stats で検算します。
### scipy.statsで答え合わせ (標本相関係数, 無相関の検定のp値)
stats.pearsonr(x_list1, x_list2)【実行結果】
（ピアソンの積率）相関係数、$${p}$$ 値です。
検算結果は一致しました。
◆ ◆ ◆
■ 無相関の検定（両側検定）例題 p.191
p.191 例）の大気汚染と水質汚濁の間の相関の有無について、無相関の検定で調べます。
仮説は次のとおりです。
・帰無仮説 $${H_0}$$ 「２つのデータは無相関である」
・対立仮説 $${H_1}$$ 「２つのデータは相関がある」
⇒両側検定：相関の有無に興味がある
① 自作関数で検定実行
まず２つのデータを設定します。
### 無相関の検定 p.191
# データ
data2 = pd.DataFrame({'市名': ['A', 'B', 'C', 'D', 'E', 'F', 'G', 'H', 'O'],
'大気汚染': [113, 64, 16, 45, 28, 19, 30, 82, 76],
'水質汚濁': [31, 5, 2, 17, 18, 2, 9, 25, 13]},
index=range(1, 10))
data2.index.name = 'No.'
data2【実行結果】
大気汚染と水質汚濁を散布図で可視化しましょう。
# 散布図の描画 p.191
sns.scatterplot(data=data2, x='大気汚染', y='水質汚濁', s=80, alpha=0.7)
plt.xlim(0, 150)
plt.ylim(0, 40);【実行結果】
右上がりの傾向が見られます。
正の相関がありそうです。
相関係数を確認します。
# pandasで相関係数を算出
data2.corr(numeric_only=True)【実行結果】
相関係数は $${0.761}$$ です。
２つのデータの間に強めの正の相関がある感じがいたします。
では有意水準 $${5\%}$$、両側検定で２つの母分散の差の検定を実行します。
# 関数利用
result = non_corr_ttest(data2['大気汚染'], data2['水質汚濁'])
result【実行結果】
検定統計量 $${3.104}$$ は、棄却限界値（上側の閾値） $${2.365}$$ より大きい（棄却域に入る）ので、有意水準 $${5\%}$$ で帰無仮説は棄却されます。
つまり、大気汚染と水質汚濁の間に 相関があると言える、です。
② 検定の可視化
$${t}$$ 分布を可視化して棄却される様子を見てみましょう。
### 図示 p.193
## 設定と準備
# 統計関連の値
N = data2.shape[0] # 標本サイズ
alpha = result['alpha'] # 有意水準
lower, upper = result['c_value'] # 棄却限界値
t_val = result['t_value'] # 検定統計量
t_dist = stats.t(df=N - 2) # t分布
# グラフ描画用設定
color = 'tab:blue' # 基本の色
x_min, x_max = -5, 5 # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001) # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101) # 下側棄却域用のx軸の値
x_val_upper = np.linspace(upper, x_max, 101) # 上側棄却域用のx軸の値
## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# t分布の確率密度関数の描画
plt.plot(x_val, t_dist.pdf(x_val), color=color,
label=f'自由度({N}-2)の$t$分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, t_dist.pdf(x_val_lower),
color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 上側棄却域の塗りつぶし
plt.fill_between(x_val_upper, 0, t_dist.pdf(x_val_upper),
color=color, alpha=0.3)
# 検定統計量の垂直線の描画
plt.axvline(t_val, color='tab:red', ls='--', label=f'検定統計量={t_val:.3f}')
# 修飾
plt.xticks([lower, 0, upper])
plt.legend();【実行結果】
ベル型の曲線は自由度 $${N -2=7}$$ の $${t}$$ 分布の確率密度関数です。
両端の青い領域が有意水準 $${5\%}$$ の棄却域です。
両側検定なので両端に棄却域が $${2.5\%}$$ 分づつあります。
赤い点線が検定統計量 $${T}$$ です。
$${T}$$ が両方の棄却域に入っているので、帰無仮説は棄却されます。
③ Python ライブラリで検定実行
ここでは、scipy.stats、pingouin を利用します。
いずれも引数に２つのサンプルデータを渡しています。
🖲️scipy.stats
# scipy.stats利用 (標本相関係数, 無相関の検定のp値)
stats.pearsonr(data2['大気汚染'], data2['水質汚濁'])【実行結果】
相関係数と $${p}$$ 値を得ました。
🖲️pingouin
# pingouin利用
pg.corr(data2['大気汚染'], data2['水質汚濁'])【実行結果】
主な出力は、r：相関係数、CI95%：相関係数の95%信頼区間、p_val：$${p}$$ 値です。
■ （参考）無相関の検定における母集団分布の仮定：ChatGPT問答
冒頭の母集団分布の仮定に関するChatGPTとのやりとりを掲載いたします。
内容の適否は確認しておりませんのでご注意ください。
母相関係数の検定
母相関係数の検定は、データの相関係数が「ある相関係数」であるかどうかの検定です。
「ある相関係数」は帰無仮説の相関係数 $${\rho_0}$$ のことです。
母相関係数の検定は何らかの母集団分布を仮定しているのでしょうか？
ちなみにChatGPTの回答は「２変量正規分布を仮定する」です。
■ 仮説
覆したい母相関係数の仮説 $${\rho_0}$$ を次のように設定します。
帰無仮説 $${H_0}$$：$${\rho = \rho_0}$$
支持したい仮説＝対立仮説は、次の３パターンのいずれかを設定します。
① 両側検定：母相関係数が $${\rho_0}$$ と異なることを検定する場合
・対立仮説 $${H_1}$$：$${\rho \neq \rho_0}$$
② 片側検定：母相関係数が $${\rho_0}$$ より小さいことを検定する場合
・対立仮説 $${H_1}$$：$${\rho < \rho_0}$$
③ 片側検定：母相関係数が $${\rho_0}$$ より大きいことを検定する場合
・対立仮説 $${H_1}$$：$${\rho > \rho_0}$$
◆ ◆ ◆
■ 母相関係数の検定：検定統計量の公式 p.194
母相関係数の検定における検定統計量の公式をテキストよりお借りします。
標本サイズ $${N}$$ の標本相関係数を $${r}$$ とすると、母相関係数の検定統計量
$$
T(r) = \sqrt{N - 3} \left(\cfrac{1}{2}\log\cfrac{1+r}{1-r} - \cfrac{1}{2}\log\cfrac{1+\rho_0}{1-\rho_0}\right)
$$
は標準正規分布で近似されます。
$${\rho_0}$$ は帰無仮説の相関係数です。
公式に則って「母相関係数の検定関数」を定義します。
確率計算には scipy.stats を利用します。
### 母相関係数の検定関数 p.194
# 標本相関係数算出関数
def calc_corr(x_list1, x_list2):
# 標本サイズ ※２つの標本サイズは同じ
N = len(x_list1)
# ２つの標本平均
x_bar, y_bar = sum(x_list1) / N, sum(x_list2) / N
# 分子、分母１，分母２の計算
numerator = sum((x - x_bar) * (y - y_bar) for x, y in zip(x_list1, x_list2))
denominator1 = (sum([(x - x_bar)**2 for x in x_list1]))**(1/2)
denominator2 = (sum([(y - y_bar)**2 for y in x_list2]))**(1/2)
# 標本相関係数r
r = numerator / (denominator1 * denominator2)
# 戻り値: 標本相関係数r, 標本サイズN
return r, N
# 母相関係数のz検定統計量の算出関数
def z_stat_pop_corr(x_list1, x_list2, rho0):
# 標本相関係数の算出
r, N = calc_corr(x_list1, x_list2)
# 戻り値: z検定統計量
return ((N-3)**(1/2)
* (1/2 * math.log((1+r)/(1-r)) - 1/2 * math.log((1+rho0)/(1-rho0))),
r)
# 母相関係数のz検定関数 ※確率計算はscipy.stats利用
def pop_corr_ztest(x_list1, x_list2, rho0, alpha=0.05, alternative='two-sided'):
# z検定統計量と標本相関係数rの算出
z_val, r = z_stat_pop_corr(x_list1, x_list2, rho0)
std_norm_dist = stats.norm(loc=0, scale=1)
match alternative:
case 'two-sided':
c_val = std_norm_dist.ppf(q=1 - alpha/2)
c_val = -c_val, c_val
p_val = std_norm_dist.sf(x=abs(z_val)) * 2
case 'less':
c_val = std_norm_dist.ppf(q=alpha)
p_val = std_norm_dist.cdf(x=z_val)
case 'greater':
c_val = std_norm_dist.ppf(q=1 - alpha)
p_val = std_norm_dist.sf(x=z_val)
return {'z_value': z_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val,
'r': r, 'rho0': rho0}次のテストデータで母相関係数の検定をやってみましょう。
$$
\begin{array}{cc}
変数1 & 変数2 \\
\hline
1 & 0.8 \\
2 & 0.9 \\
3 & 1.1 \\
4 & 1.5
\end{array}
$$
検定の条件です。
帰無仮説 $${H_0}$$ 「母相関係数$${\rho=0.9}$$」
有意水準 5%（$${\alpha=0.05}$$）
両側検定
# テスト
x_list1, x_list2 = [1, 2, 3, 4], [0.8, 0.9, 1.1, 1.5]
pop_corr_ztest(x_list1, x_list2, rho0=0.9)【実行結果】
「有意水準 5% で有意とは言えず、帰無仮説を棄却できない」です。
検定の結果、「母相関係数は $${0.9}$$」を受容します。
データの標本相関係数は $${0.959}$$ でした。
出力内容を簡単に説明します。
z_value：検定統計量 $${T}$$ の実現値
c_value：棄却限界値
帰無仮説を棄却できる検定統計量 $${T}$$ のしきい値
$${T}$$ の実現値が棄却限界値の外のとき、帰無仮説を棄却できるalpha：有意水準 $${\alpha}$$
p_value：$${p}$$ 値
帰無仮説のもとで検定統計量が実現値になる確率
$${p}$$ 値が有意水準以下（または未満）のとき、帰無仮説を棄却できるr：標本相関係数
rho0：帰無仮説の母相関係数
この流れで例題に取り組みたいところですが…
「理解度チェック」までお待ちください！
◆ ◆ ◆
💡 Tips ：フィッシャーの $${z}$$ 変換
母相関係数の検定の公式に含まれる
$$
\cfrac{1}{2}\log\cfrac{1+r}{1-r}
$$
は「フィッシャーの $${z}$$ 変換」という操作だそうです。
ChatGPTにフィッシャーの $${z}$$ 変換を聞きました。
✅ フィッシャーのz変換とは？
標本相関係数 $${r}$$ は、元々は $${[−1,1]}$$ の範囲にあり、分布が偏っています。
そこで、次の変換を行うことで、ほぼ正規分布に近づけます：
$$
z = \frac{1}{2} \log \left(\frac{1 + r}{1 - r}\right)
$$
この変換後の量 $${z}$$ は母相関係数 $${\rho}$$ のもとでほぼ正規分布 $${\mathcal{N}(\zeta, \frac{1}{N - 3})}$$ に従うとされます。
ここで：
$$
\zeta = \frac{1}{2} \log \left(\frac{1 + \rho}{1 - \rho}\right)
$$
フィッシャーの $${\boldsymbol{z}}$$ 変換で検定統計量を標準正規分布に近似させていたのですね！
実は以前、別の書籍の写経記事で「相関係数をフィッシャーの $${z}$$ 変換にかけて正規分布近似させるシミュレーション」を実践しました。
よかったら寄り道してくださいね。
理解度チェック 母相関係数の検定 p.195
国別のカロリー摂取量とある病気の患者数の間の母相関係数が $${\rho_0=-0.3}$$ かどうかを、母相関係数の検定で確かめます。
有意水準 $${5\%}$$、両側検定で、コードで淡々と書きます。
標本データを設定します。
### 母相関係数の検定 カロリー摂取とある病気の患者数 p.195
# データ
data3 = pd.DataFrame(
{'カロリー': [2750, 2956, 2675, 3198, 1816, 2233, 2375, 2288, 1932, 2036,
2183, 2882],
'患者数': [249, 713, 1136, 575, 5654, 2107, 915, 4193, 7225, 3730, 472, 291]},
index=range(1, 13))
data3【実行結果】
散布図でデータの関係を確認しましょう。
# 散布図の描画
sns.scatterplot(data=data3, x='カロリー', y='患者数', s=80, alpha=0.7);【実行結果】
カロリー2500くらいまでは右下がりの傾向が見られます。
標本相関係数を計算しましょう。
# 標本相関係数
data3.corr()【実行結果】
相関係数 $${-0.76655}$$ は強めの負の相関です。
帰無仮説の $${\rho_0=-0.3}$$ と乖離している感じがします。
では自作関数で検定を実行します。
# 母相関係数の検定 関数利用 ★テキスト解答のz検定統計量は-0.915となっている
result = pop_corr_ztest(data3['カロリー'], data3['患者数'], rho0=-0.3)
result【実行結果】
帰無仮説は棄却されます。
つまり「母相関係数は $${-0.3}$$ とは言えない」です。
標本相関係数は $${-0.767}$$ でした。
標準正規分布を可視化して棄却できる様子を見てみましょう。
### 図示
## 設定と準備
# 統計関連の値
alpha = result['alpha'] # 有意水準
lower, upper = result['c_value'] # 棄却限界値
z_val = result['z_value'] # 検定統計量
std_norm_dist = stats.norm(loc=0, scale=1) # 標準正規分布
# グラフ描画用設定
color = 'tab:blue' # 基本の色
x_min, x_max = -3, 3 # x軸の最小値、最大値
x_val = np.linspace(x_min, x_max, 1001) # 確率密度関数用のx軸の値
x_val_lower = np.linspace(x_min, lower, 101) # 下側棄却域用のx軸の値
x_val_upper = np.linspace(upper, x_max, 101) # 下側棄却域用のx軸の値
## 描画処理
# 描画領域の設定
plt.figure(figsize=(7, 3))
# 標準正規分布の確率密度関数の描画
plt.plot(x_val, std_norm_dist.pdf(x_val), color=color, label=f'標準正規分布')
# 下側棄却域の塗りつぶし
plt.fill_between(x_val_lower, 0, std_norm_dist.pdf(x_val_lower),
color=color, alpha=0.3, label=f'有意水準{alpha}の棄却域')
# 上側棄却域の塗りつぶし
plt.fill_between(x_val_upper, 0, std_norm_dist.pdf(x_val_upper),
color=color, alpha=0.3)
# 検定統計量の垂直線の描画
plt.axvline(z_val, color='tab:red', ls='--', label=f'検定統計量={z_val:.3f}')
# 修飾
plt.xticks([lower, 0, upper])
plt.legend();【実行結果】
下側（左側）の青い棄却域に検定統計量は位置しています。
帰無仮説は棄却されるのです。
（注）この記事の検定統計量の値は、テキストの解答例と相違しています。相違する理由は不明です。
２つの母相関係数の差の検定
差の検定シリーズに戻ってきました！
２つの母相関係数の差の検定は、テキストの「アイヤ しばらく！」コーナーでひっそりと紹介されています。
おそらく「対応のない母相関係数の差の検定」です。
ぜひ Python 実装いたしましょう！
■ 仮説
覆したい母相関係数の差 $${\rho_1 - \rho_2}$$ の仮説を次のように設定します。
帰無仮説 $${H_0}$$：$${\rho_1 - \rho_2 = 0}$$
支持したい仮説＝対立仮説は、次の３パターンのいずれかを設定します。
① 両側検定：母相関係数の差が ０ と異なることを検定する場合
・対立仮説 $${H_1}$$：$${\rho_1 - \rho_2 \neq 0}$$
② 片側検定：母相関係数の差が マイナス となることを検定する場合
・対立仮説 $${H_1}$$：$${\rho_1 - \rho_2 < 0}$$
③ 片側検定：母相関係数の差が プラス となることを検定する場合
・対立仮説 $${H_1}$$：$${\rho_1 - \rho_2 > 0}$$
◆ ◆ ◆
■ ２つの母相関係数の差の検定：検定統計量の公式 p.193
２つの母相関係数の差の検定における検定統計量の公式をテキストよりお借りします。
２つの母相関係数の差の検定統計量
$$
T = \cfrac{z_1 - z_2}{\sqrt{\cfrac{1}{N_1 - 3} + \cfrac{1}{N_2 - 3}}}
$$
は標準正規分布で近似されます。
ただし：
$$
z_1 = \cfrac{1}{2}\log\cfrac{1+r_1}{1-r_1},\ z_2 = \cfrac{1}{2}\log\cfrac{1+r_2}{1-r_2}
$$
公式に則って「２つの母相関係数の差の検定関数」を定義します。
確率計算には scipy.stats を利用します。
### ２つの母相関係数の差の検定 p.193
# ２つの母相関係数の差のz検定統計量の算出関数
def z_stat_pop_corr_diff(r1, r2, N1, N2):
# z1, z2の算出
z1 = 1/2 * math.log((1 + r1)/(1 - r1))
z2 = 1/2 * math.log((1 + r2)/(1 - r2))
# 戻り値: z検定統計量
return (z1 - z2) / math.sqrt(1/(N1 - 3) + 1/(N2 - 3))
# ２つの母相関係数差のz検定関数 ※確率計算はscipy.stats利用
def pop_corr_diff_ztest(r1, r2, N1, N2, alpha=0.05, alternative='two-sided'):
# z検定統計量と標本相関係数の算出
z_val = z_stat_pop_corr_diff(r1, r2, N1, N2)
# 標準正規分布の設定
std_norm_dist = stats.norm(loc=0, scale=1)
# 棄却限界値c_valとp値p_valの算出
match alternative:
case 'two-sided':
c_val = std_norm_dist.ppf(q=1 - alpha/2)
c_val = -c_val, c_val
p_val = std_norm_dist.sf(x=abs(z_val)) * 2
case 'less':
c_val = std_norm_dist.ppf(q=alpha)
p_val = std_norm_dist.cdf(x=z_val)
case 'greater':
c_val = std_norm_dist.ppf(q=1 - alpha)
p_val = std_norm_dist.sf(x=z_val)
return {'z_value': z_val, 'c_value': c_val, 'alpha': alpha, 'p_value': p_val}次のテストデータで２つの母相関係数の差の検定をやってみましょう。
$$
\begin{array}{cc}
& 標本サイズ & 標本相関係数 \\
\hline
標本1 & N_1 = 100 & r_1 = 0.8 \\
標本2 & N_2 = 120 & r_2 = 0.6 \\
\end{array}
$$
検定の条件です。
帰無仮説 $${H_0}$$ 「母相関係数の差 $${\rho_1 = \rho_2}$$」
有意水準 5%（$${\alpha=0.05}$$）
両側検定
# テスト
r1, r2 = 0.8, 0.6
N1, N2 = 100, 120
pop_corr_diff_ztest(r1, r2, N1, N2)【実行結果】
「有意水準 5% で有意であり、帰無仮説は棄却される」です。
つまり「２つの母相関係数は $${\rho_1 = \rho_2}$$ とは言えない」です。
出力内容を簡単に説明します。
z_value：検定統計量 $${T}$$ の実現値
c_value：棄却限界値
帰無仮説を棄却できる検定統計量 $${T}$$ のしきい値
$${T}$$ の実現値が棄却限界値の外のとき、帰無仮説を棄却できるalpha：有意水準 $${\alpha}$$
p_value：$${p}$$ 値
帰無仮説のもとで検定統計量が実現値になる確率
$${p}$$ 値が有意水準以下（または未満）のとき、帰無仮説を棄却できる
記事の最後をChatGPTに締めくくってもらいましょう！
📘 ChatGPTのひとこと：
「この2つ、なんか関係ありそう？」
そんな直感を、データの力でちゃんと確かめてみたくなったら、それはもう、相関との出会いのサインです。
相関係数の検定は、“つながりがありそう”を“たしかにありそう”に変えてくれる道具。
データとデータの会話に耳をすませたあなたは、きっともう統計のことばを少しずつ話し始めています。
今回の写経は以上です。
シリーズの記事
次の記事
前の記事
目次
ブログの紹介
note で７つのシリーズ記事を書いています。
ぜひ覗いていってくださいね！
１．のんびり統計
統計検定２級の問題集を手がかりにして、確率・統計をざっくり掘り下げるブログです。
雑談感覚で大丈夫です。ぜひ覗いていってくださいね。
統計検定２級公式問題集CBT対応版に対応しています。
Python、EXCELのサンプルコードの配布もあります。
２．実験！たのしいベイズモデリング１＆２をPyMC Ver.5で
書籍「たのしいベイズモデリング」・「たのしいベイズモデリング２」の心理学研究に用いられたベイズモデルを PyMC Ver.5で描いて分析します。
この書籍をはじめ、多くのベイズモデルは Ｒ言語＋Stanで書かれています。
PyMCの可能性を探り出し、手軽にベイズモデリングを実践できるように努めます。
身近なテーマ、イメージしやすいテーマですので、ぜひぜひPyMCで動かして、一緒に楽しみましょう！
３．実験！岩波データサイエンス1のベイズモデリングをPyMC Ver.5で
書籍「実験！岩波データサイエンスvol.1」の４人のベイジアンによるベイズモデルを PyMC Ver.5で描いて分析します。
この書籍はベイズプログラミングのイロハをざっくりと学ぶことができる良書です。
楽しくPyMCモデルを動かして、ベイズと仲良しになれた気がします。
みなさんもぜひぜひPyMCで動かして、一緒に遊んで学びましょう！
４．楽しい写経 ベイズ・Python等
ベイズ、Python、その他の「書籍の写経活動」の成果をブログにします。
主にPythonへの翻訳に取り組んでいます。
写経に取り組むお仲間さんのサンプルコードになれば幸いです🍀
５．RとStanではじめる心理学のための時系列分析入門 を PythonとPyMC Ver.5 で
書籍「RとStanではじめる心理学のための時系列分析入門」の時系列分析をPythonとPyMC Ver.5 で実践します。
この書籍には時系列分析のテーマが盛りだくさん！
時系列分析の懐の深さを実感いたしました。
大好きなPythonで楽しく時系列分析を学びます。
６．データサイエンスっぽいことを綴る
統計、データ分析、AI、機械学習、Pythonのコラムを不定期に綴っています。
統計・データサイエンス書籍にまつわる記事が多いです。
「統計」「Python」「数学とPython」「R」のシリーズが生まれています。
７．Python機械学習プログラミング実践記
書籍「Python機械学習プログラミング PyTorch & scikit-learn編」を学んだときのさまざまな思いを記事にしました。
この書籍は、scikit-learnとPyTorchの教科書です。
よかったらぜひ、お試しくださいませ。
最後までお読みいただきまして、ありがとうございました。
いいなと思ったら応援しよう！
応援ありがとうございます。これからもがんばって記事を作成します！