こちらのコンテンツは最新の情報ではない可能性があります。無料で学べる最新のコンテンツは Python&機械学習入門コース や 脱ブラックボックスコース をぜひご利用ください。

機械学習 実践（教師あり学習：回帰）

本章では、主に重回帰分析のアルゴリズムを用いて回帰を実装していきます。重回帰分析の理解が不十分な方は Chainer チュートリアルを先に読んでください。

本章で特に重要な論点として過学習という概念があります。教師あり学習はコンピュータに大量のデータから法則性を学習させるための手法であることはお伝えしましたが、学習に使用したデータに対して過剰に適合してしまうと新たな入力データに対して予測性能がうまく発揮できない現象が起きます。その学習済みモデルが未知のデータに対してどの程度予測がうまく行くのかの性能を汎化性能といい、もしもモデルが未知のデータに対してもうまく予測ができているようであれば、汎化性能の高いモデルであると呼ぶことができます。

また、機械学習の基礎では多くのアルゴリズムの中身についてお伝えいたしませんが、今後公開予定のテーブルデータ特化編ではより深い実践的な内容を紹介していきます。そのための準備として知識を付けていきましょう。

本章の構成

  • 回帰
  • 重回帰分析の復習
  • 重回帰分析の実装
  • 線形回帰の過学習を抑制する手法
  • 相関関係と多重共線性問題

回帰

回帰とはどういった問題でしょうか。

16_1

上記の図のように、手元にある実測値に対してアルゴリズムを選択し、実測値を近似できるように当てはめることです。特に図のような直線で表現することを線形回帰といいます。有名なアルゴリズムとしては、

  • 単回帰分析
  • 重回帰分析
  • Ridge 回帰（リッジ回帰）
  • Lasso 回帰（ラッソ回帰）

があり、本章ではこれらを実装と共に紹介していきます。

しかし、データによっては直線で表現しきれないこともあります。

16_2

上記の図をご覧ください。点線で実測値を表現しようとすると誤差が大きいことが直感的にわかります。そのような場合に、実線で表すような直線ではない形で表現するアルゴリズムを非線形回帰と呼びます。有名なアルゴリズムとしては、

  • 決定木（回帰木）
  • ランダムフォレスト
  • ニューラルネットワーク

などがあり、これらは分類でも使用されることが多いアルゴリズムです。回帰と分類で同じアルゴリズムが使われているということも覚えておきましょう。

本章で紹介するアルゴリズムの他にも、Kaggle などの世界的な分析プラットフォームでは勾配ブースティングと呼ばれる手法が多く使用されており、有名な例として XGBoost や LightGBM があります。これらはこの先のテーブルデータ特化編で詳しく紹介する予定ですので、本章では扱いません。

重回帰分析の復習

簡単に重回帰分析の復習を行いましょう。重回帰分析とは、2 つ以上の変数から回帰分析を行う多変量解析の一つです。

例えばあなたは新しく住む物件を探していて、その物件の家賃を予測することを考えているとしましょう。あなたであれば、どのような特徴を考えるでしょうか。例えば、

  • 部屋の広さ
  • 駅からの距離
  • 犯罪発生率

上記のような変数を考慮することがあると思います。今回は、部屋の広さを x1x_1、駅からの距離を x2x_2、⋯\cdots、犯罪発生率 xMx_M といった風に MM 個の入力変数を扱うことにしましょう。このとき重回帰分析のモデルは、パラメータ wmw_m と入力変数 xmx_m、バイアス bb を用いて、

y=w1x1+w2x2+⋯+wMxM+b=∑m=1Mwmxm+b\begin{aligned} y &= w_1x_1 + w_2x_2 + \cdots + w_Mx_M + b \\ &= \sum_{m=1}^{M} w_mx_m + b \end{aligned}

と定義できます。また、入力変数が一つの場合が単回帰分析だったことも覚えておきましょう。

ここでバイアス bb の扱い方を考えます。今回は w0=bw_0=b、x0=1x_0=1 とおくことで次のように bb を総和の内側の項に含めて、簡潔に表記できるようにします。

y=w1x1+w2x2+⋯+wMxM+b=w1x1+w2x2+⋯+wMxM+w0x0=w0x0+w1x1+⋯+wMxM=∑m=0Mwmxm\begin{aligned} y &= w_1x_1 + w_2x_2 + \cdots + w_Mx_M + b \\ &= w_1x_1 + w_2x_2 + \cdots + w_Mx_M + w_0x_0 \\ &= w_0x_0 + w_1x_1 + \cdots + w_Mx_M \\ &= \sum_{m=0}^{M} w_m x_m \end{aligned}

さらに線形代数で学んだベクトルの内積で表記するとシンプルに記述できます。

