支持向量机SVM

一、线性可分支持向量机

1、针对问题： 线性可分的二分类数据
2、思想： 寻找几何间隔最大的分离超平面，不像PLA，该超平面是唯一的。（硬间隔最大化）
$$\begin{gathered} \underset{w,b}{\max }\gamma \\ s.t.y_i\left(\frac{w}{\Vert w\Vert }\cdot x_i+\frac{b}{\Vert w\Vert }\right)\ge \gamma ,i=1,2,...,N\text{\hspace{1em}(1)} \end{gathered}$$
上式中$\gamma$表示几何间隔，由于函数间隔$\hat{\gamma }$与几何间隔$\gamma$间的关系：$\gamma =\frac{\hat{\gamma }}{\Vert w\Vert }$，可得到下式：
$$\begin{gathered} \underset{w,b}{\max }\frac{\hat{\gamma }}{\Vert w\Vert } \\ s.t.y_i\left(w\cdot x_i+b\right)\ge \hat{\gamma },i=1,2,...,N\text{\hspace{1em}(2)} \end{gathered}$$
令$\hat{w}=\frac{w}{\hat{\gamma }}$, $\hat{b}=\frac{b}{\hat{\gamma }}$，得到下式，由于$sign(w\cdot x+b)=sign(\hat{\gamma }\hat{w}\cdot x+\hat{\gamma }\hat{b})=sign(\hat{w}\cdot x+\hat{b})$，所以这样变换不会影响预测结果。
$$\begin{gathered} \underset{\hat{w},\hat{b}}{\max }\frac{1}{\Vert \hat{w}\Vert } \\ s.t.y_i\left(\hat{w}\cdot x_i+\hat{b}\right)\ge 1,i=1,2,...,N\text{\hspace{1em}(3)} \end{gathered}$$
由于$\underset{\hat{w},\hat{b}}{\max }\frac{1}{\Vert \hat{w}\Vert }$ 和$\underset{\hat{w},\hat{b}}{\min }\frac{1}{2}\Vert \hat{w}{\Vert }^2$等价，因此(3)式可等价于下式，称其为线性可分支持向量机的原始问题。
$$\begin{gathered} \underset{\hat{w},\hat{b}}{\min }\frac{1}{2}\Vert \hat{w}{\Vert }^2 \\ s.t.y_i\left(\hat{w}\cdot x_i+\hat{b}\right)\ge 1,i=1,2,...,N\text{\hspace{1em}(4)} \end{gathered}$$
求解(4)式得到$\widehat{w^{\ast }},\widehat{b^{\ast }}$，由此可得分离超平面和决策函数：
$$\widehat{w^{\ast }}\cdot x+\widehat{b^{\ast }}=0\text{\hspace{1em}(5)}$$
$$f(x)=sign(\widehat{w^{\ast }}\cdot x+\widehat{b^{\ast }})\text{\hspace{1em}(6)}$$
3、支持向量： 训练数据集的样本点中与分离超平面距离最近的样本点，即使(3)式约束条件等号成立的样本点。对于$y_i=1$的样本点，即满足$H_1$；对于$y_i=-1$的样本点，即满足$H_2$。
$$\begin{gathered} y_i({\hat{w}}^{\ast }\cdot x_i+{\hat{b}}^{\ast })=1 \\ {H}_1:{\hat{w}}^{\ast }\cdot x_i+{\hat{b}}^{\ast }=1 \\ {H}_2:{\hat{w}}^{\ast }\cdot x_i+{\hat{b}}^{\ast }=-1 \end{gathered}$$
4、间隔边界： $H_1$和$H_2$之间的距离，为$\frac{2}{\Vert {\hat{w}}^{\ast }\Vert }$。
5、学习的对偶算法： 对于线性可分支持向量机来说，可以直接求解原始问题(4)式得到结果，在这里为了后面的内容做铺垫，给出原始问题的对偶算法。（为了方便，把上面的 $\hat{w}$ 和 $\hat{b}$ 统一写为 $w$ 和 $b$）
① 构建极小化问题(4)的拉格朗日函数：
$$L(w,b,\alpha )=\frac{1}{2}\Vert w{\Vert }^2-\sum _{i=1}^N{\alpha }_i{y}_i\left(w\cdot x_i+b\right)+\sum _{i=1}^N{\alpha }_i\text{\hspace{1em}(7)}$$
则该极小化问题等价于：
$$\underset{w,b}{\min }\underset{\alpha }{\max }L(w,b,\alpha )\text{\hspace{1em}(8)}$$
根据拉格朗日对偶性，该问题的对偶问题为：
$$\underset{\alpha }{\max }\underset{w,b}{\min }L(w,b,\alpha )\text{\hspace{1em}(9)}$$
当满足KKT条件时，(8)和(9)等价。由KKT条件中的：
$$\begin{gathered} {\nabla }_{w}L(w,b,\alpha )=w-\sum _{i=1}^N{\alpha }_i{y}_i{x}_i=0 \\ {\nabla }_{b}L(w,b,\alpha )=-\sum _{i=1}^N{\alpha }_i{y}_i=0 \end{gathered}$$
得到：
$$\begin{gathered} w=\sum _{i=1}^N{\alpha }_i{y}_i{x}_i \\ \sum _{i=1}^N{\alpha }_i{y}_i=0\text{\hspace{1em}(10)} \end{gathered}$$
代入(7)，考虑对偶问题(9)，再将极大转化成求极小，得到与原始问题(4)等价的对偶最优化问题：
$$\begin{gathered} \underset{\alpha }{\min }\frac{1}{2}\sum _{i=1}^N\sum _{j=1}^N{\alpha }_i{\alpha }_j{y}_i{y}_j(x_i\cdot x_j)-\sum _{i=1}^N{\alpha }_i \\ s.t.\sum _{i=1}^N{\alpha }_i{y}_i=0 \\ {\alpha }_i\ge 0,i=1,2,......,N\text{\hspace{1em}(11)} \end{gathered}$$
求解该问题，得到最优的${\alpha }^{\ast }$，代入下式得到最优的$w^{\ast },b^{\ast }$，从而求得决策函数和分离超平面。
$$\begin{gathered} w^{\ast }=\sum _{i=1}^N{\alpha }_i^{\ast }{y}_i{x}_i \\ {b}^{\ast }=y_j-\sum _{i=1}^N{\alpha }_i^{\ast }{y}_i(x_i\cdot x_j)\text{\hspace{1em}(12)} \end{gathered}$$

