Python scikit-learnによる重回帰分析(データの標準化/正規化含む)

 重回帰分析と言ってもソルバーは複数ある。本記事では、sklearn(サイキットラーン)のLinearRegressionとSGDRegressorを用いた2つについて記載した。得られる回帰式は文字通り、Y=f(Xi)=a1x1+a2x2+…+aixi+eで表し、共にlinear_modelをインポートして利用する。

分析データは、sklearnに含まれているデータセットbostonを使用した。目的変数Yはボストン近郊の住宅価格で、説明変数Xiは部屋数や犯罪率など13因子ある。

下図は縦軸に住宅価格、横軸に部屋数とした場合のグラフである。青点は生データを表し、赤点はLinearRegressionを用いた重回帰分析により得られた回帰式Y=f(Xi)に生データのXiを代入してプロットしたものである。部屋数との相関が高いことがわかる。

f:id:HK29:20180514012334p:plain

そして、重回帰分析により得られた上図の赤点の回帰式の係数,切片,決定係数を下図に示す。決定係数R^2は74%でありそこそこ高い結果。

f:id:HK29:20180514014034p:plain

 

▼本プログラム

重回帰分析の実行例として3つの方法を記載した。1. 普通にした場合、2. 説明変数Xを標準化した場合、3. 説明変数Xを正規化した場合、且つscikit-learn内の他のソルバーSGDを使用した場合。ちなみに、SGDの場合はオプション設定が多々あり、例えば、本コード中にある試行回数を変更すると結果は若干変わる。

# -*- coding: utf-8 -*-
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

from sklearn import datasets
from sklearn import linear_model
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import MinMaxScaler

############################################ load and set sklearn dataset
# load boston
boston = datasets.load_boston()
print(boston)
print(boston.DESCR)
print(boston.data)
print(boston.feature_names)
print(boston.target)

# set df
boston_df=pd.DataFrame(boston.data)
boston_df.columns = boston.feature_names
boston_df['PRICE'] = pd.DataFrame(boston.target)
# 今回使用しないがcsvファイルへ出力する方法
boston_df.to_csv('boston_df.csv', sep=',', index=True, encoding='utf-8')
boston_df.to_csv('boston.csv', sep=',', index=False, encoding='utf-8')

############################################ input
# set X
X = boston_df.drop("PRICE", axis=1)
print("##### X ###############################")
print("X.describe() ->")
print(X.describe())

# set X2 :Xを標準化する場合 Standardization
scaler = StandardScaler()
scaler.fit(X)
print("##### X2 ##############################")
print("scaler.fit(X) ->" + str(scaler.fit(X)))
print("scaler.mean_ ->" + str(scaler.mean_))
print("scaler.transform ->" + str(scaler.transform(X)))
X2 = scaler.transform(X)
X2 = pd.DataFrame(X2)
X2.columns = boston.feature_names
print("X2.describe() ->")
print(X2.describe())

# set X3 :Xを正規化する場合 Normalization
mscaler = MinMaxScaler()
mscaler.fit(X)
print("##### X3 ##############################")
print("mscaler.fit(X) ->" + str(mscaler.fit(X)))
print("mscaler.data_max_ ->" + str(mscaler.data_max_))
print("mscaler.data_range_ ->" + str(mscaler.data_range_))
print("mscaler.transform ->" + str(mscaler.transform(X)))
X3 = mscaler.transform(X)
X3 = pd.DataFrame(X3)
X3.columns = boston.feature_names
print("X3.describe() ->")
print(X3.describe())

# set Y
Y = boston_df.PRICE 

############################################ Execute
# 重回帰分析
clf = linear_model.LinearRegression(fit_intercept=True, normalize=False, copy_X=True, n_jobs=1)
clf.fit(X, Y)
# Xを標準化したデータX2での重回帰分析
clf2 = linear_model.LinearRegression(fit_intercept=True, normalize=False, copy_X=True, n_jobs=1)
clf2.fit(X2, Y)
# Xを正規化したデータX3で且つ他のアルゴリズム(SGDRegressor)を利用した場合の重回帰分析
# 本家のマニュアルをみればわかるが、細かい設定が出来る
clf3_SGD = linear_model.SGDRegressor(max_iter=500)
clf3_SGD.fit(X3, Y)