y=w0x0+w1x1+w2x2+⋯+wMxM=[w0w1⋯wM][x0x1⋮xM]=wTx\begin{aligned} y &= w_0x_0 + w_1x_1 + w_2x_2 + \cdots + w_Mx_M \\\\ &= \begin{bmatrix} w_0 & w_1 & \cdots & w_M \end{bmatrix} \begin{bmatrix} x_0 \\ x_1 \\ \vdots \\\\ x_M \end{bmatrix}\\\\ &= \mathbf{w}^\mathbf{T} \mathbf{x} \end{aligned}

このモデルが持つパラメータは M+1M+1 個あり、M+1M+1 次元のベクトル w\mathbf{w} で表されています。重回帰分析では、この w\mathbf{w} のすべての要素について最適な値を求めるのでした。

目的関数の最適化

目的関数は単回帰分析と同様に二乗和誤差を用います。目標値を tnt_n、予測値を yny_n、nn はデータ数として、二乗和誤差をベクトルの内積を用いると、

L=(t1−y1)2+(t2−y2)2+⋯+(tN−yN)2=[t1−y1t2−y2⋯tN−yN][t1−y1t2−y2⋮tN−yN]=(t−y)T(t−y)\begin{aligned} \mathcal L &= (t_1 - y_1)^2 + (t_2 - y_2)^2 + \cdots + (t_N - y_N)^2 \\\\ &= \begin{bmatrix} t_1-y_1 & t_2-y_2 & \cdots & t_N-y_N \end{bmatrix} \begin{bmatrix} t_1-y_1 \\\\ t_2-y_2 \\\\ \vdots \\\\ t_N-y_N \end{bmatrix} \\\\ &= (\mathbf{t} - \mathbf{y})^\mathbf{T}(\mathbf{t} - \mathbf{y}) \end{aligned}

とできます。この目的関数 L\mathcal L を最小化するパラメータベクトル w\mathbf{w} を求めることが重回帰分析の目的です。途中の式変換に関しては Chainer チュートリアルで説明されているためこちらで紹介しませんが、

w=(XTX)−1XTt\begin{array}{c} \mathbf{w} = (\mathbf{X^T X})^{-1} \mathbf{X^T t} \end{array}

が導かれます。この式を正規方程式 (normal equation) と呼びます。

上式は、与えられたデータを並べた行列 X\mathbf{X} と、各データの目標値を並べたベクトル t\mathbf{t} から最適なパラメータ w\mathbf{w} を計算しており、新しい入力値 x\mathbf{x} =[x1,⋯,xM]T=[x_1, \cdots, x_M]^T に対して予測値 yy を得るためには、ここで求めたパラメータ w\mathbf{w} を用いて、

y=wTx\begin{array}{c} y=\mathbf{w^T x} \end{array}

のように計算すると良いことがわかります。

重回帰分析の実装

まずは scikit-learn の復習も兼ねて、重回帰分析の実装を行います。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

データセットの準備

まず、scikit-learn の中にあるデータセットを用いて、重回帰分析の実装を行います。sklearn.datasets モジュール内の load_boston メソッドを使用してボストン近郊の住宅データを取り込みます。

Pandas の DataFrame 型で表示したい方は下記の一連のコードを実行することで作成することができます。

# データセットの読み込み
from sklearn.datasets import load_boston
dataset = load_boston()
x, t = dataset.data, dataset.target
columns = dataset.feature_names

scikit-learn で用意されているデータセットは、確認すると NumPy の ndarray で保存されていることがわかります。

type(x), x.shape, type(t), t.shape
(numpy.ndarray, (506, 13), numpy.ndarray, (506,))

Pandas の DataFrame に変換を行います。また、引数で列名を指定しています。

df = pd.DataFrame(x, columns=columns)
df.head(3)
CRIM ZN INDUS CHAS NOX RM AGE DIS RAD TAX PTRATIO B LSTAT
0 0.00632 18.0 2.31 0.0 0.538 6.575 65.2 4.0900 1.0 296.0 15.3 396.90 4.98
1 0.02731 0.0 7.07 0.0 0.469 6.421 78.9 4.9671 2.0 242.0 17.8 396.90 9.14
2 0.02729 0.0 7.07 0.0 0.469 7.185 61.1 4.9671 2.0 242.0 17.8 392.83 4.03

DataFrame に先程読み込んだ目標値が格納されている t を Target という名前で列の後方に追加します。

# 目標値を追加
df['Target'] = t
df.head(3)
CRIM ZN INDUS CHAS NOX RM AGE DIS RAD TAX PTRATIO B LSTAT Target
0 0.00632 18.0 2.31 0.0 0.538 6.575 65.2 4.0900 1.0 296.0 15.3 396.90 4.98 24.0
1 0.02731 0.0 7.07 0.0 0.469 6.421 78.9 4.9671 2.0 242.0 17.8 396.90 9.14 21.6
2 0.02729 0.0 7.07 0.0 0.469 7.185 61.1 4.9671 2.0 242.0 17.8 392.83 4.03 34.7

