相関係数のシミュレーション
母集団の相関係数に対して、サンプルサイズが十分に大きくないとき、サンプルの相関係数はかなり変動する。この様子を見るために、シミュレーションプログラムを作成してみた。リスト１にプログラムを示す。プログラムファイルは、Windows用は圧縮ファイルsimu_r_file.zipとしても用意した。Linux用はsimu_r_cos.zipである。それぞれダウンロード展開すれば自由に利用できる。
PyGameは、コマンド
pip install
pygame
の実行により、簡単にインストールできた。
リスト１のプログラムを実行すると、母集団の相関係数とサンプルサイズの設定が求められる。
rho = 0.7
n = 50
これらの値を設定すると、シミュレーションが始まる。
以下の説明は、Windowsの場合である。UbuntuやCentOS Streamなどの場合は、このサイトの後半で説明している。
リスト1では、表示ウィンドウを前面に出すためにモジュールwin32guiを用いているが、これはpywin32のインストールによりimportできるようになる。pywin32のインストールは、まず次のコマンドを実行する。
pip install
pywin32
続いて、次のコマンドをPythonのインストールされているフォルダに移って実行する。
C:\Pythonのインストールされているフォルダ\python.exe
Scripts\pywin32_postinstall.py -install
Pythonのインストールされているフォルダに移るのは、フォルダScriptsがPythonのインストールされているフォルダにあるからである。
これで、再起動後、次のようにwin32guiがimportできる。
import
win32gui
ただし、Anacondaの場合は、上のpipによるインストールで完了である。\python.exe Scripts\pywin32_postinstall.py -installによるインストール作業は不要である。
rhoの値とnの値を設定すると、図１のウィンドウが表示される。
図１
ウィンドウ内を２回クリックすると、シミュレーションが始まる。この図１の画面の操作は、シミュレーションの表示画面のプログラムがキーボードの入力を受け取れるようにするためである。
シミュレーションの結果は、PyGameにより図２のように表示される。
図２
上段に、母集団における相関係数とサンプルサイズが表示され、その下に、サンプルの相関係数、ｐ値、95％信頼区間が表示されている。
最下段に示されているように、次のサンプリングを行うときはNキーを押す。Nキーを押してサンプリングを繰り返した例を以下に示す。
終了するときはEキーを押せばよい。
なお、Pythonによるデータ分析の入門書＜岡本安晴「いまさら聞けないPythonでデータ分析」丸善出版＞を上梓した。
リスト１ 相関係数のシミュレーションプログラム（simu_r.py）
import pygame as pg
import scipy.stats as ss
import numpy as np
from win32gui import SetWindowPos, SetFocus
import win32con
rho = float(input('rho = '))
n = int(input('n = '))
pg.init()
wt = int(3000 / n)
pg.display.set_caption('Simulation of
Sample Correlation')
screen =
pg.display.set_mode(size=(700, 700))
SetWindowPos(pg.display.get_wm_info()['window'],
-1, 0, 0, 0, 0,
win32con.SWP_ASYNCWINDOWPOS | win32con.SWP_NOMOVE |
win32con.SWP_NOSIZE | win32con.SWP_SHOWWINDOW)
screen.fill((255,255,255))
pg.display.flip()
#
#
母数の表示用文字列
#
my_font = pg.font.SysFont(None, 32)
img_str1 = my_font.render(f'rho =
{rho}, n = {n}', True,
(0,0,0), (255,255,255))
my_font_L = pg.font.SysFont(None, 40)
img_str_L = my_font_L.render('Next
-> key N, End ->
key E', True,
(0,0,0), (255,255,255))
rect_str_L = img_str_L.get_rect()
my_font_L1 = pg.font.SysFont(None, 60)
img_str_mouse =
my_font_L1.render('Click here, twice', True,
(0,0,0), (255,255,255))
rect_str_mouse =
img_str_mouse.get_rect()
screen.blit(img_str_mouse, (350 -
rect_str_mouse.width//2,
350 - rect_str_mouse.height//2))
pg.display.flip()
ck1 = True
while ck1:
for e in pg.event.get():
if
e.type == pg.MOUSEBUTTONDOWN:
ck1 = False
ck = True
while ck:
screen.fill((255,255,255))
#
#
データのサンプリング
#
data =
ss.multivariate_normal.rvs(mean=[0, 0], cov=[[1, rho],[rho, 1]],
size=n)
#
#
相関係数の計算
#
results =
ss.pearsonr(data.T[0], data.T[1])
#
#
信頼区間の計算
#
vL, vU =
results.confidence_interval(0.95)
v_max = np.max(data)
v_min = np.min(data)
#
#
点の描画用スケーリング関数
#
def scale_v(v):
return
100 + (700 - 200) * (v - v_min) / (v_max - v_min)
#
#
散布図
#
for d in data:
pg.draw.circle(screen, (0,0,250), (scale_v(d[0]), 710 - scale_v(d[1])),
7, 3)
pg.display.flip()
pg.time.wait(wt)
# 各点の描画時間間隔
#
#
母数の表示
#
screen.blit(img_str1, (50,
10))
#
#
相関係数、ｐ値、区間推定の表示
#
img_str2 =
my_font.render(f'r = {results.statistic:.3f}, p = {results.pvalue:.5f}' +
f', 95%CI = [{vL:.3f},
{vU:.3f}]',
True, (0,0,0), (255,255,255))
screen.blit(img_str2, (50,
40))
#
#
キー操作のメッセージ表示
#
screen.blit(img_str_L, (350
- (rect_str_L.width//2),
700 - (rect_str_L.height*2)) )
pg.display.flip()
ck1 = True
while ck1:
for
e in pg.event.get(): # キーイベントのチェック
if e.type == pg.KEYDOWN:
if e.key == pg.K_n:
ck1 = False
elif e.key == pg.K_e:
ck = False
ck1 = False
screen.fill((255,255,255)) # ホワイトアウト
pg.display.flip()
pg.quit()
サンプルスクリプトをリスト２に示す。スクリプトファイルは圧縮ファイルsimu_r_cos.zipとして用意した。
以下の実行例は、CentOS Stream 9で行った。実行環境は、Python 3.12.3/Anacondaである。ライブラリPyGameは、コマンド
pip install
pygame
によりインストールした。
リスト２のスクリプトを実行すると、母集団の相関係数とサンプルサイズの設定が求められる。
rho = 0.5
n = 20
相関係数とサンプルサイズを設定するとシミュレーションが始まり、下図のように散布図が表示される。
上部１行目に母集団パラメータ、２行目に統計量が表示されている。
下部に表示されているように、Nキーを押すと次のシミュレーションが始まり、Eキーを押すとプログラムの終了である。
上の実行は、CentOS Stream 9で行ったが、Ubuntu 24/WSL/Windows 10でもPython 3.12.3/Anacondaの仮想環境で実行できた。
リスト２ Linux（Ubuntu, CentOS Stream 9など）の場合
import pygame as pg
import scipy.stats as ss
import numpy as np
rho = float(input('rho = '))
n = int(input('n = '))
pg.init()
wt = int(3000 / n)
pg.display.set_caption('Simulation of
Sample Correlation')
screen =
pg.display.set_mode(size=(700, 700))
screen.fill((255,255,255))
pg.display.flip()
#
#
母数の表示用文字列
#
my_font = pg.font.SysFont(None, 32)
img_str1 = my_font.render(f'rho =
{rho}, n = {n}', True,
(0,0,0), (255,255,255))
my_font_L = pg.font.SysFont(None, 40)
img_str_L = my_font_L.render('Next
-> key N, End ->
key E', True,
(0,0,0),
(255,255,255))
rect_str_L = img_str_L.get_rect()
ck = True
while ck:
screen.fill((255,255,255))
#
#
データのサンプリング
#
data =
ss.multivariate_normal.rvs(mean=[0, 0],
cov=[[1, rho],[rho, 1]], size=n)
#
#
相関係数の計算
#
results =
ss.pearsonr(data.T[0], data.T[1])
#
#
信頼区間の計算
#
vL, vU =
results.confidence_interval(0.95)
v_max = np.max(data)
v_min = np.min(data)
#
#
点の描画用スケーリング関数
#
def scale_v(v):
return 100 + (700 - 200) * (v - v_min) / (v_max - v_min)
#
#
散布図
#
for d in data:
pg.draw.circle(screen, (0,0,250), (scale_v(d[0]), 710 - scale_v(d[1])),
7, 3)
pg.display.flip()
pg.time.wait(wt)
# 各点の描画時間間隔
#
#
母数の表示
#
screen.blit(img_str1, (50,
10))
#
#
相関係数、ｐ値、区間推定の表示
#
img_str2 =
my_font.render(f'r = {results.statistic:.3f}, p = {results.pvalue:.5f}' +
f', 95%CI = [{vL:.3f},
{vU:.3f}]',
True, (0,0,0), (255,255,255))
screen.blit(img_str2, (50,
40))
#
#
キー操作のメッセージ表示
#
screen.blit(img_str_L, (350
- (rect_str_L.width//2),
700 - (rect_str_L.height*2)) )
pg.display.flip()
ck1 = True
while ck1:
for
e in pg.event.get(): # キーイベントのチェック
if e.type == pg.KEYDOWN:
if e.key == pg.K_n:
ck1 = False
elif e.key == pg.K_e:
ck = False
ck1 = False
screen.fill((255,255,255)) # ホワイトアウト
pg.display.flip()
pg.quit()