二、线性支持向量机

1、针对问题： 线性不可分的二分类数据
2、思想： 允许一部分的数据点在两条分界线以内（软间隔最大化）。对每个样本点$(x_i,y_i)$引入一个松弛变量${\xi }_i\ge 0$，相应的目标函数需要支付一个代价${\xi }_i$， $C>0$为惩罚参数。得到如下凸二次规划问题（原始问题）。$w$的解是唯一的，但是$b$的解可能不唯一。
$$\begin{gathered} \underset{\hat{w},\hat{b}}{\min }\frac{1}{2}\Vert \hat{w}{\Vert }^2+C\sum _{i=1}^N{\xi }_i \\ s.t.y_i\left(\hat{w}\cdot x_i+\hat{b}\right)\ge 1-{\xi }_i,i=1,2,...,N \\ {\xi }_i\ge 0,i=1,2,...,N\text{\hspace{1em}(13)} \end{gathered}$$
该问题的对偶最优化问题为：
$$\begin{gathered} \underset{\alpha }{\min }\frac{1}{2}\sum _{i=1}^N\sum _{j=1}^N{\alpha }_i{\alpha }_j{y}_i{y}_j(x_i\cdot x_j)-\sum _{i=1}^N{\alpha }_i \\ s.t.\sum _{i=1}^N{\alpha }_i{y}_i=0 \\ 0\le {\alpha }_i\le C,i=1,2,......,N\text{\hspace{1em}(14)} \end{gathered}$$
选择合适的惩罚参数$C>0$（一般用留一法或交叉验证法对误差求和，从中选出误差最小的$C$）求解该问题，得到最优的${\alpha }^{\ast }$。代入(12)式得到最优的$w^{\ast }$；选择${\alpha }^{\ast }$的一个分量满足$0<{\alpha }_j^{\ast }<C$计算$b^{\ast }$，从而求得决策函数和分离超平面。
3. Python代码： 用交叉验证法选择合适的$C$，然后求对偶最优化问题得到最优解。

# -*- coding: utf-8 -*-
"""
Created on Sun Oct 31 15:24:29 2021

@author: 22905
"""
import numpy as np
import random
from cvxopt import matrix, solvers
import matplotlib.pyplot as plt
import csv
import time

data = np.loadtxt(open('data.csv', 'r', encoding = 'utf8'), delimiter=",",skiprows=0)
data = data.tolist() #列表
#random.shuffle(data)

x_data = np.array(data)[...,[0,1]].tolist() #如果是多维，需要增加列数，比如取前三列，括号里改为0,1,2
y_data = np.array(data)[...,[2]].tolist()

# 超参数
C_SET = [2**(i) for i in range(-15, 15, 2)]
    
def sign(v):
    if v>=0: return 1
    else: return -1
    
