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

はじめに

　ロジスティック回帰はいうまでもなく線形分類器です。なので、非線形の分類問題は本来解けません。

ロジスティック回帰が線形分離不可能な分類問題を解けないことの説明 - 静かなる名辞

　しかし、特徴量を非線形変換したり、交互作用項を入れたりして使えば、非線形の分類問題にも十分使えます。

参考：
交互作用項を入れればロジスティック回帰でも非線形分離可能になることもある - 六本木で働くデータサイエンティストのブログ

　どれくらいの威力があるのでしょうか？　やってみましょう。

準備

　便利なmain関数を作っておきましょう。

def main(X, y, model, figname):
    model.fit(X, y)
    
    cm_bright = ListedColormap(['#FF0000', '#0000FF'])
    plt.scatter(X[:, 0], X[:, 1], c=y, 
                cmap=cm_bright, edgecolors='k')

    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.1
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                         np.arange(y_min, y_max, h))
    Z = model.predict_proba(np.c_[xx.ravel(), yy.ravel()])[:,1]
    Z = Z.reshape(xx.shape)

    cm = plt.cm.RdBu
    plt.contourf(xx, yy, Z, cmap=cm, alpha=.2)
    plt.savefig(figname)

　一応下でもコード全体を示しながら説明しますが、基本的にはこの中身は変えないでやっていく予定です。

多項式は面白い

　さて、データの非線形変換といえば、多項式変換でしょう。sklearnではPolynomialFeaturesが使えます。

sklearn.preprocessing.PolynomialFeatures — scikit-learn 0.21.3 documentation

　こいつは面白くて、ドキュメントにはこういう記述があります。

For example, if an input sample is two dimensional and of the form [a, b], the degree-2 polynomial features are [1, a, b, a^2, ab, b^2]

　abも入るんだ。知らなかった。いわゆる交互作用項ですね。

　詳しい使い方の記事はこちらです。

scikit-learnのPolynomialFeaturesで多項式と交互作用項の特徴量を作る - 静かなる名辞

moons

　三日月型のグループが2つあるようなデータがmoonsです。moonsくらいならすぐできます。

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from sklearn.datasets import make_moons
from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LogisticRegression
from sklearn.pipeline import Pipeline

def main(X, y, model, figname):
    model.fit(X, y)
    
    cm_bright = ListedColormap(['#FF0000', '#0000FF'])
    plt.scatter(X[:, 0], X[:, 1], c=y, 
                cmap=cm_bright, edgecolors='k')

    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.1
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                         np.arange(y_min, y_max, h))
    Z = model.predict_proba(np.c_[xx.ravel(), yy.ravel()])[:,1]
    Z = Z.reshape(xx.shape)

    cm = plt.cm.RdBu
    plt.contourf(xx, yy, Z, cmap=cm, alpha=.2)
    plt.savefig(figname)

if __name__ == "__main__":
    X, y = make_moons(noise=0.3, random_state=0)

    pf = PolynomialFeatures(degree=4, include_bias=False)
    lr = LogisticRegression(solver="lbfgs")
    model = Pipeline([("pf", pf), ("lr", lr)])    

    main(X, y, model, "fig1.png")

　4次までです。なんとなく怪しいし、よく見ると対称になっているべき部分がぜんぜん対称じゃないなど微妙な部分もあるのですが、それでもできます。

XOR

　2次ですでにabが入るので、できて当然。

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LogisticRegression
from sklearn.pipeline import Pipeline

def main(X, y, model, figname):
    model.fit(X, y)
    
    cm_bright = ListedColormap(['#FF0000', '#0000FF'])
    plt.scatter(X[:, 0], X[:, 1], c=y, 
                cmap=cm_bright, edgecolors='k')

    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.1
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                         np.arange(y_min, y_max, h))
    Z = model.predict_proba(np.c_[xx.ravel(), yy.ravel()])[:,1]
    Z = Z.reshape(xx.shape)

    cm = plt.cm.RdBu
    plt.contourf(xx, yy, Z, cmap=cm, alpha=.2)
    plt.savefig(figname)

