[ML] 機器學習技法：第五講 Kernel Logistic Regression

ML：基礎技法學習
Package：scikit-learn
課程：機器學習技法
簡介：第五講 Kernel Logistic Regression

SVM as Regularized Model

$$\underset{b,\mathbf{w}}{min}\frac{1}{2}{\mathbf{w}}^T\mathbf{w}+C\sum _{n=1}^{N}max(1-y_n({\mathbf{w}}^T{\mathbf{z}}_n+b),0)$$

soft-margin SVM primal
$$\begin{array}{rl}\underset{b,\mathbf{w},\mathit{\xi }}{min}& \frac{1}{2}{\mathbf{w}}^T\mathbf{w}+C\cdot \sum _{n=1}^N{\xi }_n\\ subjectto& y_n({\mathbf{w}}^T{\mathbf{x}}_n+b)\ge 1-{\xi }_{n}and{\xi }_n\ge 0foralln\end{array}$$ 當 $({\mathbf{x}}_n,y_n)$ 違反時，${\xi }_n=1-y_n({\mathbf{w}}^T{\mathbf{z}}_n+b)$
當 $({\mathbf{x}}_n,y_n)$ 合法時，${\xi }_n=0$
故可改寫為 ${\xi }_n=max(1-y_n({\mathbf{w}}^T{\mathbf{z}}_n+b),0)$
得
$$\underset{b,\mathbf{w}}{min}\frac{1}{2}{\mathbf{w}}^T\mathbf{w}+C\sum _{n=1}^{N}max(1-y_n({\mathbf{w}}^T{\mathbf{z}}_n+b),0)$$ 仔細觀察，可發現這就是 L2 regularization 的變形而已

• large margin $\iff$ fewer hyperplanes $\iff$ L2 regularization of short W
• soft margin $\iff$ special $\widehat{err}$
• larger C or C $\iff$ smaller $\lambda$ $\iff$ less regularization

regularized LogReg $\Rightarrow$ approximate SVM

• linear score $s={\mathbf{w}}^T{\mathbf{z}}_n+b$
• $er{r}_{0/1}(s,y)=[ys\le 0]$
• ${\widehat{err}}_{SVM}(s,y)=max(1-ys,0)$
• $er{r}_{SCE}(s,y)=lo{g}_2(1+e^{-ys})$

從 Linear Models for Classification 得知
若 $err\ge er{r}_{0/1}$，當縮小 $err$ 時，相當於也是在做分類
而從圖中看來
$$\begin{array}{ccccc}-\mathrm{\infty }& \leftarrow & ys& \to & +\mathrm{\infty }\\ \approx -ys& & er{r}_{SVM}(s,y)& & \approx 0\\ \approx -ys& & (ln2)\cdot er{r}_{SCE}(s,y)& & \approx 0\end{array}$$ 故 L2-regularized logistic regression 其實就相當於在做 soft margin SVM

Platt's Model

$$\begin{array}{rl}g(x)& =\theta (A\cdot ({\mathbf{w}}_{SVM}^T\mathbf{\Phi }(\mathbf{x})+b_{SVM})+B)\\ & =\theta (A\cdot (\sum _{SV}{\alpha }_n{y}_{n}K({\mathbf{x}}_n,\mathbf{x})+b_{SVM})+B)\\ & =\theta (\sum _{SV}A{\alpha }_n{y}_{n}K({\mathbf{x}}_n,\mathbf{x})+A{b}_{SVM}+B)\\ \underset{A,B}{min}& \frac{1}{N}\sum _{n=1}^{N}log(1+e^{-y_n(A\cdot (\underset{{\mathbf{\Phi }}_{SVM}({\mathbf{x}}_n)}{\underbrace{{\mathbf{w}}_{SVM}^T\mathbf{\Phi }(\mathbf{x})+b_{SVM}}})+B)})\end{array}$$ ${\mathbf{\Phi }}_{SVM}({\mathbf{x}}_n)$ 可視需求代入，像是 Kernel SVM 如上推導

該怎麼做，才能將 SVM 轉換為機率的表示呢？
簡單的想法如下圖

• 若像左邊一樣，最後再將結果用 logistic function 轉換，但這就喪失了 LogReg 的 maximum likelihood 的意義了
• 若像右邊一樣，先用 SVM 跑出初始值，再代入 LogReg，但這跟直接跑 LogReg 所得到的解是差不多的，所以 SVM 幾乎無任何作用，non-linear 也無法使用，像 kernel SVM 也無法應用在此

