線形単回帰・重回帰分析に基づき予測モデルを作成|機械学習・Pythonによるscikit-learnモデル開発・評価

こんにちは、DXCEL WAVEの運営者(@dxcelwave)です!

こんな方におすすめ!
  • Pythonで単回帰および重回帰モデルを作成・評価できるようになりたい
  • scikit-learnで回帰モデルを構築する方法を知りたい
  • 回帰モデルの評価指標について学習したい
目次

単回帰分析・重回帰分析とは

単回帰分析とは、単一の説明変数xに基づいて,連続値をとる目的変数Yを予測する分析手法です。説明変数が1つだけの単回帰モデルの方程式は上図のように示されます。

重回帰分析とは、複数の説明変数x1~xnに基づいて,連続値をとる目的変数Yを予測する分析を指します。

本記事では単回帰および重回帰分析に基づく予測モデルをPythonとscikit-learnを用いてで構築する方法を解説します。

回帰分析の理論的な概要が知りたい方はこちらの記事で解説していますため、合わせてご覧下さい。

【事前準備】回帰モデルを作成に必要なデータを理解する

実際にPythonのScikit-learnライブラリを用いて回帰モデルを構築していきましょう!具体的に下記の手順でプログラムを構築していきます。

  • データセットの説明
  • データの準備
  • データの理解
  • 回帰モデル学習
  • モデル性能評価

データセットの説明

データセットには、Boston Housingというデータセットを活用します。

1970年代後半におけるボストンの住宅情報とその地域における環境をまとめたデータセットであり、住宅価格の予測を目的とした回帰モデル作成のチュートリアルによく利用されます。

Boston Housingデータセットの説明変数および目的変数の概要はそれぞれ以下になります。

説明変数一覧

特徴量名概要
CRIM(地域人工毎の)犯罪発生率
ZN25,000平方フィート以上の住宅区画の割合
INDUS(地域人工毎の)非小売業の土地面積の割合
CHASチャールズ川沿いに立地しているかどうか(該当の場合は1、そうでない場合は0)
NOX窒素酸化物の濃度(単位:pphm)
RM平均部屋数/一戸
AGE1940年よりも古い家の割合
DIS5つのボストン雇用センターまでの重み付き距離
RAD主要な高速道路へのアクセス指数
TAX10,000ドルあたりの所得税率
PTRATIO(地域人工毎の)学校教師1人あたりの生徒数
B(地域人工毎の)アフリカ系アメリカ人居住者の割合
LSTAT低所得者の割合

目的変数一覧

特徴量名概要
MEDV住宅価格の中央値(単位 $1,000)

データの準備

前述したBoston Housingのデータセットを準備します。目的変数と説明変数をPandasのデータフレーム形式で取り扱うために、下記のコードを実行しましょう。

# データ前処理
import numpy as np
import pandas as pd

# データ可視化
import matplotlib.pyplot as plt
import seaborn as sns
plt.style.use("ggplot")
# %matplotlib inline

# グラフの日本語表記対応
from matplotlib import rcParams
rcParams["font.family"]     = "sans-serif"
rcParams["font.sans-serif"] = "Hiragino Maru Gothic Pro"

# データセット読込
from sklearn.datasets import load_boston
boston = load_boston()

# DataFrame作成
df = pd.DataFrame(boston.data)
df.columns = boston.feature_names
df["MEDV"] = boston.target

データの理解

回帰モデル作成に用いる説明変数の特徴および説明変数と目的変数間の関係性を理解するために、グラフや表によってデータを可視化します。

相関係数をヒートマップとして可視化する

numpyのcorrcoef関数を用いて相関行列を作成し、serbonのheatmapメソッドを用いて可視化します。

ここでの相関行列は、標準化された特徴量から算出した共分散行列と等しくなります。

# 相関係数を算出
corr = np.corrcoef(df.values.T)

# ヒートマップとして可視化
hm   = sns.heatmap(
                 corr,                         # データ
                 annot=True,                   # セルに値入力
                 fmt='.2f',                    # 出力フォーマット
                 annot_kws={'size': 8},        # セル入力値のサイズ
                 yticklabels=list(df.columns), # 列名を出力
                 xticklabels=list(df.columns)) # x軸を出力