それぞれの変数ごとの説明は下記の表を参考にしてください。

列名 説明
CRIM 人口 1 人あたりの犯罪発生率
ZN 25,000 平方フィート以上の住宅区画が占める割合
INDUS 小売業以外の商業が占める面積の割合
CHAS チャールズ川の川沿いかどうか（1: 川の周辺, 0: それ以外）
NOX 窒素酸化物の濃度
RM 住居の平均部屋数
AGE 1940 年より前に建てられた持ち主が住んでいる物件の割合
DIS 5 つのボストン雇用施設からの重み付き距離
RAD 環状高速道路へのアクセシビリティ指標
TAX $10,000 あたりの固定資産税率
PTRATIO 町ごとにみた教師 1 人あたりの生徒数
B 町ごとにみた黒人の比率を Bk としたときの (Bk−0.63)2(Bk - 0.63)^2 の値
LSTAT 給与の低い職業に従事する人口の割合

scikit-learn を用いて機械学習アルゴリズムを実装を行う際には NumPy の ndarray に変換する必要があり、現在の Pandas の DataFrame から .values の属性で取得することができます。

※ 厳密には Pandas の DataFrame でも動きますが、基本は NumPy の ndarray なので変換する癖をつけておきましょう。

また、同時にデータフレーム内の入力変数と目的変数の切り分けも行います。データフレーム内の Target が目的変数に該当し、それ以外が入力変数に該当します。目的変数以外のすべてを選択する方法はいくつかありますが、最もポピュラーな方法は drop() を使用することです。引数に設定した任意の列、行を削除できます。

  • labels : （行、列）ラベルを指定
  • axis：行方向 (axis=0) または列方向 (axis=1) を指定

今回は列を削除したいので axis=1 としましょう。

# 入力変数と目的変数の切り分け
t = df['Target'].values
x = df.drop(labels=['Target'], axis=1).values

学習用データセットとテスト用データセットへ分割

ここで、このデータセットを 2 つに分割します。それは、モデルの学習に用いるためのデータと、学習後のモデルのパフォーマンスをテストするために用いるデータは、異なるものになっている必要があるためです。これは、冒頭で少し触れた汎化性能というものに関わる重要なことです。

例え話を使ってなぜデータセットを分割する必要があるかを説明します。

例えば、大学受験の準備のために 10 年分の過去問を購入し、一部を勉強のために、一部を勉強の成果をはかるために使用したいとします。10 年分という限られた数の問題を使って、結果にある程度の信頼のおけるような方法で実力をチェックするには、下記の 2 つのうちどちらの方法がより良いでしょうか。

  • 10 年分の過去問全てを使って勉強したあと、もう一度同じ問題を使って実力をはかる
  • 5 年分の過去問だけを使って勉強し、残りの 5 年分の未だ見たことがない問題を使って実力をはかる

一度勉強した問題を再び解くことができると確認できても、大学受験の当日に未知の問題が出たときにどの程度対処できるかを事前にチェックするには不十分です。よって、後者のような方法で数限られた問題を活用する方が、本当の実力をはかるには有効でしょう。

これは機械学習におけるモデルの学習と検証でも同様に言えることです。

実力をつけるための勉強に使うデータの集まりを、学習用データセット (training dataset) といい、実力をはかるために使うデータの集まりを、テスト用データセット (test dataset) と言います。このとき、学習用データセットとテスト用データセットに含まれるデータの間には、重複がないようにします。

早速、用意した x と t を、学習用データセットとテスト用データセットに分割しましょう。どのように分けるかには色々な方法がありますが、単純に全体の何割かを学習用データセットとし、残りをテスト用データセットとする、といった分割を行う方法はホールドアウト法 (holdout method) と呼ばれます。scikit-learn では、データセットから指定された割合（もしくは個数）のデータをランダムに抽出して学習用データセットを作成し、残りをテスト用データセットとする処理を行う関数 sklearn.model_selection.train_test_split() が提供されています。

from sklearn.model_selection import train_test_split
      
# 2 つに分割
x_train, x_test, t_train, t_test = train_test_split(x, t, test_size=0.3, random_state=0)

ここで、train_test_split() の test_size という引数に 0.3 を与えています。これはテスト用データセットを全体の 30% のデータを用いて作成することを意味しています。自動的に残りの 70% は学習用データセットとなります。上のコードはデータ全体からランダムに 70% を学習用データセットとして抽出し、残った 30% をテスト用データセットとして返します。

