『Rによる計量経済学』第8章「多重共線性」をPythonで実行する。
テキスト付属データセット(「k0801.csv」等)については出版社サイトよりダウンロードしてください。
また、以下の説明は本書の一部を要約したものですので、より詳しい説明は本書を参照してください。
次のような重回帰モデルを推定。
$Y_{i} = \alpha + \beta X_{i} + \gamma Z_{i} + u_{i} (i = 1, 2, ..., n)$
回帰係数$\hat{\beta}$の分散である$s_{\hat{\beta}}^{2}$は次の3つの部分からなる。
$s_{\hat{\beta}}^{2} = \frac{s^{2}}{\Sigma (X_{i} - X)^{2}(1 - r_{XZ}^{2})}$
$= \frac{(A)}{(B)[1-(C)]}$
ここで(C)が理由で回帰係数$\hat{\beta}$が有意にならないことを多重共線性(multi collinearity)の問題があると呼ぶ。
① 他の説明変数を除いた場合に説明力があるかを確認する。
② (A)の部分の計算
③ (B)の部分の計算
④ (C)の部分の計算
⑤ ①〜④の結果から④が原因であるかどうか判断する。
「k0801.csv」を分析。
In [1]:
%matplotlib inline
In [2]:
# -*- coding:utf-8 -*-
from __future__ import print_function
import numpy as np
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings('ignore')
In [3]:
# データ読み込み
data = pd.read_csv('example/k0801.csv')
data
Out[3]:
In [4]:
# 説明変数設定
X = data[['X', 'Z']]
X = sm.add_constant(X)
X
Out[4]:
In [5]:
# 被説明変数設定
Y = data['Y']
Y
Out[5]:
In [6]:
# OLSの実行(Ordinary Least Squares: 最小二乗法)
model = sm.OLS(Y,X)
results = model.fit()
print(results.summary())
これから、Xの係数 $\hat{\beta}$ のP値が小さく有意水準5%でも有意ですが、Zの係数 $\hat{\gamma}$ のP値が大きく明らかに有意でないことがわかります。
In [7]:
# 説明変数設定
X = data['X']
X = sm.add_constant(X)
X
Out[7]:
In [8]:
# 説明変数をXのみにして推定
# OLSの実行(Ordinary Least Squares: 最小二乗法)
model = sm.OLS(Y,X)
results = model.fit()
print(results.summary())
In [9]:
# 説明変数設定
X = data['Z']
X = sm.add_constant(X)
X
Out[9]:
In [10]:
# 説明変数をZのみにして推定
# OLSの実行(Ordinary Least Squares: 最小二乗法)
model = sm.OLS(Y,X)
results = model.fit()
print(results.summary())
In [11]:
# 標準偏差
data.std()
Out[11]:
In [12]:
# 平均
data.mean()
Out[12]:
In [13]:
# 相関係数
data.corr()
Out[13]:
以上の結果から、重回帰分析において $Z_{i}$ は有意ではなかったが、多重共線性によるものだったとうかがえる。
自由度調整済決定係数 $\overline{R}^{2}$ と同様に、回帰係数の実質的な説明力を表す。
しかし、Statsmodelsにはステップワイズ法を使って変数選択をしてくれるメソッドがないらしいので、ここのページからforward_selected()を借りて変数選択を行う。
In [14]:
# データ読み込み
data2 = pd.read_csv('example/k0802.csv')
data2
Out[14]:
In [15]:
import statsmodels.formula.api as smf
def forward_selected(data, response):
"""Linear model designed by forward selection.
Parameters:
-----------
data : pandas DataFrame with all possible predictors and response
response: string, name of response column in data
Returns:
--------
model: an "optimal" fitted statsmodels linear model
with an intercept
selected by forward selection
evaluated by adjusted R-squared
"""
remaining = set(data.columns)
remaining.remove(response)
selected = []
current_score, best_new_score = 0.0, 0.0
while remaining and current_score == best_new_score:
scores_with_candidates = []
for candidate in remaining:
formula = "{} ~ {} + 1".format(response,
' + '.join(selected + [candidate]))
score = smf.ols(formula, data).fit().rsquared_adj
scores_with_candidates.append((score, candidate))
scores_with_candidates.sort()
best_new_score, best_candidate = scores_with_candidates.pop()
if current_score < best_new_score:
remaining.remove(best_candidate)
selected.append(best_candidate)
current_score = best_new_score
formula = "{} ~ {} + 1".format(response,
' + '.join(selected))
model = smf.ols(formula, data).fit()
return model
In [16]:
# 第一引数にdata、第二引数に被説明変数のコラム名
model = forward_selected(data2, 'Y')
print(model.model.formula)
In [17]:
# OLS実行
X = data2[['X5', 'X6', 'X8']]
X = sm.add_constant(X)
Y = data2['Y']
model = sm.OLS(Y,X)
results = model.fit()
print(results.summary())