plt.tight_layout()
plt.show()

今回の重回帰モデルの説明変数は、目的変数(MEDV)と比較的高い相関性(0.45以上)を示したINDUS, RM, TAX, PTRATIO, LSTATを用いることとします。

ここで、説明変数同士の多重共線性が懸念されますが、本記事では回帰モデルの作成と評価方法を中心に解説するため、内容を割愛しています。

説明変数と目的変数の関係性を散布図として可視化(データの分布・外れ値の確認)

前述で指定した説明変数(INDUS, RM, TAX, PTRATIO, LSTAT、MEDV)および目的変数(MEDV)を散布図として可視化します。散布図はseabornのpairplotメソッドを用いて表現できます。

以下散布図をもとにデータの分布および外れ値の有無を確認してみましょう。

# 目的変数(MEDV)と相関が強い上位5つの変数を散布図として可視化
sns.pairplot(df[['INDUS', 'RM', 'TAX', 'PTRATIO','LSTAT', 'MEDV']])
plt.tight_layout()
plt.show()

目的変数(MEDV)は6行6列目のヒストグラムを確認すると、正規分布に従っているように見えるものの、いくつか外れ値が含まれるように見えます。

説明変数(RM)と目的変数(MEDV)の関係を示した6行2列目の散布図を見ると、やや外れ値が気になるものの、2変数の関係は線形として表現されているように見受けられます。

よって、単回帰モデルの説明変数にはRMを用いることとします。

RMvsMEDVの散布図として可視化

参考として、RMおよびMEDVによる回帰直線については、sebornのlmplotメソッドを用いても表現できます。

# 目的変数(MEDV)と相関が強い変数としてRM(部屋数)を選択し可視化
sns.lmplot("RM","MEDV",data=df)
plt.tight_layout()
plt.show()

目的変数のヒストグラム

こちらも参考として、目的変数(MEDV)のヒストグラムの可視化方法を示します。ヒストグラムはmatplotlibを用いて表現可能です。

# 目的変数(MEDV)のヒストグラムを表示
plt.hist(df.MEDV,bins=60)
plt.xlabel("MEDV(住宅価格の中央値)")
plt.ylabel("住宅数")
plt.tight_layout()
plt.show()

【実践】線形単回帰モデルの作成・評価

それでは実際に回帰モデルの作成・評価を行っていきましょう!モデルは単回帰モデルと重回帰モデルをそれぞれ作成し、評価します。

まず単回帰モデルから作成していきます。以下モデル作成に用いる説明変数と目的変数を再掲します。

説明変数RM
目的変数MEDV

【モデル学習】単回帰モデルを作成

単回帰モデルを作成するためのコードを下記に示します。

from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split

# 変数定義
X = df[['RM']].values # 説明変数
y = df['MEDV'].values # 目的変数

# 学習・テストデータ分割
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)

# 回帰モデルのインスタンス
model = LinearRegression(
                         fit_intercept=True, # 切片の計算
                         normalize=False,    # 正規化
                         copy_X=True,        # メモリ内でのデータ複製
                         n_jobs=1,           # 計算に用いるジョブの数
                         positive=False
                        )

# モデル学習
model.fit(X_train, y_train)

# 予測値(Train)
y_train_pred = model.predict(X_train)

# 予測値(Test)
y_test_pred = model.predict(X_test)

【補足】線形回帰モデル作成時のメソッドについて

線形単回帰・重回帰モデルはどちらもscikit-learnのLinearRegression()メソッドを用いて作成可能です。このメソッドに用いる引数についても下記に補足します。

# 回帰モデルのインスタンス
model = LinearRegression(
                         fit_intercept=True,
                         normalize=False,
                         copy_X=True, 
                         n_jobs=1,
                         positive=True
                        )