例えば、データセット中の目標値が 1 のデータが 10 個、2 のデータが 8 個、3 のデータが 12 個 ⋯\cdots というように、カテゴリごとにまとめられて並んでいることがあります。そのとき、データセットの先頭から 18 個目のところで学習用データセットとテスト用データセットに分割したとすると、学習用データセットには目標値が 3 のデータが 1 つも含まれないこととなります。

そこで、ランダムにデータセットを分割する方法が採用されています。random_state という引数に毎回同じ整数を与えることで、乱数のシード値を固定することができ、実行の度に結果が変わることを防いでいます。こちらの値は今回は 0 としていますが、1234 でも 42 でも問題ないので好きな数字で固定しましょう。

それでは、分割後の学習用データセットを用いてモデルの学習、精度の検証を行いましょう。

モデルの学習・検証

scikit-learn を用いて学習を行う際には基本的に下記のステップに沿って実装を行います。

  • Step 1：モデルの定義
  • Step 2：モデルの学習
  • Step 3：モデルの検証

それでは実装していきましょう。

Step 1：モデルの定義

モデルの定義ではどの機械学習アルゴリズムを使用してモデルの構築を行うのかを定義します。scikit-learn で重回帰分析を行う場合は、LinearRegression クラスを使用します。sklearn.linear_model 以下にある LinearRegression クラスを読み込んで、インスタンスを作成しましょう。

# Step 1：モデルの定義
from sklearn.linear_model import LinearRegression
model = LinearRegression()

Step 2：モデルの学習

次にモデルの学習を行います。scikit-learn に用意されている手法の多くは、共通して fit() と呼ばれるメソッドを持っており、 使いやすいインターフェースとなっています。

model を用いて学習を実行するには、fit() の引数に入力値 x と目標値 t を与えます。今回学習に使用するデータセットは先程切り分けを行った、学習用のデータセットになります。

# Step 2：モデルの学習
model.fit(x_train, t_train)
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None, normalize=False)

Step 3：モデルの検証

モデルの学習が完了しました。求まったパラメータの値を確認してみましょう。重回帰分析では、重み w とバイアス b の ２ つがパラメータでした。求まった重み w の値は model.coef_ 属性に、バイアス b の値は model.intercept_ 属性に格納されています。

# 学習後のパラメータ w 
model.coef_
array([-1.21310401e-01, 4.44664254e-02, 1.13416945e-02, 2.51124642e+00, -1.62312529e+01, 3.85906801e+00, -9.98516565e-03, -1.50026956e+00, 2.42143466e-01, -1.10716124e-02, -1.01775264e+00, 6.81446545e-03, -4.86738066e-01])
# パラメータの分布をヒストグラムで可視化
plt.figure(figsize=(10, 7))
plt.bar(x=columns, height=model.coef_);
<Figure size 720x504 with 1 Axes>
# 学習後のバイアス b
model.intercept_
37.93710774183309

繰り返しますが、モデルの学習が完了したら精度の検証を行います。LinearRegression クラスは score() メソッドを提供しており、入力変数と目的変数を与えると学習済みのモデルを用いて計算した決定係数 (coefficient of determination) という指標を返します。

これは、使用するデータセットのサンプルサイズを NN、nn 個目の入力値に対する予測値を yny_n、目標値を tnt_n、そしてそのデータセット内の全ての目標値の平均値を t¯\bar t としたとき、

R2=1−∑n=1N(tn−yn)2∑n=1N(tn−t¯)2\begin{array}{c} R^2 = 1 - \frac{\sum_{n=1}^{N} (t_n - y_n)^2}{\sum_{n=1}^{N} (t_n - \bar t)^2} \end{array}

で表される指標です。

決定係数の最大値は 1 であり、値が大きいほど（1 に近いほど）モデルが与えられたデータに当てはまっていることを表します。数式の意味を最初から深く理解する必要はありませんので 1 に近ければ良い、と現時点では覚えましょう。

# Step 3：モデルの検証
print('train score : ', model.score(x_train, t_train))
print('test score : ', model.score(x_test, t_test))
train score : 0.7645451026942549 test score : 0.6733825506400171

学習済みモデルを用いて学習用データセットで計算した決定係数は、およそ 0.765 でした。テスト用データセットは およそ 0.673 と学習用データセットに対しての決定係数よりも値が小さくなりました。

教師あり学習の目的は、学習時には見たことがない新しいデータ、ここではテスト用データセットに含まれているデータに対して、またさらなる未知なデータに対しても高い性能を発揮するようにモデルを学習することです。逆に、上記の例のように学習時に用いたデータに対してはよく当てはまっていても、学習時に用いなかったデータに対しては予測値と目標値の差異が大きくなってしまう現象を、過学習 (overfitting) と言います。

