線形回帰: ダミー変数 + 標準化の注意

scikit-learnで線形回帰をしていたら、切片がとんでもない大きさ（2e+11）になった。
どうも、特徴量に大量にあったダミー変数（業種区分のようなもの）と、数値型の標準化が悪さをしていたらしい。
それの再現と対策を書く。

線形回帰はモデルの解釈がしやすい、という点で有用な一方で、（予測単体ならともかく）まともに解釈ができる状態にもっていくまでに色々な落とし穴がある。

ダミー変数は１つで、OneHot encodingすると出てくるものとする。

方針：何をするのか

ボストン住宅価格データ + 適当なダミー変数（１変数）で、上記の問題を再現する。

コード：トラブルの再現

import pandas as pd
import numpy as np

from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import StandardScaler

from sklearn.datasets import load_boston

どれか1つでフラグが立つようなダミー変数が、フィッティングパラメータに与える影響をみたいので

1. CHAS列を削除
2. dummy変数の列を用意
3. それをOneHot enncoding

data = load_boston()
df = pd.DataFrame(data['data'], columns=data['feature_names'])
df['target'] = data['target']

df.drop('CHAS', axis=1, inplace=True)

str_list = [f'A{i}' for i in range(20)]
df['dummy'] = np.random.choice(str_list, len(df))
df = pd.get_dummies(df)

あとはダミー変数以外の入力と、ターゲットをそれぞれで標準化して、LinearRegressionでフィッティングする。

num_fe_cols = df.columns.to_list()
ex_fe = [f'dummy_{i}' for i in str_list]
ex_fe += ['target']

num_fe_cols = [fe for fe in num_fe_cols if fe not in ex_fe]

# 標準化
x_std_scaler = StandardScaler()
y_std_scaler = StandardScaler()

inputs = x_std_scaler.fit_transform(df[num_fe_cols])
targets = y_std_scaler.fit_transform(df[['target']])

# 念のため、DataFrameを分けてから入力作成（ダミー変数付き）
df_std = df.copy()
df_std.loc[:, num_fe_cols] = inputs
df_std.loc[:, 'target'] = targets

feature_names = df_std.columns.to_list()
feature_names.remove('target')

# フィッティング
lm = LinearRegression()
lm.fit(df_std.loc[:, feature_names], targets)

すると、係数はこのような結果……。ダミー変数に対応する部分の係数が"-2.87597484e+11"というとんでもない値になっている。

lm.coef_

array([[-1.23998168e-01,  1.11775353e-01,  3.20912206e-02,
        -2.20312936e-01,  2.86927225e-01, -4.96117784e-03,
        -3.37975459e-01,  3.19861003e-01, -2.37376673e-01,
        -2.38160875e-01,  9.27067463e-02, -4.29321436e-01,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11, -2.87597484e+11,
        -2.87597484e+11, -2.87597484e+11]])

切片を見ると

lm.intercept_

array([2.87597484e+11])

係数とちょうど逆符号になっている。何が起こったのだろう。

原因：どうしてこうなったのか

この現象が起きる原因は

• 目的変数は標準化されて平均が0になっている
• ダミー変数が全て0になることはない
• そのため、両辺で平均を取ったときに、$0 = 切片 + w_n$($w_n$をn番目のダミー変数に対応する係数)という関係になっている
• 切片とダミー変数用のパラメータは、逆符号かつ値が決まらないので、非常に大きくなることがある

おそらく、このような理由だと考えている。
実際、何回か計算をさせると、係数が変わったり、一見まともな値になることもあった。
切片を一意に決められないことから来る計算の不安定さが、こういった現象の原因だろう。

対策：どうすればいいのか

値が決まらないなら、切片を決めてしまえばよい、ということで切片=0にしてしまう。

lm = LinearRegression(fit_intercept=False) # デフォルトはTrue

これで解決する。

lm.coef_

array([[-0.12400473,  0.1117352 ,  0.03217197, -0.22048052,  0.2869507 ,
        -0.00498836, -0.33796915,  0.319808  , -0.23725577, -0.23817404,
         0.09274132, -0.42924623,  0.08038865,  0.0622342 ,  0.09501951,
        -0.08434313,  0.11698025,  0.02396338,  0.08164036,  0.0452252 ,
        -0.19803361, -0.10092905, -0.07053491,  0.03973364, -0.04118562,
        -0.20115119, -0.0753999 ,  0.16343603, -0.08174667,  0.22527764,
        -0.10547045, -0.18457662]])

lm.intercept_

0.0

他には

• pd.get_dummies(df, drop_first=True)にする（オールゼロのダミー変数が作られるので、上記の問題が生じない）
• 目的変数が平均がゼロにならないような変換にする（MinMaxScalerなど）

などの対策がある（が、本質的ではないと思う）。

今回のデータセットだと、ダミー変数が10以上でこの現象が起きやすくなり、それ未満だとあまり起こらないことが多かった。
エラーが出ないからといって、案外見過ごしているかもしれないので注意が必要。

多重共線性とは似て非なるというか、数値計算特有のトラブルであるといった印象を持った。