ラビットチャレンジ - 機械学習 第6章 サポートベクターマシン

# 0．概要
本記事は日本ディープラーニング協会認定の講座プログラムである「ラビット・チャレンジ」が提供している科目の1つである機械学習のレポートである。
記事タイトルに記載のとおり、第6章のサポートベクターマシンについて以下にまとめる。

# 1．サポートベクターマシン(SVM)
回帰または分類に適用できる教師あり学習の機械学習手法(回帰ならサポートベクター回帰)。

ここでは2クラス分類に適用する場合について記載する。
分類の場合、線形モデルの正負で2値分類する。

線形判別関数($w^{T}x+b=0$)と最も近いデータ点(サポートベクトル)との距離(マージン)が最大となる線形判別関数を求める。

平均二乗誤差($MSE$)の残差はグラフで表すと予測値からデータ点へ縦に真っ直ぐ差を取るが、マージンは線形判別関数からデータ点へ90°の線を引いた差。

決定境界に近いデータの値が変わったりすることで間違った予測されることも起きうるが、マージンを大きく取る(マージン最大化)によって間違った予測されることを防ぐことができる。

間違った予測を防ぐためにはサポートベクトルのみ考慮すれば良い(どのクラスに分類されるか微妙なサポートベクトルと違って他のデータは多少変わっても予測結果が変わることはないため)。

https://logics-of-blue.com/svm-concept/

## 1.1．目的関数

各データ点と決定境界の距離

```math
\frac{|w^{T}x+b|}{||w||}=\frac{t_{i}(w^{T}x+b)}{||w||}
```

点と直線との距離の公式
https://mathtrain.jp/tentotyokusen

マージンは決定境界と最も距離の近いデータ点との距離

```math
min_{i}=\frac{t_{i}(w^{T}x+b)}{||w||}
```

目的はマージンを最大化することなため

```math
max_{w,b}[min_{i}=\frac{t_{i}(w^{T}x+b)}{||w||}]
```

マージン上の点において、$t_{i}(w^{T}x+b)=1$が成り立つため、全ての点に対して、$t_{i}(w^{T}x+b)\geq1$が成り立つ。

```math
max_{w,b}=\frac{1}{||w||}
```

制約条件：$t_{i}(w^{T}x+b)\geq1$ $\ \ (i=1,2,...,n)$

$x$：説明変数
$w$：パラメータ
$b$：切片
$t$：分類

## 1.2．ラグランジュの未定乗数法
制約付き最適化問題を解く手法のこと。

制約$g_{i}(x)\geq0(i=1,2,...,n)$のもとで、$f(x)$が最小となる$x$は、変数$\lambda_{i}\geq0(i=1,2,...,n)$を用いて新たに定義したラグランジュ関数$L(x, \lambda)=f(x)-\sum_{i=1}^{n}{\lambda_{i}g(x)}$において$\frac{\partial L}{\partial x_{j}}=0(j=1,2,...,m)$を満たす。

https://mathtrain.jp/mlm
https://www.yunabe.jp/docs/lagrange_multiplier.html

### 1.2.1．KKT条件
制約付き最適化問題において最適解が満たす条件のこと。

制約$g_{i}(x)\geq0(i=1,2,...,n)$のもとで、$f(x)$が最小となる$x^*$は以下を満たす。

```math
\frac{\partial L}{\partial x_{j}}|x^*=0\ (j=1,2,...,m)
```

```math
g_{i}(x^*)\le0,\lambda_{i}\geq0,\lambda_{i}g_{i}(x^*)=0\ (i=1,2,...,n)
```

$L$：ラグランジュ関数
$\lambda$：ラグランジュ乗数

### 1.2.2．ラグランジュ関数

```math
L(w,b,a)=\frac{1}{2}||w||^2-\sum_{i=1}^{n}{a_{i}(t_{i}(w^{T}x+b)-1)}
```

$a_{i}\geq0\ (i=1,2,...,n)$：ラグランジュ乗数


最小となる$w,b$は以下を満たす。

```math
\frac{\partial L}{\partial w}=w-\sum_{i=1}^{n}{a_{i}t_{i}x_{i}}=0
```

```math
\frac{\partial L}{\partial b}=-\sum_{i=1}^{n}{a_{i}t_{i}}=0
```

↓

```math
w=\sum_{i=1}^{n}{a_{i}t_{i}x_{i}}
```

```math
\sum_{i=1}^{n}{a_{i}t_{i}}=0
```

■双対問題
ラグランジュ関数から$w,b$を消去したもの

https://takehr.hatenadiary.org/entry/20120608/1339138209

講義資料には載ってましたが、これは正直よく分からなかったので一旦置いておいてあとから見直す…

## 1.3．ソフトマージンSVM
データを線形分離できない時、誤差$\xi_{i}$を許容して誤差に対してペナルティを与える。