所以若改為先跑完 SVM，再運用 LogReg 修正一下呢？
$$g(x)=\theta (A\cdot ({\mathbf{w}}_{SVM}^T\mathbf{\Phi }(\mathbf{x})+b_{SVM})+B)$$ 看起來可行，而且其實就是利用 A 調整 ${\mathbf{w}}_{SVM}$，利用 A & B 調整 $b_{SVM}$
或是將 SVM 視為一種特別的轉換 ${\mathbf{\Phi }}_{SVM}(\mathbf{x})$，從多維轉換至一維，再進行處理
骨子裡仍可使用 SVM 的各種方法，也就是 ${\mathbf{\Phi }}_{SVM}(\mathbf{x})$ 可視需求代入，像是 Kernel SVM

• 這不是 Z 空間中最好的 LogReg 解，只是還不錯的解，需用 KLR 得最好解
• 若 ${\mathbf{w}}_{SVM}$ 夠好，通常 $A>0$
  假如 $A<0$，表示 SVM 在亂做，不然怎麼會完全反向
• 若 $b_{SVM}$ 夠好，通常 $B\approx 0$
• 與原本 SVM 的邊界不一定一樣，因為 $B$ 做了平移
• 如何解 LogReg
• GD
• SGD
• 甚至更簡單的方法，因為只有兩個變數

Kernel Logistic Regression (KLR)

Representer Theorem
claim: for any L2-regularized linear model
$$\underset{\mathbf{w}}{min}\frac{\lambda }{N}{\mathbf{w}}^T\mathbf{w}+\frac{1}{N}\sum _{n=1}^{N}err(y_n,{\mathbf{w}}^T{\mathbf{z}}_n)$$ optimal ${\mathbf{w}}_{\ast }=\sum _{n=1}^N{\beta }_n{\mathbf{z}}_n$
最佳化 ${\mathbf{w}}_{\ast }$ 必是 ${\mathbf{z}}_n$ 的線性組合

令 ${\mathbf{w}}_{\ast }={\mathbf{w}}_{\parallel }+{\mathbf{w}}_{\perp }$
${\mathbf{w}}_{\parallel }\in span({\mathbf{z}}_n)$ & ${\mathbf{w}}_{\perp }\perp span({\mathbf{z}}_n)$
希望證得 ${\mathbf{w}}_{\perp }=0$
反證法
假設存在 ${\mathbf{w}}_{\perp }$
$$\begin{array}{rl}err(y_n,{\mathbf{w}}_{\ast }^T{\mathbf{z}}_n)& =err(y_n,({\mathbf{w}}_{\parallel }+{\mathbf{w}}_{\perp }{)}^T{\mathbf{z}}_n)\\ & =err(y_n,{\mathbf{w}}_{\parallel }^T{\mathbf{z}}_n+\underset{0}{\underbrace{{\mathbf{w}}_{\perp }^T{\mathbf{z}}_n}})\\ & =err(y_n,{\mathbf{w}}_{\parallel }^T{\mathbf{z}}_n)\end{array}$$ err 的值，不因 ${\mathbf{w}}_{\perp }$ 而變
$$\begin{array}{rl}{\mathbf{w}}_{\ast }^T{\mathbf{w}}_{\ast }& =({\mathbf{w}}_{\parallel }+{\mathbf{w}}_{\perp }{)}^T({\mathbf{w}}_{\parallel }+{\mathbf{w}}_{\perp })\\ & ={\mathbf{w}}_{\parallel }^T{\mathbf{w}}_{\parallel }+2{\mathbf{w}}_{\parallel }^T{\mathbf{w}}_{\perp }+{\mathbf{w}}_{\perp }^T{\mathbf{w}}_{\perp }\\ & ={\mathbf{w}}_{\parallel }^T{\mathbf{w}}_{\parallel }+{\mathbf{w}}_{\perp }^T{\mathbf{w}}_{\perp }\\ & >{\mathbf{w}}_{\parallel }^T{\mathbf{w}}_{\parallel }\end{array}$$ 矛盾，因 err 固定，而 ${\mathbf{w}}_{\ast }$ 為最佳解，那麼在 ${\mathbf{w}}^T\mathbf{w}$ 為何 ${\mathbf{w}}_{\parallel }$ 還可以比 ${\mathbf{w}}_{\ast }$ 更小
故最佳化 ${\mathbf{w}}_{\ast }$ 必是 ${\mathbf{z}}_n$ 的線性組合