def find_alpha(x_train, y_train, C):
    n = len(x_train)
    X = np.array(x_train).dot(np.transpose(np.array(x_train)))
    Y = np.array(y_train).dot(np.transpose(np.array(y_train)))
    P = matrix(X*Y)
    q1 = np.array([[-1.] for i in range(n)]) #cvx只接受双精度
    q = matrix(q1)
    A = matrix(np.transpose(np.array(y_train)))
    b = matrix(np.array([[0]]))
    G1 = np.diag([-1. for i in range(n)])
    G2 = np.diag([1. for i in range(n)])
    G = matrix(np.vstack((G1,G2)))
    h1 = [[0.] for i in range(n)]
    h2 = [[C] for i in range(n)]
    h = matrix(np.vstack((h1, h2)))
    solvers.options['show_progress'] = False
    sol = solvers.qp(P, q, G, h)
    alpha = sol['x'] #n*1的
    alpha = np.array(alpha) # (n, 1)的矩阵
    return alpha

def find_weight_bias(x_train, y_train, alpha):
    weight = []
    for i in range(np.array(x_train).shape[1]):
        temp = alpha*np.array(y_train)*np.array(x_train)[...,[i]] # [...,[i]]表示索引矩阵所有行，第i列
        weight.append(temp.sum())#temp内所有元素之和
    for i in range(alpha.shape[1]):
        if alpha[i][0] > 0 and alpha[i][0] < C:
            j = i
            break
    sigmod = 0
    alpha = alpha.tolist()

    for i in range(len(x_train)):
        sigmod += y_train[i][0]*alpha[i][0]*float(np.dot(x_train[i],x_train[j]))
    bias = y_train[j][0] - sigmod
    return weight, bias

def find_false_number(weight, bias, x_val, y_val):
    false_number = 0
    for i in range(len(x_val)):
        y_predict = sign(np.array(weight).dot(np.array(x_val[i])) + bias)
        if y_predict != y_val[i][0]:
            false_number += 1
    return false_number

if __name__ == '__main__':
    
    time_start = time.time()
    false_num = [] #存放错误数量
    interval = len(data)//5 #索引间隔
    # 交叉验证
    for C in C_SET:
        num = 0
        for i in range(5):
            start, end = i*interval, (i+1)*interval
            x_val = x_data[start: end]
            y_val = y_data[start: end]
            x_train = x_data[0:start] + x_data[end:len(x_data)]
            y_train = y_data[0:start] + y_data[end:len(y_data)]
            alpha = find_alpha(x_train, y_train, C)
            weight, bias = find_weight_bias(x_train, y_train, alpha)
            num += find_false_number(weight, bias, x_val, y_val)
        false_num.append(num)
        
    fine_C = C_SET[false_num.index(min(false_num))]
    print("the fine C is:", fine_C)
    fine_alpha = find_alpha(x_data, y_data, fine_C)
    fine_weight, fine_bias = find_weight_bias(x_data, y_data, fine_alpha)
    fine_false_number = find_false_number(fine_weight, fine_bias, x_data, y_data)
    print("right rate is:", (len(x_data)-fine_false_number)/len(x_data))
    print("the fine weight and bias are:", fine_weight, bias)
    time_end = time.time()
    print('running time is:', time_end - time_start)
    
    for i in range(len(x_data)):
        if y_data[i][0] == 1:
            plt.plot(x_data[i][0], x_data[i][1], 'ro')
        else:
            plt.plot(x_data[i][0], x_data[i][1], 'bo')
    # 二维空间中的线
    X = np.arange(round(min(np.array(x_data)[...,[0]])[0])-2, round(max(np.array(x_data)[...,[0]])[0])+2)
    Y = (-fine_weight[0]*X-fine_bias)/fine_weight[1]
    plt.plot(X, Y)
    plt.show()

分类结果：

三、非线性支持向量机

对于线性不可分数据，可以寻找一个正定核函数$K(x_i,x_j)$，将原数据映射到另一个特征空间中使得数据线性可分。可将(14)写为：
$$\begin{gathered} \underset{\alpha }{\min }\frac{1}{2}\sum _{i=1}^N\sum _{j=1}^N{\alpha }_i{\alpha }_j{y}_i{y}_{j}K(x_i,x_j)-\sum _{i=1}^N{\alpha }_i \\ s.t.\sum _{i=1}^N{\alpha }_i{y}_i=0 \\ 0\le {\alpha }_i\le C,i=1,2,......,N\text{\hspace{1em}(15)} \end{gathered}$$
求解该问题，得到最优的${\alpha }^{\ast }$。代入(12)式得到最优的$w^{\ast }$；选择${\alpha }^{\ast }$的一个分量满足$0<{\alpha }_j^{\ast }<C$计算$b^{\ast }$，从而求得决策函数和分离超平面。
$$b^{\ast }=y_j-\sum _{i=1}^N{\alpha }_i^{\ast }{y}_{i}K(x_i,x_j)\text{\hspace{1em}(16)}$$
构造决策函数：
$$f(x)=sign\left(\sum _{i=1}^N{\alpha }_i^{\ast }{y}_{i}K(x,x_i)+b^{\ast }\right)\text{\hspace{1em}(17)}$$