```math
\xi_{i}=1-t_{i}(w^{T}x+b)
```

線形分離できない場合でも対応でき、パラメータ$C$の大小で決定境界が変化する。
パラメータ$C$が大きいと誤差をあまり許容せず、小さいと誤差をより許容する。

### 1.3.1．目的関数

```math
\frac{1}{2}||w||^2+C\sum_{i=1}^{n}{\xi_{i}}
```

制約条件：
$t_{i}(w^{T}x+b)\geq1-\xi_{i}$ $\ \ \xi_{i}\geq0$
$0\le a_{i}\le C$

$\frac{1}{2}||w||^2$：マージンの大きさ
$C$：誤差に対するペナルティを制御するパラメータ
$\sum_{i=1}^{n}{\xi_{i}}$：誤差に対するペナルティ

## 1.4．非線形分離
線形分離できない時、特徴空間に写像してその空間で線形に分離する。

特徴空間上では線形に分離、元の空間上では非線形な分離となる。

## 1.5．カーネルトリック
非線形分離を可能とする手法の1つ。
ソフトマージンでは必ずしも良い性能となるとは限らないが、カーネルトリックは特徴ベクトルを非線形変換してその空間上で線形の識別を行うことで性能向上をしている。

カーネル関数を使い、高次元ベクトルの内積をスカラ関数で表現する。
RBFカーネル、ガウシアンカーネルを用いる。

特徴空間が高次元でも計算コストを抑えられる。

https://home.hiroshima-u.ac.jp/tkurita/lecture/svm/node6.html

# 2．ハンズオン

## 2.1．演習

### 2.1.1．ライブラリのインポート

```Python:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
```

### 2.1.2．線形カーネル

#### 2.1.2.1．線形分離可能なデータの生成

```Python:
def gen_data():
    x0 = np.random.normal(size=50).reshape(-1, 2) - 2.
    x1 = np.random.normal(size=50).reshape(-1, 2) + 2.
    X_train = np.concatenate([x0, x1])
    ys_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)
    return X_train, ys_train

X_train, ys_train = gen_data()
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
```

![liner1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/ecc14409-b30c-1c2d-a20e-3c577be17e94.png)

第5章の時と同じような感じで乱数に-2と+2を加算することで-2と+2の周辺に散らばる乱数を生成。

#### 2.1.2.2．最急降下法(バッチ勾配降下法)

```Python:
t = np.where(ys_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# 線形カーネル
K = X_train.dot(X_train.T)

eta1 = 0.01
eta2 = 0.001
n_iter = 500

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
    # 1−Ha
    grad = 1 - H.dot(a)
    # a+η_{1}(1−Ha)
    a += eta1 * grad
    # a−η_{2}(a^{T}t)t
    a -= eta2 * a.dot(t) * t
    # 0より大きければa、でなければ0
    a = np.where(a > 0, a, 0)
```

線形カーネル(線形分離)を使用する。
一部コメントで補足。

https://kenyu-life.com/2019/08/24/iris_svm_kernel/

#### 2.1.2.3．予測

```Python:
index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

# 総和 a_{i}t_{i}k(x_{i},x_{s})
term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
# t_{s} - ↑の項 の平均
b = (support_vector_t - term2).mean()

xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        # 総和 a_{i}t_{i}k(x_{i},x_{s})
        # 切片は初期化のときに掛けてる
        y_project[i] += a * sv_t * sv.dot(X_test[i])
y_pred = np.sign(y_project)

# 訓練データを可視化
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])

# マージンと決定境界を可視化
plt.quiver(0, 0, 0.1, 0.35, width=0.01, scale=1, color='pink')
```

![liner2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/8e99a93f-e23c-a4df-fe7d-84a9e801d760.png)

切片等を取得し、予測している。
一部コメントで補足。

マージンと決定境界を可視化すると分離できていることが確認できる。
丸が付いているデータがサポートベクトル。

### 2.1.3．RBFカーネル(ガウシアンカーネル)

#### 2.1.3.1．線形分離不可能なデータの生成

```Python:
factor = .2
n_samples = 50
linspace = np.linspace(0, 2 * np.pi, n_samples // 2 + 1)[:-1]
outer_circ_x = np.cos(linspace)
outer_circ_y = np.sin(linspace)
inner_circ_x = outer_circ_x * factor
inner_circ_y = outer_circ_y * factor

X = np.vstack((np.append(outer_circ_x, inner_circ_x),
               np.append(outer_circ_y, inner_circ_y))).T
y = np.hstack([np.zeros(n_samples // 2, dtype=np.intp),
               np.ones(n_samples // 2, dtype=np.intp)])
X += np.random.normal(scale=0.15, size=X.shape)
x_train = X
y_train = y

plt.scatter(x_train[:,0], x_train[:,1], c=y_train)
```

