scikit-learnのPolynomialFeaturesで多項式と交互作用項の特徴量を作る

はじめに

　回帰などで非線形の効果を取り扱いたいとき、多項式回帰は定番の方法です。また、交互作用項も使うと有用なときがあります。

　pythonユーザはいきなりSVRやランダムフォレスト回帰などの非線形回帰を使うことが多い気もしますが、線形モデルでも特徴量を非線形変換すればできます。scikit-learnでやるのであれば、PolynomialFeaturesを使います*1。

使い方

　リファレンスはここです。

sklearn.preprocessing.PolynomialFeatures — scikit-learn 0.21.3 documentation

　以下のような引数を取ります。

class sklearn.preprocessing.PolynomialFeatures(degree=2, interaction_only=False, include_bias=True, order=’C’)

degree

　多項式の次数です。議論の余地なし。

interaction_only

　Trueにすると交互作用項だけ出力します。交互作用項というのは、つまりX[:,0]*X[:,1]のようなものです。

include_bias

　デフォルトではTrueで、バイアス項（=すべて1の列）を入れてくれます。sklearnの他のモデルで切片を推定してくれないものは少ないので、ほとんどのケースではFalseにしても差し支えないでしょう。

order

　出力の列の順番が変わります。普通は意識しなくてもいいでしょう。

　さて、面白そうなのはinteraction_onlyです。ちょっと設定を変えて動作をチェックしてみましょう。

　interaction_only=False (default)の場合

>>> from sklearn.preprocessing import PolynomialFeatures
>>> import numpy as np
>>> a = np.arange(6).reshape(3, 2)
>>> pf = PolynomialFeatures(degree=3, include_bias=False)
>>> pf.fit_transform(a)
array([[  0.,   1.,   0.,   0.,   1.,   0.,   0.,   0.,   1.],
       [  2.,   3.,   4.,   6.,   9.,   8.,  12.,  18.,  27.],
       [  4.,   5.,  16.,  20.,  25.,  64.,  80., 100., 125.]])
# 簡略のためX[:,0]をa, X[:,1]をbと表記すると、
# a, b, a**2, a*b, b**2, a**3, a**2 * b, a * b**2, b**3
# つまり1次以上3次以下の項がすべて出てくる

　interaction_only=Trueの場合

>>> pf = PolynomialFeatures(degree=3, interaction_only=True, \
...                         include_bias=False)
>>> pf.fit_transform(a)
array([[ 0.,  1.,  0.],
       [ 2.,  3.,  6.],
       [ 4.,  5., 20.]])
# 上と同様に表記すると、
# a, b, a*b
# のみが現れる

　ドキュメントに書いてある

If true, only interaction features are produced: features that are products of at most degree distinct input features (so not x[1] ** 2, x[0] * x[2] ** 3, etc.).

　という記述はこのことを示しています。

　交互作用項を使わない、というオプションはないようです。ちょっと残念かも。

試してみる

　3次多項式で表せるデータを生成し、多項式回帰してみます。

import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LinearRegression
from sklearn.pipeline import Pipeline

def main():
    np.random.seed(0)
    x = np.linspace(-3, 3, 40)
    y = 2 + 2 * x + 2 * x ** 2 - 2 * x ** 3 +\
        np.random.normal(scale=10, size=x.shape)
    X = x.reshape(-1, 1)

    pf = PolynomialFeatures(degree=3, include_bias=False)
    lr = LinearRegression()
    pl = Pipeline([("PF", pf), ("LR", lr)])
    pl.fit(X, y)

    pred_y = pl.predict(X)
    plt.scatter(x, y, c="b", alpha=0.2)
    plt.plot(x, pred_y, c="b")
    plt.title(    
        "predicted model: "
        "{0:.2f} + {1:.2f}x + {2:.2f}x^2 + "
        "{3:.2f}x^3".format(
            pl.named_steps.LR.intercept_,
            *pl.named_steps.LR.coef_))
    plt.savefig("result.png")

if __name__ == "__main__":
    main()

　ノイズが大きめなせいもありますが、あまり係数はあてにならないことがわかります。でも予測線自体は（この区間なら）いい感じで引けているので、内挿的な用途ならこれで十分です。外挿するのにはあまり向きません。

結論

　多項式はけっこう強力なので、使えるときは使うと良いと思います。今回は回帰で試しましたが、分類でも使えるときがあります。

非線形がなんだ！ロジスティック回帰＋多項式でやってやる！ - 静かなる名辞

　ただし交互作用項が出てくるので、次数が増えると組み合わせ爆発的に特徴量次元数が大きくなるのが厄介なところです。逆に、変数のn乗にあまり価値を見出さなければ、interaction_onlyで使うという手もあります。あるいは、2次にして変数のn乗+2次の交互作用だけというのも使いやすいと思います。

*1:なお、タイトルの「多項式の特徴量」という表現は自分でも気に入らないと思っていますが、他に言いようがないのでこう書いています。いい案があったら教えてください