斜めから機械学習！多重共線性

次々と新しい技術が登場するITの世界。それを追いかけるだけでも時間が足りないのに機械学習もとなると大変です。疑問に感じたことを深堀りしている暇はない。でも気になる。そんな機械学習の暗闇に斜めからスポットライトをあてたい。

今回のテーマは「多重共線性」です。

## 多重共線性
入力変数間に強い相関関係が存在する場合、誤差の拡大効果によって回帰係数の推定結果を信頼できなくなります。

どういうことなのか。極端な例で考えてみます。

|腹囲|身長|体重|
|-:|-:|-:|
|1|1|1|
|2|2|2|
|3|3|3|

このデータは、腹囲と身長と体重が完全に相関しています。このデータを用いて、身長と体重から腹囲を予測する線形回帰モデルを構築すると回帰係数はどうなるか。

次のいずれのモデルでも腹囲を完全に表現できます。

（腹囲） = 0.2（身長） ＋ 0.8（体重）
（腹囲） = 0.5（身長） ＋ 0.5（体重）
（腹囲） = 0.7（身長） ＋ 0.3（体重）

なので解は一つに定まりません。
身長と体重、どちらがどれだけ寄与しているのか判断できる情報がデータに含まれていないからです。

次に、このデータに僅かなノイズを加えます。
腹囲と体重の１行目にノイズを加えます。

|腹囲|身長|体重|
|-:|-:|-:|
|0.99|1|1.01|
|2|2|2|
|3|3|3|

そうすると解は一つに定まりますが ↓

（腹囲） = 2.0（身長）－ 1.0（体重）

体重が増えると腹囲が減るという不自然な結果が得られます。これが「強い相関関係が存在する場合に誤差が拡大解釈される」ということです。

実務では、多重共線性はいくつかの理由で十分起こりえるので注意が必要です。少量のサンプルしか手に入らない、データに何らかの偏りがある、説明変数が多い（次元の呪い）など。

回帰係数が極端に大きい場合に多重共線性を疑います。標準化で説明変数や目的変数のスケールを合わせておくと回帰係数が飛びぬけて大きいことが際立ち多重共線性に気づきやすくなります。


## Python で試してみる
このような説明をすると苦情？が来ることがあります。完全に相関しているデータは解が不定のハズなのに、scikit-learn の LinearRegression で試すと解が得られるからです。

```py
import numpy as np
from sklearn.linear_model import LinearRegression

y = np.array([1.,2.,3.])
X = np.array([
    [1., 1.],
    [2., 2.],
    [3., 3.],
])

lr = LinearRegression()
lr.fit(X, y)
lr.coef_
```
実行結果（回帰係数）
```
array([0.5, 0.5])
```

試しに正規方程式を用いた最小二乗法で解を計算してみます。
正規方程式
$$
X^{\mathrm T} X \boldsymbol {w} - X ^{\mathrm T}{\boldsymbol {y}} = 0\\
$$
より、回帰係数は
$$
\boldsymbol {w} = \left(X^{\mathrm T} X \right)^{-1} X ^{\mathrm T}{\boldsymbol {y}}
$$

```py
# 中心化
dev_X = X - np.mean(X, axis=0)
dev_y = y - np.mean(y, axis=0)

# 回帰係数
w = np.linalg.inv(dev_X.T.dot(dev_X)).dot(dev_X.T.dot(dev_y))
w
```
実行結果（抜粋）
```
LinAlgError: Singular matrix
```

同じデータですがエラーになります。逆行列を求められないからです。

scikit-learn の LinearRegression は、これとは違う方法で解を求めているということです。[ドキュメント](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LinearRegression.html)を読むと解を求めるソルバーに SciPy パッケージの [scipy.linalg.lstsq](https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.lstsq.html) を使っているようです。[ソースコード](https://github.com/scikit-learn/scikit-learn/blob/1.1.2/sklearn/linear_model/_base.py#L736
)もそのように実装されています。詳しくは巻末補足を参照ください。

次にノイズを加えたデータで回帰係数を求めます。

正規方程式を用いた場合
```py
# ノイズを付加したデータ
noise_y = np.array([0.99, 2., 3.])
noise_X = np.array([
    [1.01, 1.],
    [2., 2.],
    [3., 3.],
])

# 中心化
dev_X = noise_X - np.mean(noise_X, axis=0)
dev_y = noise_y - np.mean(noise_y, axis=0)

# 回帰係数
w = np.linalg.inv(dev_X.T.dot(dev_X)).dot(dev_X.T.dot(dev_y))
w
```
実行結果（回帰係数）
```
array([ 2., -1.])
```

scikit-learn の LinearRegression を用いた場合
```py
lr = LinearRegression()
lr.fit(noise_X, noise_y)
lr.coef_
```
実行結果（回帰係数）
```
array([ 2., -1.])
```

## 対処法
多重共線性の対処法の１つに正則化があります。正則化は、回帰係数が大きくなりすぎることを抑制します。ですが対処法は正則化だけではありません。次回は、他のアプローチを紹介する予定です。



## 補足
[scipy.linalg.lstsq](https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.lstsq.html) は、最小二乗問題を解くために LAPACK（Linear Algebra PACKage）を使っています。LAPACK は、線形代数演算ライブラリです。対応ライブラリに Intel MKL（Math Kernel Library）や OpenBLAS などがあります。

scipy.linalg.lstsq に `lapack_driver` というオプションがあり、デフォルトは `gelsd` となっています。詳細は、[gelsd のリファレンス](https://netlib.org/lapack/explore-html/d7/d3b/group__double_g_esolve_ga94bd4a63a6dacf523e25ff617719f752.html#ga94bd4a63a6dacf523e25ff617719f752) を参照ください。特異値分解 (SVD) を用いて最小正規解を求めていると書かれています。


使用している LAPACK ライブラリの確認方法
```py
import scipy
scipy.show_config()
```
出力例（抜粋）
```
lapack_opt_info:
    libraries = ['openblas', 'openblas']
    library_dirs = ['/usr/local/lib']
    language = c
    define_macros = [('HAVE_CBLAS', None)]
    runtime_library_dirs = ['/usr/local/lib']
```