スクロールできます
引数名概要デフォルト
fit_interceptFalse の場合、切片を計算しない。
目的変数が原点を必ず通る場合に有効。
True
normalizeTrue の場合、説明変数を事前に正規化。
変数間の尺度が異なるときに有効。
False
copy_XTrueの場合、メモリ内でデータを複製して実行。True
n_jobs計算に使うジョブの数。
-1 にすると全ての CPU を用いて計算。
1
positiveTrueに設定すると係数が強制的に正の値となる。
(高次元の配列でのみ利用)
False

【モデル学習】作成したモデルを可視化

前述で作成した回帰モデル(回帰直線)を可視化してみましょう。

# 単回帰直線を可視化
plt.scatter(X_train, y_train, c='blue', s=30)       # 散布図
plt.plot(X_train, y_train_pred, color='red', lw=4)  # 直線

# グラフの書式設定
plt.xlabel('RM(平均部屋数/一戸)')
plt.ylabel('MEDV(住宅価格の中央値)')
plt.tight_layout()
plt.show()

# 回帰係数情報出力
print('偏回帰係数: %.2f' % model.coef_[0])
print('切片: %.2f' % model.intercept_)

# 出力結果
# 偏回帰係数: 9.31
# 切片: -35.99

【評価】RMSEおよび決定係数R2を用いた性能評価

RMSE(平均平方二乗誤差)と決定係数R2を用いて回帰モデルの性能を評価します。下記のコードを実行しましょう!

from sklearn.metrics import r2_score            # 決定係数
from sklearn.metrics import mean_squared_error  # RMSE

# 平均平方二乗誤差(RMSE)
print('RMSE 学習: %.2f, テスト: %.2f' % (
        mean_squared_error(y_train, y_train_pred, squared=False), # 学習
        mean_squared_error(y_test, y_test_pred, squared=False)    # テスト
      ))

# 決定係数(R^2)
print('R^2 学習: %.2f, テスト: %.2f' % (
        r2_score(y_train, y_train_pred), # 学習
        r2_score(y_test, y_test_pred)    # テスト
      ))

# 出力結果
# RMSE 学習: 6.49, テスト: 6.86
# R^2  学習: 0.50, テスト: 0.44

上記決定係数R2に着目すると、学習データを使ったR2は0.5、テストデータを使ったR2は0.44であり、さほど性能は良くない結果となりました。

ここで回帰モデルの評価指標について詳しく知りたい方はこちらの記事で解説しています。

【評価】残差プロットによる可視化

単回帰モデルの予測結果を残差プロットとして表現する方法を下記に示します。良い回帰モデルは、サンプルが残差の中央値0近辺にランダムの分布しますが、今回作成した回帰モデルは、誤差の大きいプロットが目立つように見受けられれます。

# 予測値と残差をプロット(学習データ)
plt.scatter(y_train_pred,             # グラフのx値(予測値)  
            y_train_pred - y_train,   # グラフのy値(予測値と学習値の差)
            c='blue',                 # プロットの色
            marker='o',               # マーカーの種類
            s=40,                     # マーカーサイズ
            alpha=0.7,                # 透明度
            label='学習データ')         # ラベルの文字


# 予測値と残差をプロット(テストデータ)
plt.scatter(y_test_pred,            
            y_test_pred - y_test, 
            c='red',
            marker='o', 
            s=40,
            alpha=0.7,
            label='テストデータ')

# グラフの書式設定
plt.xlabel('予測値')
plt.ylabel('残差')
plt.legend(loc='upper left')
plt.hlines(y=0, xmin=-20, xmax=60, lw=2, color='black')
plt.xlim([-10, 50])
plt.ylim([-40, 20])
plt.tight_layout()
plt.show()

【実践】重回帰モデルの作成・評価

続いて、重回帰モデルを作成していきます。以下モデル作成に用いる説明変数と目的変数を再掲します。

説明変数INDUS, RM, TAX, PTRATIO, LSTAT
目的変数MEDV

【モデル学習】重回帰モデルを作成

重回帰モデルを作成するためのコードを下記に示します。重回帰モデルは、単回帰モデル同様に、scikit-learnのLinearRegression()メソッドを用いて作成できます。