推論

学習が終わったモデルに、新たな入力値を与えて予測値を計算させるには、predict() メソッドを用います。学習済みのモデルを使ったこのような計算は、推論 (inference) と呼ばれます。

本来であれば、まったく別のデータセットを与えるのですが、テスト用データセットからサンプルを 1 つ取り出し、推論を行ってみましょう。このとき、predict() メソッドに与える入力値の ndarray の形が (サンプルサイズ, 各サンプルの次元数) となっている必要があることに気をつけてください。

# 推論
y = model.predict(x_test)
      
print('予測値: ', y[0])
print('目標値: ', t_test[0])
予測値: 24.935707898576915 目標値: 22.6

線形回帰の過学習を抑制する手法

重回帰分析の結果を確認すると、学習用データセットとテスト用データセットそれぞれの結果に乖離があることに疑問を持たれたのではないでしょうか。

この乖離が大きい場合、学習したモデルは学習用データセットには適合できているが、新しいデータに対しては良い精度で予測ができていない事になります。この現象を過学習と呼びます。過学習を抑制するには下記のアプローチが有効です。

  • データセットのサンプル数を増やす
  • ハイパーパラメータを調整する
  • 他のアルゴリズムを使用する

データセットのサンプル数を増やすについては、チュートリアルの範囲で扱えるデータセットのサンプル数に限りがあるため実施できないのですが、知識として重要なので覚えておきましょう。また、ハイパーパラメータの調整については後章で詳しく紹介します。

本節では、重回帰分析以外のアルゴリズムを用いて過学習を抑制できるか確認します。どのアルゴリズムが適切かを検証するためには手法を実装し、精度を比較する方法が最も簡単です。今回は下記に挙げる 2 つのアルゴリズムを実装します。

  • Ridge 回帰（リッジ回帰）
  • Lasso 回帰（ラッソ回帰）

Ridge 回帰

過学習を抑制するための手法でよく用いられるものの 1 つに正則化 (regularization) があります。機械学習において、正則化は目的関数に正則化項と呼ばれるモデルの複雑性に罰則（ペナルティ）を科すために追加の項を導入する手法を指します。

Ridge 回帰は重回帰分析に対して重みの 2 乗で表現される L2L2 ノルムを用いて正則化を行うことで、モデルの過度な複雑さに罰則を課して過学習を抑制する手法です。

数式で表現すると、目的関数 L{\mathcal L} が下記の様に変化します。

L=∣∣Xw−t∣∣22+α∣∣w∣∣22\begin{array}{c} \mathcal L = || {\mathbf X} {\mathbf w} - {\mathbf t}||_2 ^2 + \alpha ||{\mathbf w}||_2 ^2 \end{array}

重回帰分析と比較すると第 2 項に正則化項と呼ばれる α∣∣w∣∣22\alpha ||{\mathbf w}||_2 ^2 が増えています。数式だとイメージが湧きづらいのですが、入力変数を 2 つだとすると下記の図のように可視化できます。重回帰分析の目的関数である最小二乗法と正則化項という 2 つの最小化した場所を最適解とします。

16_7

正則化項の α\alpha が罰則の強さを表し、この値を大きくするほど小さな幅でパラメータが調整されます。機械学習のモデルは学習によってパラメータの値を算出もしくは調整を行います。それに対し、この α\alpha のようにモデルの学習を実行する前に設定する値をハイパーパラメータと呼びます。

ここで新しくノルムという概念が出てきましたが、ノルムとは皆さんの馴染みやすい言葉で説明すると長さを一般化したものです。ベクトル空間で距離を与えるために数学ではノルムという概念が存在します。

入力変数の数が 2 のときのイメージ図を載せておきますので、さらに詳しく知りたい方はこちらをご覧ください。

16_4

全体の手順は重回帰分析と同じ流れです。実装を進めていきましょう。

# モデルの定義、ハイパーパラメータの値を設定
from sklearn.linear_model import Ridge
ridge = Ridge(alpha=1)
      
# モデルの学習
ridge.fit(x_train, t_train)
Ridge(alpha=1, copy_X=True, fit_intercept=True, max_iter=None, normalize=False, random_state=None, solver='auto', tol=0.001)
# モデルの検証
print('train score : ', ridge.score(x_train, t_train))
print('test score : ', ridge.score(x_test, t_test))
train score : 0.7623440182689594 test score : 0.6665819091486688

色々とハイパーパラメータ（罰則）の値も変えて試してみても、先ほどの「train : 0.765, test : 0.673」からほとんど差が無いので、今回のデータセットにおいて Ridge 回帰はあまり有効でないようです。

調整後のパラメータを確認して、どのような変化があったのか比較して確認してみます。比較する際に別々に可視化してもよいのですが、並べた方が見やすい場合もあります。そのようなときには、plt.figure.add_subplot() を使用します。