$$\begin{array}{rl}\underset{\mathit{\beta }}{min}{E}_{in}(\mathit{\beta })& =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\beta }_{m}K({\mathbf{x}}_n,{\mathbf{x}}_m)+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n\sum _{m=1}^N{\beta }_{m}K({\mathbf{x}}_m,{\mathbf{x}}_n)})\\ & =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}{\mathit{\beta }}^T\mathbf{K}\mathit{\beta }+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}})\end{array}$$ 因無限制條件，可利用 GD/SGD/... 求解 $$$$ $$\begin{array}{rl}▽E_{in}(\mathit{\beta })& =\frac{\lambda }{N}2{\mathbf{K}}^T\mathit{\beta }+\frac{1}{N}\sum _{n=1}^N\frac{-y_n{\mathbf{K}}_{(:,n)}}{1+e^{y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}}}\end{array}$$ $$$$

根據 機器學習基石：第十講 Logistic Regression & 機器學習基石：第十四講 Regularization 可得
L2-regularized logistic regression
$$\begin{array}{rl}\underset{\mathbf{w}}{min}{E}_{in}(\mathbf{w})& =\underset{\mathbf{w}}{min}\frac{\lambda }{N}{\mathbf{w}}^T\mathbf{w}+E_{in}(\mathbf{w})\\ & =\underset{\mathbf{w}}{min}\frac{\lambda }{N}{\mathbf{w}}^T\mathbf{w}+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n{\mathbf{w}}^T{\mathbf{z}}_n})\\ (\because {\mathbf{w}}_{\ast }=\sum _{n=1}^N{\beta }_n{\mathbf{z}}_n)\underset{\mathit{\beta }}{min}{E}_{in}(\mathit{\beta })& =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}(\sum _{n=1}^N{\beta }_n{\mathbf{z}}_n{)}^T\sum _{n=1}^N{\beta }_n{\mathbf{z}}_n+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n(\sum _{m=1}^N{\beta }_m{\mathbf{z}}_m{)}^T{\mathbf{z}}_n})\\ & =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\mathbf{z}}_n^T{\beta }_m{\mathbf{z}}_m+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n(\sum _{n=1}^N{\beta }_m{\mathbf{z}}_m^T{\mathbf{z}}_n)})\\ (\because {\mathbf{z}}_n^T{\mathbf{z}}_m=K({\mathbf{x}}_n,{\mathbf{x}}_m))& =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\beta }_m{\mathbf{z}}_n^T{\mathbf{z}}_m+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n(\sum _{n=1}^N{\beta }_{m}K({\mathbf{x}}_m,{\mathbf{x}}_n))})\\ & =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\beta }_{m}K({\mathbf{x}}_n,{\mathbf{x}}_m)+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n(\sum _{n=1}^N{\beta }_{m}K({\mathbf{x}}_m,{\mathbf{x}}_n))})\\ & =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\beta }_{m}K({\mathbf{x}}_n,{\mathbf{x}}_m)+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n\sum _{m=1}^N{\beta }_{m}K({\mathbf{x}}_m,{\mathbf{x}}_n)})\\ & =\underset{\mathit{\beta }}{min}\frac{\lambda }{N}{\mathit{\beta }}^T\mathbf{K}\mathit{\beta }+\frac{1}{N}\sum _{n=1}^{N}ln(1+e^{-y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}})\\ \\ ▽E_{in}(\mathit{\beta })& =\frac{\lambda }{N}2{\mathbf{K}}^T\mathit{\beta }+\frac{1}{N}\sum _{n=1}^N\frac{1}{1+e^{-y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}}}{e}^{-y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}}(-y_n{\mathbf{K}}_{(:,n)})\\ & =\frac{\lambda }{N}2{\mathbf{K}}^T\mathit{\beta }+\frac{1}{N}\sum _{n=1}^N\frac{-y_n{\mathbf{K}}_{(:,n)}}{1+e^{y_n{\mathit{\beta }}^T{\mathbf{K}}_{(:,n)}}}\end{array}$$

• 兩種不同的角度
• KLR = linear model of $\mathit{\beta }$
• $\sum _{n=1}^N\sum _{m=1}^N{\beta }_n{\beta }_{m}K({\mathbf{x}}_n,{\mathbf{x}}_m)={\mathit{\beta }}^T\mathbf{K}\mathit{\beta }$，一種特別的 regularizer
•  $\sum _{m=1}^N{\beta }_{m}K({\mathbf{x}}_m,{\mathbf{x}}_n)$ 可視同 $\mathit{\beta }$ 與 轉換過的資料 作內積
• kernel 作用在 轉換 與 regularizer
• KLR = linear model of $\mathbf{w}$
• 藏在 kernel 中，且為 L2 regularizer
• $\mathit{\beta }$
• 長度為 $N$，同資料個數
• ${\beta }_n$ 通常不為 0，與 ${\alpha }_n$ 相反

程式碼

import matplotlib.pyplot as plt
import numpy as np
 
