とても有名な「統計学入門」をpythonで実装しながら読んでいきます。
本をしっかり読むことと実装の勉強が目的ですので説明が不足していたりします。
また、実装方法に関してはベストなものではないと思いますのでご了承ください。
前回
1次元のデータ
2次元のデータ
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from IPython.display import display
import seaborn as sns
import japanize_matplotlib
d = pd.read_csv('pokemon_for_stats.csv').iloc[:,1:]
gen = ["第1世代"]*151+["第2世代"]*100+["第3世代"]*135+["第4世代"]*107+["第5世代"]*156+["第6世代"]*72+["第7世代"]*88+["第8世代"]*89
d['世代'] = gen
num_cols = ['重さ','高さ','HP','攻撃','防御','特攻','特防','素早さ','合計']
d
多次元データの統計学は、多くの変数を一括してその間の関係を扱うものであり、1次元データ分析を繰り返すことではない。
簡単のため$p=2$で、2変数$x,y$を考えてみる。
$x$と$y$の間に区別を設けず対等に見る方法を相関、$x$から$y$を見るとき回帰という。
散布図と分割表
観測対象$i(i=1,\cdots,n)$から、1組のデータ$(x_i,y_i)$が得られ、それらが両方とも量的データである場合、横軸に$x$、縦軸に$y$をとって、各観測対象を平面上にプロットした図を散布図と呼ぶ。
plt.rcParams['font.size'] = 15
sns.pairplot(d[num_cols])
二つの変数間の関係のことを、一般に相関関係と呼ぶが、統計学では直線関係に近い傾向が見られるときに「相関関係がある」ということが多い。
一方の変数の増加につれて他方の変数も増加する傾向を「正の相関関係がある」といい、
逆に一方の変数の増加が他方の変数の減少に対応している場合を「負の相関関係がある」という。
片方、または両方が質的データである場合は分割表と呼ばれる手法を用いる。
pd.crosstab(d['タイプ1'], d['世代'])
d.groupby(['世代'])[num_cols].mean()
分割表で二つのデータの関係を見るには、相対度数を用いる。
相対度数は縦・横・全体に対して3種類できる。どれを用いるかは分析の目的による。
d_agg = pd.crosstab(d['タイプ1'], d['世代'])
print('横方向の相対度数')
display((d_agg.T / d_agg.sum(axis=1)).T)
print('縦方向の相対度数')
display(d_agg / d_agg.sum(axis=0))
print('データ全体の大きさに対する相対度数')
display(d_agg / d_agg.sum())
相関係数
「攻撃」に対する「防御」「素早さ」「特攻」の散布図を示す。
fig, axes = plt.subplots(1, 3, figsize=(14,4))
axes[0].scatter(d['攻撃'],d['防御'])
axes[0].set_xlabel('攻撃')
axes[0].set_ylabel('防御')
axes[1].scatter(d['攻撃'],d['素早さ'])
axes[1].set_xlabel('攻撃')
axes[1].set_ylabel('素早さ')
axes[2].scatter(d['攻撃'],d['特攻'])
axes[2].set_xlabel('攻撃')
axes[2].set_ylabel('特攻')
「防御」「素早さ」「特攻」の順に関係がなくなっていくように見える。
このようなときは相関係数を計算して比較すればよい。
もっともよく用いられるものは、ピアソンの積率相関係数であり、単に相関係数というときはこれを指す。
データが$(x_1,y_1),(x_2,y_2),\cdots,(x_n,y_n)$で与えられた場合、変数$x$と$y$の間の相関係数は
$$
r_{xy}=\frac{\sum (x_i-\bar{x})(y_i-\bar{y})/n}{\sqrt{\sum(x_i-\bar{x})^2/n}\sqrt{\sum(y_i-\bar{y})^2/n}}
=\frac{\sum (x_i-\bar{x})(y_i-\bar{y})}{\sqrt{\sum(x_i-\bar{x})^2}\sqrt{\sum(y_i-\bar{y})^2}}
$$
で定義される。相関係数は常に$-1\leq r_{xy}\leq 1$である。
分母の$\sqrt{\sum(x_i-\bar{x})^2/n}$と$\sqrt{\sum(y_i-\bar{y})^2/n}$は、変数$x,y$のそれぞれの標準偏差$S_x,S_y$である。
分子の$C_{xy}=\sum (x_i-\bar{x})(y_i-\bar{y})/n$は共分散と呼ぶ。
def calc_r(x, y):
dx = x - np.mean(x)
dy = y - np.mean(y)
result = np.sum(dx*dy)/(np.sqrt(np.sum(dx**2))*np.sqrt(np.sum(dy**2)))
return result
print('攻撃-防御 ', calc_r(d['攻撃'],d['防御']))
print('攻撃-素早さ ', calc_r(d['攻撃'],d['素早さ']))
print('攻撃-特攻 ', calc_r(d['攻撃'],d['特攻']))
攻撃-防御 0.45443082983760247
攻撃-素早さ 0.3374489049359815
攻撃-特攻 0.3196123813673278
scipyを利用して計算する。
from scipy.stats import pearsonr
pearsonr(d['攻撃'],d['防御'])
(0.4544308298376027, 5.775292721219386e-47)
相関行列を計算して可視化したものを示す。
R = np.corrcoef(d[num_cols].T)
plt.figure(figsize=(10,10))
sns.heatmap(R, annot=True, square=True, xticklabels=num_cols, yticklabels=num_cols, vmin=-1, vmax=1, cmap='RdBu')
相関係数が最大、最小である$±1$をとるのは、各 $(x_i,y_i)$の関係が、
y_i=±\frac{S_x}{S_y}(x_i-\bar{x})+\bar{y}\ (i=1,\cdots,n)
の時に限られる。
書き換えると$y_i=bx_i+a\ (b=±S_x/S_y,a=\bar{y}-\bar{x}S_x/S_y)$のときであり、すべてのデータが直線$y=bx+a$上にあるとき、相関係数は1または-1となる。
これを正の完全相関($r_{xy}=1,b>0$)、または負の完全相関($r_{xy}=-1,b<0$)という。
一般的には強い相関関係があるからといって、必ずしもその二つのデータの間に因果関係があるということではない。
みかけ上の相関と偏相関係数
2つの事象に因果関係がないのに、共通の要因の影響などにより相関関係があるように見えてしまうことを見かけ上の相関と呼ぶ。
共通の要因(変数3)の影響を除いた後の、変数1と変数2の相関係数は偏相関係数と呼ばれ、次のように計算される。
$$
r_{12・3}=\frac{r_{12} - r_{13}r_{23}}{\sqrt{1-r_{13}^2}\sqrt{1-r_{23}^2}}
$$
def calc_pr(x, y, z):
r12 = calc_r(x, y)
r13 = calc_r(x, z)
r23 = calc_r(y, z)
return (r12 - r13*r23) /(np.sqrt(1-r13**2)*np.sqrt(1-r23**2))
calc_pr(d.重さ, d.高さ, d.HP)
0.5546807100645829
順位相関係数
2つの質的基準がある場合に、観測対象$i$の、2つの基準による順位$R_i$,$r_i'$の間の相関を示す指標である。
順位相関係数は、スピアマンによる定義のものと、ケンドールの定義によるものがしばしば用いられる。
スピアマンの順位相関係数$r_s$は、$R$と$R'$に通常の積率相関係数の式を適用したもので、
$$
r_s=1-\frac{2}{n^3-n}\sum_{i=1}^n(R_i-R_i')^2
$$
これに対しケンドールの順位相関係数$r_K$は、観測対象の対$(i,j)$を考え、次のようにカウントする。
正順、すなわち$R_i>R_j,R_i'>R_j'$あるいは$r_i<R_j,R_i'<R_j'$の場合には$+1$、
逆順、すなわち$R_i>R_j,R_i'R_j'$の場合には$-1$の値を与え、
$+1$を与えた対の数$G$と$-1$を与えた対の数$H$の大小関係をみる。
$$
r_K=\frac{G-H}{n(n-1)/2}
$$
で定義する。
まず、順位の計算を行う。
d_rank = d[['タイプ1','攻撃','防御']].groupby(['タイプ1']).mean().rank(ascending=False)
d_rank
そして、scipyを利用して各相関係数を求める。
from scipy.stats import spearmanr
correlation, pvalue = spearmanr(d_rank['攻撃'], d_rank['防御'])
correlation
0.32094943240454077
from scipy.stats import kendalltau
correlation, pvalue = kendalltau(d_rank['攻撃'], d_rank['防御'])
correlation
0.21568627450980396
直線および平面のあてはめ
最小二乗法によるあてはめを考える。説明変数$x_i$から予想される従属変数$y$の値$bx_i+a$と現実の値$y_i$が、最も小さい隔たりを持つのが、最適な直線$y=bx+a$の引き方である。
従って、二乗和
$$
L=\sum_{i=1}^n{y_i-(bx_i+a)}^2
$$
を最小にする$a,b$を求める。
$a,b$でそれぞれ偏微分して0と置くと、
$$
na+(\sum x_i)b=\sum y_i\
(\sum x_i)a+(\sum x_i^2)b=\sum x_iy_i
$$
となる。これを正規方程式と呼ぶ。
これを$a,b$について解くと
$$
b=\frac{\sum x_iy_i-n\bar{x}\bar{y}}{\sum x_i^2-n\bar{x}^2}\
a=\bar{y}-b\bar{x}
$$
のように得られる。得られた$a,b$による1次式を$y$の$x$への回帰方程式(回帰直線)と呼ぶ。
x = d['高さ']
y = d['重さ']
n = len(x)
x_mean = np.mean(x)
y_mean = np.mean(y)
b = (np.sum(x*y)-n*x_mean*y_mean) / (np.sum(x**2) - n * x_mean**2)
a = y_mean - b * x_mean
print('傾きb: ', b)
print('切片a: ', a)
傾きb: 62.19271304253602
切片a: -9.788573931292724
「高さ」が1 m 高くなると「重さ」は62.2 kg増加する傾向があるとわかる。
$b$が
$$
b=\frac{\sum(x_i-\bar{x})(y_i-\bar{y})}{\sum(x_i-\bar{x})^2}
$$
とも表されることから、
$$
b=r\frac{S_y}{S_x}
$$
という関係が成り立つことがわかる。
相関係数$r$は直線のあてはまりの良さの尺度でもある。
回帰直線による$y$のあてはめ値$\hat{y}_i=bx_i+a$は一般に実際に観測された$y_i$と一致せず、その外れの二乗$\sum d_i^2$について
$$
\sum d_i^2=\sum(y_i-\hat{y}_i)^2=(1-r^2)\sum(y_i-\bar{y})^2
$$
が成り立つ。
$r^2$が1に近いほど$y_iは\hat{y}_i$に近づき、$r^2=1$なら正確に$yy_i=bx_i+a$が成り立ち、$y$は$x$から完全に決定される。
したがって、$r^2$は$x$から$y$を決定する強弱の度合いを表し、$r^2$を決定係数と呼ぶ。
$$
SS_0=\sum(y_i-\bar{y})\
SS_E=\sum(y_i-\hat{y})\
$$
とおくと、$S_0$は被説明変数$y$のばらつき、$SS_E$は回帰をした後のばらつきであるが、
$$
SS_E=(1-r^2)SS_0
$$
であるから、$y$のばらつきは$SS_0$から$SS_E$へちょうそ割合$r^2 %$で減少したことになる。
減少高を
$$
SS_R=SS_0-SS_E
$$
とおくと、
$$
SS_R=r^2SS_0
$$
となり、$r^2$が大きいほど回帰の効果も大きいこととなる。
R2 = R**2
plt.figure(figsize=(10,10))
sns.heatmap(R2, annot=True, square=True, xticklabels=num_cols, yticklabels=num_cols, vmin=-1, vmax=1, cmap='RdBu')
説明変数が二つ以上の場合、重回帰と呼ぶ。$p=2$の場合は、
$$
y=b_1x_1+b_2x_2+a
$$
となるが、これは3次元空間における平面のあてはめの問題となる。
次回
大数の法則と中心極限定理
参考書