from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split

# 変数定義
X = df[['INDUS','RM','TAX','PTRATIO','LSTAT']].values # 説明変数
y = df['MEDV'].values                                 # 目的変数(住宅価格の中央値)

# 学習データ/テストデータ分割
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=0)

# 重回帰のインスタンス
model_multi = LinearRegression()

# モデル学習
model_multi.fit(X_train, y_train)

# 予測値(学習データ)
y_train_pred = model_multi.predict(X_train)

# 予測値(テストデータ)
y_test_pred  = model_multi.predict(X_test)

【モデル学習】作成したモデルの偏回帰係数を導出

作成した重回帰モデルの偏回帰係数を確認する方法を下記に示します。

# 偏回帰係数
df_coef = pd.DataFrame(model_multi.coef_.reshape(1,5), 
                       columns=['INDUS','RM','TAX','PTRATIO','LSTAT'])
print(df_coef)

# 切片
print('切片: %.2f' % model_multi.intercept_)


# 出力結果(偏回帰係数)
#      INDUS        RM       TAX   PTRATIO     LSTAT
#   0.056633  4.627714 -0.005232 -1.017411 -0.526502

# 出力結果(切片)
# 20.37

【評価】RMSEおよび決定係数R2を用いた性能評価

RMSEと決定係数R2を用いて重回帰モデルの性能を評価していきます。下記のコードを実行してみましょう。

from sklearn.metrics import r2_score            # 決定係数
from sklearn.metrics import mean_squared_error  # RMSE

# 平均平方二乗誤差(RMSE)
print('RMSE 学習: %.2f, テスト: %.2f' % (
        mean_squared_error(y_train, y_train_pred, squared=False), # 学習
        mean_squared_error(y_test, y_test_pred, squared=False)    # テスト
      ))

# 決定係数(R^2)
print('R^2 学習: %.2f, テスト: %.2f' % (
        r2_score(y_train, y_train_pred), # 学習
        r2_score(y_test, y_test_pred)    # テスト
      ))

# 出力結果
# RMSE 学習: 4.92, テスト: 5.83
# R^2  学習: 0.71, テスト: 0.59

重回帰モデルの決定係数R2を単回帰モデルと比較すると、学習データを使ったR2は0.5→0.71、テストデータを使ったR2は0.44→0.59に改善されました。性能としてはまだまだ良いとは言えませんが、今回のデータセットの場合、重回帰モデルを用いると目的変数の説明性が向上すると言えます。

【評価】残差プロットによる可視化

最後に重回帰モデルの残差プロットの描画方法を下記に示します。

# 予測値と残差をプロット(学習データ)
plt.scatter(y_train_pred,             # グラフのx値(予測値)  
            y_train_pred - y_train,   # グラフのy値(予測値と学習値の差)
            c='blue',                 # プロットの色
            marker='o',               # マーカーの種類
            s=40,                     # マーカーサイズ
            alpha=0.7,                # 透明度
            label='学習データ')         # ラベルの文字


# 予測値と残差をプロット(テストデータ)
plt.scatter(y_test_pred,            
            y_test_pred - y_test, 
            c='red',
            marker='o', 
            s=40,
            alpha=0.7,
            label='テストデータ')

# グラフ書式設定
plt.xlabel('予測値')
plt.ylabel('残差')
plt.legend(loc='upper left')
plt.hlines(y=0, xmin=-20, xmax=60, lw=2, color='black')
plt.xlim([-20, 60])
plt.ylim([-40, 20])
plt.tight_layout()
plt.show()

【参考】AI・機械学習における配信情報まとめ

当サイトではAI・機械学習における「基礎」から「最新のプログラミング手法」に至るまで幅広く解説しております。また「おすすめの勉強方法」をはじめ、副業・転職・フリーランスとして始める「AI・機械学習案件の探し方」についても詳しく言及しています。

【仕事探し】副業・転職・フリーランス

【教育】おすすめ勉強法

【参考】記事一覧

最後に

この記事が気に入ったら
フォローしてね!

本記事をシェア!
目次