############################################ print out
# out put
print("### out put from X data ##############################")
print(pd.DataFrame({"Name":X.columns,
                    "Coefficients":clf.coef_}).sort_values(by='Coefficients') ) # 係数でソート
print(clf.intercept_) # 切片
print(clf.score(X,Y)) # R^2

# out put2
print("### out put2 from X2 data ############################")
print(pd.DataFrame({"Name":X2.columns,
                    "Coefficients":clf2.coef_}).sort_values(by='Coefficients') )
print(clf2.intercept_)
print(clf2.score(X2,Y))

# out put3
print("### out put3 from X3 data and used SGD ###############")
print(pd.DataFrame({"Name":X3.columns,
                    "Coefficients":clf3_SGD.coef_}).sort_values(by='Coefficients') )
print(clf3_SGD.intercept_)
print(clf3_SGD.score(X3,Y))

############################################ plot
# out put
plt.figure(1)
plt.title('raw X-data')
plt.xlabel('RM (number of rooms)', fontsize=14)
plt.ylabel('PRICE (target)', fontsize=14)
plt.scatter(X["RM"], Y, c='blue', label='Raw data')
plt.scatter(X.RM, clf.predict(X), c='red', label='Multiple regression analysis')
plt.legend(loc='lower right', fontsize=12)
#plt.show()

# out put2
plt.figure(2)
plt.title('Standardization of X-data')
plt.xlabel('RM (number of rooms)', fontsize=14)
plt.ylabel('PRICE (target)', fontsize=14)
plt.scatter(X2.RM, Y, c='blue', label='Raw data')
plt.scatter(X2.RM, clf2.predict(X2), c='red', label='Multiple regression analysis')
plt.legend(loc='lower right', fontsize=12)
#plt.show()

# out put3
plt.figure(3)
plt.title('Normalization of X-data and used SGD')
plt.xlabel('RM (number of rooms)', fontsize=14)
plt.ylabel('PRICE (target)', fontsize=14)
plt.scatter(X3.RM, Y, c='blue', label='Raw data')
plt.scatter(X3.RM, clf3_SGD.predict(X3), c='red', label='Multiple regression analysis')
plt.legend(loc='lower right', fontsize=12)
plt.show()

そして、説明変数Xを標準化した場合と正規化した場合の各々の重回帰分析の結果を以降に示す。共に、目的変数Yに対する説明変数Xの全ての感度を揃えることを意図している。

▼説明変数Xを標準化した場合(本プログラム中ではX2データを指す)

標準化とは「各生データの平均値との差」を標準偏差で割ったものである。横軸の中心が0になっている。この意味は説明変数Xの全ての平均を0にして標準偏差が1になるようデータの分散を揃えている。

f:id:HK29:20180520122412p:plain

下図をみると、決定係数R^2は74%で変化はない。そして、係数はRM(部屋数)が一番感度が高いことに変化はないが、その他の説明変数の順番は微妙に変化している。

f:id:HK29:20180520123604p:plain

▼説明変数Xを正規化した場合(本プログラム中ではX3データを指す)

正規化とは「各生データの最小値との差」を「最大値と最小値の差」で割ったものを指す。横軸が0~1になっている。この意味は説明変数Xの全てのデータ範囲が0~1になるように揃えている。

f:id:HK29:20180520122429p:plain

f:id:HK29:20180520123747p:plain

どの方法が適切かはない。ただ、扱うデータを標準化もしくは正規化した方が、理屈的に妥当なのは理解できる。結局、重回帰分析の目的は回帰式を作成することであり、また得られた各説明変数の係数から感度が高い因子を見つけることである。今回扱ったbostonの住宅価格に与える感度が一番高いのは部屋数(RM)であることは言える。

▼参考リンク

sklearn.linear_model.LinearRegression — scikit-learn 0.19.1 documentation

sklearn.linear_model.SGDRegressor — scikit-learn 0.19.1 documentation

sklearn.preprocessing.StandardScaler — scikit-learn 0.19.1 documentation

sklearn.preprocessing.MinMaxScaler — scikit-learn 0.19.1 documentation

hk29.hatenablog.jp

以上

<広告>