16_3

最初に plt.figure(figsize()) で大きな箱を準備し、add_subplot(2, 1, 1) とすると 2 行 1 列に区切った 1 番目（上）にプロットするということです。使用頻度が高い機能なので覚えていきましょう。今回は箱を縦に 2 行 1 列で用意します。

# 箱を準備
fig = plt.figure(figsize=(7, 10))

# 重回帰分析
ax1 = fig.add_subplot(2, 1, 1)
ax1.bar(x=columns, height=model.coef_)
ax1.set_title('Linear Regression')

# リッジ回帰
ax2 = fig.add_subplot(2, 1, 2)
ax2.bar(x=columns, height=ridge.coef_)
ax2.set_title('Ridge Regression');
<Figure size 504x720 with 2 Axes>

重回帰分析の学習済みモデルの最も大きかったパラメータと今回の最も大きいパラメータの値を比較すると正則化が適用され、Ridge 回帰の方が値が -15 から -8 に小さくなっています。正則化項の α\alpha の値を調整することによって、パラメータが調整される幅にどのように変化があるか確かめてみてください。

予測性能に大きな変化は無かったですが、重回帰分析に比べて特定の入力変数に過剰に適合していないことがわかります。

Lasso 回帰

Lasso 回帰は重回帰分析に対して、∣∣w∣∣1||{\mathbf w}||_1で表現される L1L1 ノルムを使用して正則化を行う手法です。図から入り、Ridge 回帰との違いを見つけましょう。

16_6

先ほどの Ridge 回帰と同じく重回帰分析で用いた目的関数に正則化項を加えてペナルティを課しますが、 L1L1 ノルムを使用する点が異なります。図では丸みがなくなっています。これは L1L1 ノルム、重みの絶対値を求めているためです。

L=∣∣Xw−t∣∣22+α∣∣w∣∣1\begin{array}{c} \mathcal L = ||{\mathbf X} {\mathbf w} - {\mathbf t}||_2 ^ 2 + \alpha ||{\mathbf w}||_1 \end{array}

Lasso 回帰の特徴として不要な入力変数を特定し、該当する重み ww を 00 にする事で実質的に入力変数の種類を減らすことが出来る点あります。このとき生成される重み ww のベクトル w{\mathbf w} は 00 を多く含み、これをスパース性があると表現します。

入力変数を減らすことができる特徴から、高次元の入力変数に対して高い効果を期待できます。

Lasso 回帰も今までと同様に scikit-learn で簡単に実装が可能です。

# モデルの定義
from sklearn.linear_model import Lasso
lasso = Lasso(alpha=1)
      
# モデルの学習
lasso.fit(x_train, t_train)
Lasso(alpha=1, copy_X=True, fit_intercept=True, max_iter=1000, normalize=False, positive=False, precompute=False, random_state=None, selection='cyclic', tol=0.0001, warm_start=False)
# モデルの検証
print('train score : ', lasso.score(x_train, t_train))
print('test score : ', lasso.score(x_test, t_test))
train score : 0.7084095500978869 test score : 0.6115433359595555

スコアが下がってしまいました。先程、重みをスパースにすると説明しましたが、実際にどの程度スパースになっているのか確認してみます。重みが 0 でない (lasso.coef_ != 0) 特徴量の数 (np.sum()) を出力します。

# 0 になっていない特徴量の数
print('元の特徴量の数 : ', x.shape[1])
print('Lasso の特徴量 : ', np.sum(lasso.coef_ != 0))
元の特徴量の数 : 13 Lasso の特徴量 : 10

3 つの特徴量がの重みが 0 になっています。正則化の強さを表す α\alpha の値を小さくして、変化を確認してみましょう。

# アルファを変更してみる
lasso_005 = Lasso(alpha=0.05)
lasso_005.fit(x_train, t_train)

print('train score : ', lasso_005.score(x_train, t_train))
print('test score : ', lasso_005.score(x_test, t_test))
train score : 0.7548928631432029 test score : 0.6541502573235292 
# 0 になっていない特徴量の数
print('元の特徴量の数 : ', x.shape[1])
print('Lasso005 の特徴量 : ', np.sum(lasso_005.coef_ != 0))
元の特徴量の数 : 13 Lasso005 の特徴量 : 12 
fig = plt.figure(figsize=(7, 15))

# 重回帰分析
ax1 = fig.add_subplot(3, 1, 1)
ax1.bar(x=columns, height=model.coef_)
ax1.set_title('Linear Regression')

# lasso
ax2 = fig.add_subplot(3, 1, 2)
ax2.bar(x=columns, height=lasso.coef_)
ax2.set_title('Lasso Regression 1')