from sklearn.svm import SVC
 
 
# kernel logistic regression
class KLR:
    def __init__(self):
        pass
    
    
    # 高斯 kernel
    def kernel_rbf(self, xn, xm):
        kv = np.exp(-self.gamma * (np.linalg.norm(xn-xm)**2))
        return kv
    
    
    # 產生 kernel matrix
    def kernelMatrix(self, x):
        m = []
        for xn in x:
            row = [self.kernel_rbf(xn, xm) for xm in x]
            m.append(row)
 
        return np.array(m)
    
    
    # error function
    # lamda/N B^T K B + 1/N log(...)
    def errFuc(self, x, y):
        sum = 0
        for n in range(self.N):
            sum += np.log(1+np.exp(-y[n] * np.dot(self.beta,self.K[:,n])))
        return self.lambda_/self.N * np.dot(self.beta, self.K).dot(self.beta) + sum/self.N
        
    
    # 梯度方程式
    def gradFuc(self, x, y):
        err = np.array([0] * self.N)
        for n in range(self.N):
            val = y[n] * np.dot(self.beta, self.K[:,n])
            temp = -y[n] * self.K[:,n] / (1 + np.exp(val))
            err = np.add(err, temp)
 
        val = self.lambda_/self.N * 2 * np.dot(self.K, self.beta) + 1/self.N * err
 
        return val
    
    
    # 訓練
    def fit(self, x, y, gamma=0.5, lambda_=0.1, n=1, max_number=1000):
        # 初始化
        y[y==0] = -1
        self.N = len(y)
        self.beta = [1] * self.N
        self.gamma = gamma   
        self.lambda_ = lambda_
        self.x = x
        self.K = self.kernelMatrix(x)
        
        t = 0
        while True:     
            g = self.gradFuc(x, y)            
            self.beta = self.beta - n*g
            err = self.errFuc(x,y)
            
            # print('{:4d} n:{:.3f} err : {}'.format(t, n, err))
            t += 1     
            if t > max_number or err < 0.01:
                break
    
    
    # 預測結果
    # w = sigma(b_n z_n)
    # g = sigma(b_n kernel(x_n,x))
    def predict(self, X):
        k = []
        for x in X:
            m = []
            for xn in self.x:
                m.append(self.kernel_rbf(xn, x))
            k.append(m)
        k = np.array(k).T
        p = np.dot(np.array(self.beta), k)
        return p
        
 
# 訓練資料
X = np.c_[(.4, -.7),
          (-1.5, -1),
          (-1.4, -.9),
          (-1.3, -1.2),
          (-1.1, -.2),
          (-1.2, -.4),
          (-.5, 1.2),
          (-1.5, 2.1),
          (1, 1),
          # --
          (1.3, .8),
          (1.2, .5),
          (.2, -2),
          (.5, -2.4),
          (.2, -2.3),
          (0, -2.7),
          (1.3, 2.1)].T
Y = np.r_[[0] * 8 + [1] * 8]
 
# C 越大表示越無法容忍錯誤
C = 1.0
 
# exp(-2|x-x'|^2)
gamma = 2
classifiers = {"Platt's Model": SVC(kernel='rbf', gamma=gamma, C=C, probability=True),
               'KLR': KLR(),
               }
# 多少分類器
n_classifiers = len(classifiers)
 
# 產生 meshgrid，共 100x100 的點
XX1, XX2 = np.mgrid[-3:3:200j, -4:4:200j]
# 打平合併
Xfull = np.c_[XX1.ravel(), XX2.ravel()]
 
# 因是 dict 所以顯示不一定照順序
for index, (name, classifier) in enumerate(classifiers.items()):
    # 訓練資料
    classifier.fit(X, Y)
    
    # 得到機率
    if name == "Platt's Model":
        # 取 1 的機率值，因接近 1 則機率越大，若取 0 會相反
        probas = classifier.predict_proba(Xfull)[:, 1]
    elif name == 'KLR':
        probas = classifier.predict(Xfull)
              
    plt.subplot(n_classifiers, 1, index + 1)
    # 調整之間的空白高度
    plt.subplots_adjust(hspace=.6)
    
    # 重新排列結果，符合 XX1
    probas = probas.reshape(XX1.shape)  
    # 畫圖，contour 不填滿顏色，contourf 填滿顏色
    im = plt.contourf(XX1, XX2, probas, cmap="winter")
    # 依圖畫出 colorbar
    plt.colorbar(im)
    
    # 畫出所有點
    plt.scatter(X[:, 0], X[:, 1], c=Y, cmap=plt.cm.Paired)
        
    # 標題
    plt.title("{}".format(name))
 
plt.show()

參考

sklearn.svm.SVC
sklearn.kernel_ridge.KernelRidge