lazy diary

統計とその周辺

XORのようなデータを非線形ロジスティック回帰で分類する

参考書は相変わらずこちら

多変量解析入門――線形から非線形へ

多変量解析入門――線形から非線形へ

今回は、基底関数展開を利用した非線形のロジスティック回帰を行ないます。基底関数を求めて、データを変換、最後にロジスティック回帰を行なうだけです。今回は基底関数としてガウスカーネルを使いました。
非線形の恩恵を確かめるために、通常の線形ロジスティック回帰では分類できないXORのようなデータを今回は使います。

import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
from sklearn.cluster import KMeans
from sklearn import linear_model
from scipy.stats import multivariate_normal

cov = 0.5 * np.identity(2)
data1 = pd.DataFrame(multivariate_normal.rvs(mean = (-2.0, -2.0), cov = cov, size = 25))
data2 = pd.DataFrame(multivariate_normal.rvs(mean = (2.0, -2.0), cov = cov, size = 25))
data3 = pd.DataFrame(multivariate_normal.rvs(mean = (2.0, 2.0), cov = cov, size = 25))
data4 = pd.DataFrame(multivariate_normal.rvs(mean = (-2.0, 2.0), cov = cov, size = 25))

data1["c"] = 0
data2["c"] = 1
data3["c"] = 0
data4["c"] = 1

data = pd.concat([data1, data2, data3, data4], join = "inner")
data = data.reset_index(drop = True)
data = data.rename(columns = {0:"x", 1:"y", "c":"c"})

x = data["x"]
y = data["y"]
clusters = data["c"]

plt.scatter(x, y, c = clusters)

f:id:mt19:20181109161138p:plain
色付けをまともにしてないので見づらいですが、たしかにXORです。ここから前回と同じようにk-meansからの基底関数の構成およびデータの変換を行ないます。

#k-meansなどでクラスタリングする
dist = []
for i in range(1,30):
  km = KMeans(n_clusters = i)
  km.fit(data[["x", "y"]])
  dist.append(km.inertia_)
  
plt.plot(range(1,30), dist, marker = ".")

f:id:mt19:20181109161217p:plain
エルボー則に従って4くらいかなと考えられます。

from sklearn import linear_model
#クラスタ数を4とする

#ガウス型基底関数
def gauss_func(x, mu, sigma):
  return np.exp((np.sum(np.power(x - mu, 2))) / (-2 * sigma))

m = 4
km = KMeans(n_clusters = m)
km.fit(data[["x", "y"]])
cluster_center = km.cluster_centers_
dist = np.power(km.transform(data[["x", "y"]]), 2)
cluster_num = km.predict(data[["x", "y"]])
least_dist = np.zeros(shape = (m, 1))
for i in range(len(data)):
  for j in range(m):
    if cluster_num[i] == j:
      least_dist[j] += dist[i][j]
for i in range(m):
  least_dist[i] = least_dist[i] / np.sum(cluster_num == i)

transformed_X = pd.DataFrame(np.zeros(shape = (len(data), m)))
for i in range(len(data)):
  for j in range(m):
    transformed_X.iloc[i, j] = gauss_func(data[["x", "y"]].iloc[i,:], cluster_center[j], least_dist[j])

さて、カーネルパラメータの推定とカーネルを使ってデータの変換まで終わったのでロジスティック回帰を実行します。

#ロジスティック回帰の実行
log_reg = linear_model.LogisticRegression(penalty = "l2", solver = "saga", C = 1000.0, max_iter = 1000).fit(transformed_X, clusters)
print(log_reg.score(transformed_X, clusters)) #->1.0

とりあえず分類できました。決定境界をプロットしてみます。

seq = np.linspace(data["x"].min(), data["x"].max(), 50)
X1, X2 = np.meshgrid(seq, seq)
Z = []
x1 = X1.reshape(-1, 1)
x2 = X2.reshape(-1, 1)
testdata = np.hstack((x1, x2))
transformed_testdata = pd.DataFrame(np.zeros(shape = (len(x1), m)))
for i in range(len(x1)):
  for j in range(m):
    transformed_testdata.iloc[i, j] = gauss_func(testdata[i,:], cluster_center[j], least_dist[j])
Z = log_reg.predict_proba(transformed_testdata)

plt.title("Non Linear Logistic Regression")
plt.xlabel("x")
plt.ylabel("y")
Z = np.array(Z[:,0]).reshape(X1.shape)
plt.contour(X1, X2, Z, colors = "black", levels = [0.5])
plt.scatter(data["x"], data["y"], c = clusters)

f:id:mt19:20181109161626p:plain
ガウスカーネルっぽい丸みを帯びた決定境界です。
基底関数になるものならばどれを使ったって構いません。(うまく分類できるかは別ですが)。例えばx^n(n = 0, 1, 2,\dots)を用いても良いです。