# lasso_005
ax3 = fig.add_subplot(3, 1, 3)
ax3.bar(x=columns, height=lasso_005.coef_)
ax3.set_title('Lasso Regression 0.05');
<Figure size 504x1080 with 3 Axes>

更新後のパラメータを確認すると、いくつかの入力変数に紐づく重み ww が 00 になっています。冒頭にお伝えした Lasso 回帰の特長によって不要な入力変数を減らすことが出来ています。

相関関係と多重共線性問題

重回帰分析を使用しているとよく直面する問題として、多重共線性と呼ばれる問題があります。その問題と対処法を確認しましょう。

データセットの準備

こちらから regression_pls.csv をダウンロードしてください。リンクをクリックすると完了です。

データセットを Colab へアップロードして、データセットの中身を確認しましょう。Colab へのファイルのアップロードは GUI （画面操作）でも行うことが可能ですが、下記のコードを実行し、「ファイル選択」をクリックしても可能です。

from google.colab import files
uploaded = files.upload()
      
# データの確認
df = pd.read_csv('regression_pls.csv')
df.head(3)
Target x1 x2 x3 x4 x5 x6 x7 x8 x9 ... x187 x188 x189 x190 x191 x192 x193 x194 x195 x196
0 1.58 59.068 54.028 59.037114 24 0 0.213790 -0.369921 0.369921 0.213790 ... 0 0 0 0 0 0 0 0 0 0
1 1.34 46.073 40.025 46.053098 20 0 -0.001725 -0.271722 0.271722 0.001725 ... 0 0 0 0 0 0 0 0 0 0
2 1.22 60.052 56.020 60.021129 24 0 0.299685 -0.481433 0.481433 0.299685 ... 0 0 0 0 0 0 0 0 0 0

3 rows × 197 columns

df.shape
(1290, 197)

入力変数 x と目的変数 Target に切り分けてから scikit-learn で使えるデータ形式への変換を行います。今回の目的変数は Target の列が該当し、その他の列が入力変数に該当します。

# データの切り分け -> ndarray に変換
x = df.drop('Target', axis=1).values
t = df['Target'].values

print(x.shape, t.shape)
(1290, 196) (1290,)

入力変数の数が 196 個と先ほど扱ったデータセットと比べて大幅に数が増えていることが確認できます。

学習用データとテスト用データの分割

テスト用データセット全体を 30%、random_state を 0 とします。

x_train, x_test, t_train, t_test = train_test_split(x, t, test_size=0.3, random_state=0)
      
# モデルの定義
from sklearn.linear_model import LinearRegression
model = LinearRegression()
      
# モデルの学習
model.fit(x_train, t_train)
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None, normalize=False)
# モデルの検証
print('train : ', model.score(x_train, t_train))
print('test : ', model.score(x_test, t_test))
train : 0.9365471359950432 test : 0.024335695916876232

テスト用データセットに対しての決定係数が大幅に低くなっていることが確認できます。このように学習用データセットに対しての予測精度が高く、テスト用データセットの予測精度が低い状態を過学習 (Over fitting) と呼びます。

過学習が起こる原因は様々ですが、今回は入力変数同士の相関が強いものが含まれており、多重共線性 (Multicollinearity) という問題が起こっているためです。ここでは詳細な説明は割愛しますが、一般的に機械学習アルゴリズムは学習時に相関の高い入力変数を使用することにより、モデルの予測精度向上の妨げとなる場合があります。その対策としてはいくつかありますが、単純にデータセットの中から相関の高い入力変数のどちらかを削除する方法や使用するアルゴリズムを変更するなどが挙げられます。

多重共線性について詳しい解説はこちらを参照してください。

相関関係の関係

多重共線性の問題に取り組む前に、実際に入力変数同士の相関が高いのか確認してみましょう。.corr() を使用して相関関係を表す数値（相関係数）を確認することができます。

相関係数は 1 が最大で、その 2 つの変数が完全に正の相関があることを表しています。例えば変数 a と変数 b の相関係数が 1 の場合は、変数 a が 1 増えると変数 b も 1 増えるというような関係になります。