def make_xor():
    np.random.seed(0)
    x = np.random.uniform(-1, 1, 300)
    y = np.random.uniform(-1, 1, 300)
    target = np.logical_xor(x > 0, y > 0)
    return np.c_[x, y], target

if __name__ == "__main__":
    X, y = make_xor()

    pf = PolynomialFeatures(degree=5, include_bias=False)
    lr = LogisticRegression(solver="lbfgs")
    model = Pipeline([("pf", pf), ("lr", lr)])    

    main(X, y, model, "fig2.png")

　悪くなさそう。しれっと5次でやっているのは、この方が表現力が若干高いような気がするからです（厳密に検証していません。実用的には2次の方が良いかも。次数が少ないに越したことはないので）。

circles

　円形のデータも実はきれいに分類できます。しかもたった2次で。

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from sklearn.datasets import make_circles
from sklearn.preprocessing import PolynomialFeatures
from sklearn.linear_model import LogisticRegression
from sklearn.pipeline import Pipeline

def main(X, y, model, figname):
    model.fit(X, y)
    
    cm_bright = ListedColormap(['#FF0000', '#0000FF'])
    plt.scatter(X[:, 0], X[:, 1], c=y, 
                cmap=cm_bright, edgecolors='k')

    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.1
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                         np.arange(y_min, y_max, h))
    Z = model.predict_proba(np.c_[xx.ravel(), yy.ravel()])[:,1]
    Z = Z.reshape(xx.shape)

    cm = plt.cm.RdBu
    plt.contourf(xx, yy, Z, cmap=cm, alpha=.2)
    plt.savefig(figname)

if __name__ == "__main__":
    X, y = make_circles(noise=0.2, factor=0.5, random_state=1)

    pf = PolynomialFeatures(degree=2, include_bias=False)
    lr = LogisticRegression(solver="lbfgs")
    model = Pipeline([("pf", pf), ("lr", lr)])    

    main(X, y, model, "fig3.png")

　これについてはにわかには信じがたかったので、決定木で可視化してみました。

決定木をいろいろな方法で可視化する - 静かなる名辞

from sklearn.datasets import make_circles
from sklearn.preprocessing import PolynomialFeatures
from sklearn.tree import DecisionTreeClassifier
from dtreeviz.trees import dtreeviz

def main():
    X, y = make_circles(noise=0.2, factor=0.5, random_state=1)
    pf = PolynomialFeatures(degree=2, include_bias=False)
    X_pf = pf.fit_transform(X)
    feature_names = ["x", "y", "x^2", "xy", "y^2"]

    dtc = DecisionTreeClassifier(max_depth=4)
    dtc.fit(X_pf, y)

    viz = dtreeviz(dtc, X_train=X_pf, y_train=y, 
                   feature_names=feature_names, target_name="",
                   class_names=["outer", "center"])
    viz.save("dtc.svg")

if __name__ == "__main__":
    main()

# rsvg-convert dtc.svg --format=png --output=dtc.png

　xとyの二乗が大きければ外側、小さければ内側って感じ。まあ、そりゃそうか。

　追記：
　これは二次判別分析をやったのとだいたい同じことになるんだっけ……と思って英語版wikiを見たら、そんな感じの記述があったのでたぶんそういうこと。

Quadratic classifier - Wikipedia

　違う解釈としては、円の方程式を思い出してください。半径一定の円で切ることができるので、二次式で十分です。

まとめ

　このように、ロジスティック回帰でもデータを多項式変換することで非線形の分類問題を解くことができます。

　sklearnだと、PolynomialFeaturesは気楽に使えていいですね。ただ、次数が増えると急速に変数の数が増えていくので、無思慮には使えませんが。実用でいけるかどうかは、特徴量の大きさが現実的な範囲に収まるかどうかにかかっています。3次以上は実際問題としては厳しいので、それなりに気を使う必要があります。
（経験的には、元データが高次元であればあるほど分類器の非線形への適応度合いは問われなくなるような気がするので、なんとかなる可能性が高いです。）

　低次元（<100くらい）の非線形問題であれば、だいたいこれ＋線形分類器で解けそうな雰囲気があります。