![noliner1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/d447ff74-df96-98e4-6b3f-c05f8b9e4bbc.png)

円形になる非線形の分離となるデータであることが確認できる。

#### 2.1.3.2．RBFカーネルを用いた学習

```Python:
# RBfカーネルの定義
def rbf(u, v):
    sigma = 0.8
    return np.exp(-0.5 * ((u - v)**2).sum() / sigma**2)
    
X_train = x_train
t = np.where(y_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# RBFカーネル
K = np.zeros((n_samples, n_samples))
for i in range(n_samples):
    for j in range(n_samples):
        K[i, j] = rbf(X_train[i], X_train[j])

eta1 = 0.01
eta2 = 0.001
n_iter = 5000

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
    grad = 1 - H.dot(a)
    a += eta1 * grad
    a -= eta2 * a.dot(t) * t
    a = np.where(a > 0, a, 0)
```

元空間では線形分離できないので特徴空間上で線形分離する。
RBFカーネル使っていること以外は基本的に同じ。

RBFカーネル

```math:
K(x,x')=exp(-\gamma||x-x'||^2)
```

```math:
\gamma=\frac{1}{2\sigma^2}
```

https://qiita.com/sz_dr/items/f3d6630137b184156a67
https://upura.hatenablog.com/entry/2018/01/22/003429
http://ibisforest.org/index.php?Gaussカーネル

#### 2.1.3.3．予測

```Python:
index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()

xx0, xx1 = np.meshgrid(np.linspace(-1.5, 1.5, 100), np.linspace(-1.5, 1.5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        y_project[i] += a * sv_t * rbf(X_test[i], sv)
y_pred = np.sign(y_project)

# 訓練データを可視化
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# 領域を可視化
plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])
```

![noliner2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/6ffac4be-5068-9dd7-a1ee-c21ab845934c.png)

マージンと決定境界を可視化すると非線形で分離できていることが分かる。
丸が付いているデータがサポートベクトル。

### 2.1.4．ソフトマージンSVM

#### 2.1.4.1．重なりのあるデータの生成

```Python
x0 = np.random.normal(size=50).reshape(-1, 2) - 1.
x1 = np.random.normal(size=50).reshape(-1, 2) + 1.
x_train = np.concatenate([x0, x1])
y_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)

plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
```

![soft1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/4aed73b8-bd3d-3e81-fee0-be267fc336e0.png)

これまでのデータと違い、綺麗に分離できないデータであることが分かる。

#### 2.1.4.2．マージン内部のデータやご分類データにペナルティを与える最急降下法

```Python:
X_train = x_train
t = np.where(y_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# 線形カーネル
K = X_train.dot(X_train.T)

# 誤差に対するペナルティを制御するパラメータ
C = 1
eta1 = 0.01
eta2 = 0.001
n_iter = 1000

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
    grad = 1 - H.dot(a)
    a += eta1 * grad
    a -= eta2 * a.dot(t) * t
    # 制約条件 0 ≤ a_{i} ≤ C に合うようクリッピング
    a = np.clip(a, 0, C)
```

ペナルティのパラメータ$C$は1。
np.clip(a, 0, C)で制約条件のクリッピングを行っている。

#### 2.1.4.3．予測

```Python:
index = a > 1e-8
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()

xx0, xx1 = np.meshgrid(np.linspace(-4, 4, 100), np.linspace(-4, 4, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        y_project[i] += a * sv_t * sv.dot(X_test[i])
y_pred = np.sign(y_project)

# 訓練データを可視化
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# 領域を可視化
plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])
```

![soft2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/fbe192e7-1ecb-0705-9327-365201d45cad.png)

線形カーネルやRBFカーネルを使っていた時と違ってマージン内にデータが入ることを許容して決定境界が引かれていることが分かる。

以前、CourseraのMachine Learningを受講した時、Andrew Ng教授がSVMは複雑なので自前で実装するのではなくライブラリを使った方が良いと話していたが、このハンズオンで改めて感じた(もちろん数学的な背景を理解しておくことは大切だが)。

# X．ラビットチャレンジとは
ラビットチャレンジとは、日本ディープラーニング協会認定の講座プログラムの1つ。
E資格を受験するためにはこのラビットチャレンジ等、いずれかの講座プログラムを修了しなければならない。

ラビットチャレンジの特徴は「現場で潰しが効くディープラーニング講座」の通学講座録画ビデオを編集した教材を使用した自習スタイルであるという点。
サポートは他の講座より少なく、受け身ではなく自主的に学んでいく姿勢でなければ進められないが、その分、他の講座に比べると安価であり、手が出しやすい。
ある程度知識がある人、自力で頑張るぞというガッツのある人向けではないかと感じる。

![bnr_jdla.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/389699/bce66671-a2f5-a7f9-de6f-3f87e15c253b.png)