# 相関係数の算出
df_corr = df.corr()
df_corr.head()
Target x1 x2 x3 x4 x5 x6 x7 x8 x9 ... x187 x188 x189 x190 x191 x192 x193 x194 x195 x196
Target 1.000000 -0.642326 -0.648078 -0.640489 -0.524453 NaN 0.111829 -0.360696 0.357026 0.113189 ... -0.032287 -0.015204 0.019244 -0.047169 NaN 0.007788 NaN 0.002448 -0.113820 0.043600
x1 -0.642326 1.000000 0.997571 0.999978 0.908895 NaN 0.322508 -0.117193 0.134074 0.298204 ... 0.051291 0.186110 0.027947 -0.002219 NaN 0.001304 NaN -0.015226 -0.038657 0.027857
x2 -0.648078 0.997571 1.000000 0.997252 0.883891 NaN 0.322631 -0.097297 0.115794 0.294947 ... 0.047416 0.191792 0.029659 -0.007914 NaN 0.005912 NaN -0.015106 -0.062823 0.027773
x3 -0.640489 0.999978 0.997252 1.000000 0.910855 NaN 0.324352 -0.120477 0.137237 0.300415 ... 0.051542 0.186772 0.028046 -0.002001 NaN 0.001447 NaN -0.015093 -0.038138 0.028359
x4 -0.524453 0.908895 0.883891 0.910855 1.000000 NaN 0.385792 -0.284647 0.293981 0.382603 ... 0.044125 0.174983 0.022996 0.018780 NaN -0.010834 NaN -0.016378 0.027813 0.055553

5 rows × 197 columns

今回はヒートマップ（数値の大小を色の違いで表現する可視化方法）を用いて、可視化を行いましょう。可視化には seaborn を使用します。seaborn は Python の可視化用パッケージで、 Matplotlib をベースに作成されており、Matplotlib と比較してより見やすい可視化を行うことが可能です。

seaborn を用いてヒートマップを使用する際には、seaborn.heatmap() メソッドを使用します。引数のオプションで annot=True とすることにより、ヒートマップ中に数値を表示することができ、seaborn の読み込みの際は慣習的に as sns として読み込むことが多いです。

今回は入力変数の数が多いため、データセットの最初の 20 行、20 列のみを選択します。また、ノートブック上で表示するプロットのサイズを調整するために plt.figure() メソッドを使用しています。

plt.figure(figsize=(12, 8))
sns.heatmap(df_corr.iloc[:20, :20], annot=True);
<Figure size 864x576 with 2 Axes>

相関係数の値が高い部分がベージュ色に近く、低い部分が赤黒い色になっており直感的にわかるように可視化することができました。見てみると、入力変数同士の相関が高い部分があることが確認できます。例えば x1 と x16 は相関係数が 0.92 となっています。

2 つの変数の関係を可視化するには sns.jointplot() を使ってみましょう。

sns.jointplot(x='x1', y='x16', data=df);
<Figure size 432x432 with 3 Axes>

ヒートマップで確認した相関係数が高い箇所を可視化すると、線形で表現できるほどの関係が確認できました。

今回は重回帰分析を実装し、その予測精度が低かった原因として、この入力変数同士の相関が高い多重共線性の問題が起こっている可能性が高いと言えます。単純に相関の高い入力変数をデータセットから地道に取り除くアプローチも考えられますが、今回はモデル側でこの問題を解決できないか試してみます。

Partial Least Squares (PLS)

Partial Least Squares (PLS) は上記のような多重共線性の問題の対処法として有力なアルゴリズムです。中身は下記のステップで実装されています。

  1. 入力値と目標値の共分散が最大になるように主成分を抽出
  2. 抽出された主成分に対して重回帰分析を用いてモデルの学習を行う

16_5

「主成分を抽出」とありますが、こちらは教師なし学習で紹介する主成分分析と似た概念になります。主成分分析は、100 個の特徴量があった場合にその特徴量の数を分散を用いて 10 や 20 個などに削減する手法（次元削減）です。PLS との違いは教師あり学習と教師なし学習で、主成分を抽出するときに目標値の情報を使うか使わないかの違いがあります。

scikit-learn に用意されているモジュール PLSRegression を使用すると、これらのステップを意識せずに実装できますが、裏側の挙動として理解しておいてください。引数に n_components があり、ここで何次元に落とし込むかを決定します。何次元に落とすかどうかは正解がなく、試行錯誤していくことになりますが、今回は 7 として実装してみましょう。

# モデルの定義（ n_components:7 とする）
from sklearn.cross_decomposition import PLSRegression
pls = PLSRegression(n_components=7)
      
# モデルの学習
pls.fit(x_train, t_train)
PLSRegression(copy=True, max_iter=500, n_components=7, scale=True, tol=1e-06)
# モデルの検証
print('train score : ', pls.score(x_train, t_train))
print('test score : ', pls.score(x_test, t_test))
train score : 0.906376310202351 test score : 0.7387281471807299

PLS を使用することで多重共線性の解消を緩和できることがわかりました。

このように scikit-learn を用いると様々な機械学習アルゴリズムを簡単に実装できます。細かい理論も重要なのですが、最初はそれぞれのアルゴリズムの特徴を大きく捉え、どのような現場で活用できるのかといった点を抑えることを意識して学習を進めることが挫折しないコツです。

次章では分類を紹介します。

shareアイコン