ここでは、３つ以上の変数をあつかう線形の重回帰分析を学びます。
アメリカのマサチューセッツ州北東部にある大都市、ボストン市の住宅価格をさまざまな説明変数を用いて分析してみます。
####⑴ ライブラリをインポートする
# 数値計算に必要なライブラリ
import numpy as np
import pandas as pd
# グラフを描画するパッケージ
import matplotlib.pyplot as plt
# 機械学習ライブラリscikit-learnの線形モデル
from sklearn import linear_model
####⑵ データを読み込む
scikit-learnに付属しているデータセットの一つ「Boston house prices dataset」を読み込み、変数bostonとして格納します。
from sklearn.datasets import load_boston
boston = load_boston()
####⑶ データの中身を確認し、分析用に加工する
print(boston)
まず中身を表示して、データを概観します。
print(boston.DESCR)
ボストン・データセットの解説を確認してみます。
データセットには説明変数として13項目、また目的変数として1項目が用意されています。
| 変数名 | もとの語 | 定義 |
|---|---|---|
| CRIM | crime rate | 町別の一人当たり犯罪率 |
| ZN | zone | 25,000平方フィート超に区画された住宅地の割合 |
| INDUS | industry | 町あたりの非小売業（製造業）の面積の割合 |
| CHAS | Charles River | チャールズ川ダミー変数（地域が川に接している場合は1、そうでない場合は0） |
| NOX | nitric oxides | 一酸化窒素濃度（1000万分の1） |
| RM | rooms | 一戸当たり平均室数 |
| AGE | ages | 1940年以前に建てられた持ち家の割合 |
| DIS | distances | ボストンの5つの雇用センターまでの加重距離 |
| RAD | radial highways | 幹線道路へのアクセス性指数 |
| TAX | tax rate | 10,000ドル当たりの固定資産税率 |
| PTRATIO | pupil-teacher ratio | 町別の生徒・教師比率 |
| B | blacks | 町別の黒人の割合 |
| LSTAT | lower status | 人口に占める下層階級の割合 |
| MEDV | median value | 1000ドル台の持家の中央値 |
####⑷ データをPandasのデータフレームに変換する
boston_df = pd.DataFrame(boston.data)
print(boston_df)
[506 rows x 13 columns] とあり、データの形状は506行×13列、すなわち13変数について標本数は506となっています。
####⑸ カラム名を指定し、目的変数を追加する
# カラム名を指定
boston_df.columns = boston.feature_names
print(boston_df)
boston_dfのcolumnsとして、bostonのfeature_namesを指定します。
# 目的変数を追加
boston_df['PRICE'] = pd.DataFrame(boston.target)
print(boston_df)
boston_dfにカラム名PRICEとして、bostonのデータtargetをPandasのデータフレームに変換して格納します。最右端の列に目的変数となる「PRICE」が加わります。
##scikit-learnの基本的な文法
これからscikit-learnを利用して重回帰分析を行ないますが、その手順を整理しておきます。
①モデルの元（インスタンス）を作成
モデルの変数名 = LinearRegression()
②説明変数Xと目的変数Yをもとにモデルを作成
モデルの変数名.fit(X, Y)
③作成したモデルを使って回帰係数を算出
モデルの変数名.coef_
④作成したモデルを使って切片を算出
モデルの変数名.intercept_
⑤決定係数を算出してモデルの精度を求める
モデルの変数名.score(X, Y)
####⑹ 説明変数、目的変数を作成する
# 目的変数のみ削除して変数Xに格納
X = boston_df.drop("PRICE", axis=1)
# 目的変数のみ抽出して変数Yに格納
Y = boston_df["PRICE"]
print(X)
print(Y)
データフレームの列を、drop関数に引数 ("カラム名", axis=1) を指定して削除する。
さらに、データフレームの列をデータフレーム名["カラム名"]で指定して抽出し、変数Yに格納する。
####⑺ インスタンスを作成する①
model = linear_model.LinearRegression()
####⑻ モデルを作成する②
model.fit(X,Y)
####⑼ 回帰係数を算出する③
model.coef_
coef_ は係数（coefficient）の意味です。
いささか係数がわかりにくいので、データフレームに変換して見やすくします。
# 変数coefficientに係数の値を格納
coefficient = model.coef_
# データフレームに変換し、カラム名とインデックス名を指定
df_coefficient = pd.DataFrame(coefficient,
columns=["係数"],
index=["犯罪率", "宅地率", "製造業比率", "チャールズ川", "一酸化窒素濃度",
"平均室数", "持ち家率", "雇用センター", "幹線道路", "固定資産税率",
"生徒･教師比率", "黒人比率", "下層階級比率"])
df_coefficient
####⑽ 切片を算出する④
model.intercept_
####⑾ 決定係数を算出して精度を確認する⑤
model.score(X, Y)
##クロスバリデーション
モデルの精度、すなわち解析そのものの妥当性を検証する方法として、交差検証（クロスバリデーション cross-validation）という考え方があります。
多くは、サンプルをランダムに２つの群に分けて、まず１群をつかってモデルを作成し、残りの１群をつかってそのモデルをテストするというものです。
前者をトレーニングデータ、後者をテストデータと呼びます。
さて、scikit-learnには、トレーニングデータとテストデータの分割を行うメソッドが用意されています。
sklearn.model_selection.train_test_split
トレーニングとテストをどういう比率で分けるかという決まったルールはありませんが、train_test_split の引数に指定をしなければ４分の１すなわち全サンプルの25%分がテストデータとして分割されます。
####➊ サンプルを分割する
# sklearnのtrain_test_splitメソッドをインポート
from sklearn.model_selection import train_test_split
# 変数X, Yをそれぞれトレーニング用とテスト用に分割
X_train, X_test, Y_train, Y_test = train_test_split(X,Y)
変数Xのトレーニング用・テスト用に分割されたそれぞれの中身を確認してみます。
print("変数xのトレーニング用データ:", X_train, sep="\n")
print("変数xのテスト用データ:", X_test, sep="\n")
####➋ インスタンスを作成する①
変数名multi_lregは、multiple linear regression analysis（重回帰分析）の略です。
multi_lreg = linear_model.LinearRegression()
####➌ トレーニングデータでモデルを作成する②
multi_lreg.fit(X_train, Y_train)
####➍ トレーニングデータでの決定係数を算出する⑤
multi_lreg.score(X_train, Y_train)
####➎ テストデータでの決定係数を算出する⑤
multi_lreg.score(X_test,Y_test)
トレーニングとテストという一対ですが、これは「既知のデータ」と「未知のデータ」と置き換えることができます。
既に持っているデータでモデルを作成し、そのモデルを、新たに得られるデータに適用した場合どの程度それが通用するのか。
常に、分析対象となるデータは全体の一部分であり、過去形もしくは現在完了形のものです。ただ、なにも現状確認のためだけに分析をしているわけではなく、そこから「先を読む」「今後を占う」という狙いがあるはずです。
その意味で、決定係数は、低いとなればモデル自体の精度が疑われるわけですが、かといって高ければ高いほど良いともいえません。むしろ高過ぎることは問題です。
決定係数が高い、すなわち残差（誤差）が小さいということは、それだけモデルが分析対象としたデータに適合しているということです。それが適合し過ぎているとなると、今度新たに得られたデータに適用したら決定係数はガタ落ちということになりかねません。いわゆる「過学習」という問題です。
次に、重回帰分析の計算のしくみを理解するために、敢えて便利な scikit-learn を利用しないで重回帰分析を行ないます。