単回帰モデル分析
＜相関係数の検定＞にあげてあるプログラムに回帰分析の機能を付け加えた。回帰分析は次の単回帰モデル
を用いているが、詳しくは（岡本安晴「データ分析のための統計学入門」、おうふう）を参照されたい。また、回帰分析一般については（岡本安晴「心理学データ分析と測定」勁草書房；第１章〜第４章）で説明している。
プログラムは、相関係数の検定のプログラムに単回帰分析のためのコードを追加したもので、以下のようになっている。
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ simpleRegMain.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
import matplotlib.pyplot as plt
from simpleRegSub import Ccalc
pairs = [[70, 168], [68, 165], [62,
159], [56, 156], [70, 170],
[64, 162], [59, 162], [70, 169], [68, 168], [69, 170],
[67, 168], [63, 166], [58, 163], [74, 170], [65, 168],
[60, 160], [65, 165], [72, 175], [61, 164], [77, 174],
[68, 164], [64, 161], [70, 166], [47, 147], [66, 165],
[68, 170], [73, 174], [74, 172], [63, 164], [83, 182],
[58, 159], [64, 166], [67, 167], [75, 174], [67, 163],
[68, 166], [68, 167], [66, 165], [63, 165], [55, 154],
[61, 163], [60, 160], [62, 164], [62, 160], [59, 160],
[51, 155], [59, 163], [68, 170], [56, 157], [79, 180]]
calc_p = Ccalc(pairs)
t, df, p, r, a, b = calc_p.calc()
print("\nt = %f df = %d"%(t, df))
print("Prob(x < %f) = %f Prob(x
> %f) = %f df
= %d"%(t, p, t, 1.0 - p, df))
print("\nY = %f * X +
%f"%(a, b))
min_x = max_x = pairs[0][0]
min_y = max_y = pairs[0][1]
for v in pairs:
if min_x > v[0]:
min_x = v[0]
if max_x < v[0]:
max_x = v[0]
if min_y > v[1]:
min_y = v[1]
if max_y < v[1]:
max_y = v[1]
reg_min = a * min_x + b
reg_max = a * max_x + b
xdata = []
ydata = []
for v in pairs:
xdata.append(v[0])
ydata.append(v[1])
plt.xlabel("X")
plt.ylabel("Y")
plt.title("Y = %f * X + %f r = %.3f"%(a, b,
r))
plt.plot(xdata, ydata, 'go')
plt.plot([min_x, max_x], [reg_min,
reg_max])
plt.show()
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
データ対
をリスト
として表し、これらのリストを要素とするリストpairsとして入力データを表している。データを表すリストpairsを実引数としてクラス型Ccalcのオブジェクトcalc_pを作成して、メソッドcalc_p.calc()によりｔ値などを算出している。クラス型CcalcはプログラムsimpleRegSub.pyにおいて以下のように宣言されている。
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ simpleRegSub.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
import math
from calc_p_from_t
import p_value
def sqr(x):
return x * x
class Ccalc:
def __init__(self, pairs):
self.pairs = pairs
def calc(self):
n =
len(self.pairs)
print("n = %d"%n)
sum_x = 0.0
sum_y = 0.0
for
i in range(0, n):
sum_x += self.pairs[i][0]
sum_y += self.pairs[i][1]
mean_x = sum_x / n
mean_y = sum_y / n
print("Mean_x = %f"%mean_x)
print("Mean_y = %f"%mean_y)
ssum_x = 0.0
ssum_y = 0.0
sum_xy = 0.0
for
i in range(0, n):
ssum_x += sqr(self.pairs[i][0] - mean_x)
ssum_y += sqr(self.pairs[i][1] - mean_y)
sum_xy += (self.pairs[i][0] - mean_x) * (self.pairs[i][1] - mean_y)
var_x = ssum_x / n
var_y = ssum_y / n
print("Var_x = %f"%var_x)
print("Var_y = %f"%var_y)
sd_x
= math.sqrt(var_x)
sd_y
= math.sqrt(var_y)
print("sd_x = %f"%sd_x)
print("sd_y = %f"%sd_y)
cor
= sum_xy / math.sqrt(ssum_x * ssum_y)
print("r = %f"%cor)
t =
cor * math.sqrt(n - 2.0) / math.sqrt(1.0 - (cor ** 2.0))
df =
n - 2
p =
p_value(t, df)
a =
cor * sd_y / sd_x
b =
mean_y - a * mean_x
return t, df, p, cor, a, b
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
単回帰モデルの係数は次式
によって求めている。
上のプログラムを実行すると、図１のように散布図に回帰直線が描かれたWindowsが表示される。
図１
右上のＸ印をクリックしてWindowを閉じるとプログラムの実行終了となる。実行終了後のIDLEのコンソール画面は図２のようになっている。
図２
表示されている確率は帰無仮説「相関係数＝０」すなわち「a=0」に対する確率である。
プログラムのソースコードファイルは圧縮ファイルSimpleRegPrg.zipにまとめた。ファイル名SimpleRegPrg.zipをクリックしてダウンロード展開すれば利用できる。