LIBSVM使用：二维数据SVM支持向量、决策边界与分类间隔可视化（动图）

| 图源

  LIBSVM是台湾大学林智仁教授等开发的一个简单易用、快速有效的SVM、SVR的开源软件包，有各种语言接口，包括python、matlab和java等。现在有很多的机器学习软件包的SVM，SVR都是基于LIBSVM开发的。软件包里超参优化方法是网格法，即按照超参的范围和步长，遍历所有超参的组合可能，最后在训练集上平均准确率最高的那一组参数作为最佳参数。当然，本文不涉及寻优，只是涉及LIBSVM数据转换，模型训练，返回值解析和可视化。这是笔者上学期模式识别课程的一个作业，写了很详细注释，也算是一个存货，放在这里，供大家参考。

  博客发布时间：每周二上午10点

copyright© 意疏：https://blog.csdn.net/sinat_35907936/article/details/115110787

---

LibSVM

---

• 数据格式

  LIBSVM的数据格式如下，一行是一个样本，行首是标签（类别编码，二分类时是-1或者1），然后是特征编码：特征的值一组一组的依次排列，组与组之间用空格分开。这样的数据格式的好处是，某些值为0的特征在数据储存时可以省略，因为有特征编号的存在，读数据的时候知道哪些值应该补全。

label < space > 1: feature 1 < space > 2:Feature 2 < space >...... < space > n: feature n

• 模型参数解析

  SVM，支持向量机，它最本质的含义其实应该是最优分类超平面，表达式可以由式（1）给出，其中$W^{\ast }$和$b^{\ast }$为其参数。由于在优化完成之后，最优分类超平面（决策面）可以由少量向量完全支持住（确定），这些向量就被称为支持向量，所以这个算法才叫支持向量机。

  支持向量也不是什么新的东西，它就是输入样本特征组成的向量$X$本身，只是要满足一个关系，二分类且不含松弛因子时的支持向量满足的关系由式（2）给出，$y_i,X_i$为第$i$个样本的标签和特征组成的向量。

$$W^{\ast }X+b^{\ast }=0\text{\hspace{1em}(1)}$$

$$y_i(W^{\ast }{X}_i+b^{\ast })-1=0\text{\hspace{1em}(2)}$$

$$W^{\ast }=\sum _{i=1}^N{\alpha }_i^{\ast }{y}_i{X}_i(\alpha =0,X_i\notin SV)\text{\hspace{1em}(3)}$$

  我们要可视化支持向量、决策面，由式（1）确定、和分类间隔，由式（2）得出，我们就先要获取到支持向量SV，决策面的参数$W^{\ast }$和$b^{\ast }$。

  $W^{\ast }$由式（3）确定，${\alpha }_i$为优化时条件约束的拉格朗日乘子，当样本$X_i$不在SV中时，${\alpha }_i=0$,当样本$X_i$在SV中时，${\alpha }_i$ 不为0。

  相关参数都存在model中，我们只需要调用以下接口即可获得。

  SV=model.get_ SV()

  $y{\alpha }^{\ast }$=model.get_ sv_ coef()

  $b^{\ast }$=model.rho

  sv_idx = model.get_sv_indices()

  注：sv_indices是指属于支持向量的样本在训练集中的序号，用于找到该支持向量的标签，在可视化时需要根据标签配色。

---

二维数据点产生与格式转换

---

• 二维线性可分与不可分数据点产生

  用sklearn.datasets的make_blobs工具分别产生120个线性可分与120个线性不可分的二维数据点，来作为我们的训练数据，是否线性可分由方差std控制。

import sklearn.datasets
import matplotlib.pyplot as plt


def tod_data(n, std):
    X, Y = sklearn.datasets.make_blobs(n_samples=n, n_features=2, centers=2, cluster_std=std, random_state=1)
    return X, Y

if __name__ == '__main__':
    X1, Y1 = tod_data(120, 2)
    X2, Y2 = tod_data(120, 4)
    plt.figure(0)
    plt.scatter(X1[:, 0], X1[:, 1], marker='h', c=Y1)
    plt.figure(1)
    plt.scatter(X2[:, 0], X2[:, 1], marker='h', c=Y2)

    plt.show()
    print(X2.shape, Y2.shape)
    print(X1.shape, Y1.shape)
    # (120, 2) (120,)
    # (120, 2) (120,)

• 数组转LIBSVM格式

  上述过程产生的数据中，数据是一个二维数组，标签是一个一维向量，标签的值为0或1。我们把它按照LIBSVM的数据格式要求将它们转换到一个文件中，便于使用LIBSVM的svm_read_problem()函数进行读取。

import os

# X是一个二维数组，Y是一个一维向量，测试集可以随意给一个值如x

# libsvm数据整理
def libsvm_data_arrange(filename, X, Y):
    temp_out = []

    if os.path.exists('libsvm_data') == 0:
        os.mkdir('libsvm_data')
    flag = 1  # 行首标签标志
    with open('libsvm_data//' + filename, mode='w+') as f:
        for i in range(len(X)):
            for j in range(len(X[0])):
                if flag:
                    if Y[i] == 0:
                        Y[i] = -1
                    temp_out.append(str(Y[i]))
                    temp_out.append(str(j + 1) + ':' + str('%.6f' % X[i][j]))
                    flag = 0
                else:
                    temp_out.append(str(j + 1) + ':' + str('%.6f' % X[i][j]))
            # print(temp_out)
            flag = 1
            f.write(' '.join(temp_out) + '\n')
            temp_out = []
    f.close()

  转换后的数据如下，我们是做二分类，所以标签是-1或者1，make_blobs产生的数据标签是0或1所以代码中也做了修改。二维数据只有两个特征，即x轴坐标和y轴坐标，他们也是按照特征编码：特征的值的方式存储的。

---

可视化

---

• 训练

  我们需要通过训练来找最优分类超平面，调用LIBSVM的svm_train()即可开始训练，非常简单，返回值是SVM的model。model在保存之后其实就是一个文本，里面包含一些参数和所有的支持向量（远小于样本数），如下图。

model1 = svm_train(y1, x1, '-t 0 -c ' + str(c))

   -t： kernel_type，核函数类型

      0：线性函数，相当于没有变换
      2：RBF(default)，最常用

   -c： cost，惩罚系数或者松弛因子。

      1：default

  核函数技术，对输入数据进行非线性变换的一种简化手段，目的是希望原线性不可分的数据，在变换后尽量线性可分。在转换后的数据空间，决策面是线性的，投影到原始数据空间决策面就变成了非线性。非线性的决策面和分类间隔不好画，所以只可视化了线性的，即-t 0。

  松弛因子C为线性不可分的样本设计的超参数，它的存在，允许了错误分类和分类间隔内样本的出现，分类间隔是指式（3）确定的两条直线中间的部分。松弛因子的大小用于控制这两种样本的数量。C越大分类间隔越小，这两种样本的数量越小，越容易过拟合。C越小分类间隔越大，这两种样本的数量越多，越容易欠拟合，实际应用中需要用寻优方法找到合适的C值。笔者采用了多个C值来可视化，便于直观的感受C与分类间隔之间的关系。

  由于分类间隔内的样本点，也是支持向量，故，含有松弛因子时支持向量的条件：

   边界上的SV：Margin SV
$$y_i(W^{\ast }{X}_i+b^{\ast })-1=0\text{\hspace{1em}(4)}$$

   分类间隔内的SV：Other SV，$(|{\alpha }_i|=C)$
$$y_i(W^{\ast }{X}_i+b^{\ast })-1<0\text{\hspace{1em}(5)}$$

• 可视化

   训练完成之后依次取出需要的值，然后用matplotlib.pyplot可视化。每一个C值都对应一张图片，用ScreenToGif可以将其转换成动图。图中边界上的SV用□标识，分类间隔里的SV用△标识。

   分类间隔：

$$dis=\frac{2}{||W^{\ast }|{|}_2}\text{\hspace{1em}(6)}$$

import os
import sklearn.datasets
from svmutil import *
import numpy as np
import math
import matplotlib.pyplot as plt


# 数据点产生
def tod_data(n, std):
    X, Y = sklearn.datasets.make_blobs(n_samples=n, n_features=2, centers=2, cluster_std=std, random_state=1)
    return X, Y


# X是一个二维数组，Y是一个一维向量，测试集没有Y，可以随意赋一个值如'x'

# libsvm数据整理
def libsvm_data_arrange(filename, X, Y):
    temp_out = []

    if os.path.exists('libsvm_data') == 0:
        os.mkdir('libsvm_data')
    flag = 1
    with open('libsvm_data//' + filename, mode='w+') as f:
        for i in range(len(X)):
            for j in range(len(X[0])):
                if flag:
                    if Y[i] == 0:
                        Y[i] = -1
                    temp_out.append(str(Y[i]))
                    temp_out.append(str(j + 1) + ':' + str('%.6f' % X[i][j]))
                    flag = 0
                else:
                    temp_out.append(str(j + 1) + ':' + str('%.6f' % X[i][j]))
            # print(temp_out)
            flag = 1
            f.write(' '.join(temp_out) + '\n')
            temp_out = []
    f.close()


graph_title = ''
model_name = ''
x1 = []
y1 = []
X1 = []
Y1 = []

Linearly_separable_data = True
Linearly_separable_data = False

# 选择线性可分数据，和线性不可分数据
if Linearly_separable_data:
    X1, Y1 = tod_data(120, 2)

    print(X1, Y1)

    # 将产生的数据点转换成libsvm格式
    libsvm_data_arrange('LSD_data', X1, Y1)
    graph_title = 'Linearly separable data'
    model_name = 'LSD_model'
    # 读数据
    y1, x1 = svm_read_problem(r'libsvm_data//LSD_data')
else:
    X1, Y1 = tod_data(120, 4)
    libsvm_data_arrange('LID_data', X1, Y1)
    graph_title = 'Linearly inseparable data'
    model_name = 'LID_model'
    y1, x1 = svm_read_problem(r'libsvm_data//LID_data')

# 惩罚因子c
C = [1, 0.5, 0.2, 0.1, 0.08, 0.06, 0.04, 0.02, 0.008, 0.005, 0.003]

for n in range(len(C)):

    c = C[n]
    # 训练，选线性模型
    model1 = svm_train(y1, x1, '-t 0 -c ' + str(c))
    print(model1)
    sv_temp = []

    # 创建绘图窗口
    plt.figure(2)
    ax = plt.subplot(1, 1, 1)
    ax.set_xlabel('x1')
    ax.set_ylabel('x2')
    ax.set_title(graph_title)

    sv_x1 = []
    sv_y1 = []
    # 找支持向量SV的坐标
    sv1 = model1.get_SV()
    # 找支持向量SV在训练集中的序号，从1开始，注意与列表序号区别，需要减一，y2[sv_idx[i]-1]
    sv_idx1 = model1.get_sv_indices()
    print(sv_idx1)

    for i in range(len(sv_idx1)):
        sv_y1.append(int(y1[sv_idx1[i] - 1]))  # 通过序号找到标签
        for j in range(len(sv1[0])):
            sv_temp.append(sv1[i][j + 1])  # 特征序号从1，开始，故需j + 1
        # print(x1[sv_idx1[i]-1])
        sv_x1.append(sv_temp)
        sv_temp = []
    sv_x1 = np.array(sv_x1).reshape(i + 1, j + 1)
    # print(sv1)
    # print(sv_y1)

    # 系数alpha*y
    alpha1 = model1.get_sv_coef()

    # 求决策面系数w，由alpha和支持向量坐标得到，画支持向量

    w = 0.0
    for i in range(len(alpha1)):
        # w = w + abs(alpha1[i][0]) * sv_y1[i] * sv_x1[i]  # 求w的公式
        w = w + alpha1[i][0] * sv_x1[i]  # 求w的公式
        # print(alpha1[i][0] * sv_y1[i])

        # 画支持向量散点图
        if abs(alpha1[i][0]) < C[n]:
            plt.scatter(sv_x1[i, 0], sv_x1[i, 1], marker='s', color="none", linewidths=1, s=100, edgecolors='red')
        else:
            plt.scatter(sv_x1[i, 0], sv_x1[i, 1], marker="^", color="none", linewidths=1, s=100, edgecolors='blue')

    # 决策面偏置b
    b = -model1.rho[0]

    # 画数据散点图
    ax.scatter(X1[:, 0], X1[:, 1], marker='o', s=20, c=Y1)
    ax.scatter(sv_x1[:, 0], sv_x1[:, 1], marker='x', s=20, c=sv_y1)

    # 画决策边界与决策间隔

    xx = 15
    yy = 15
    dis_x = np.arange(-xx, xx, 0.1)
    dis_y = np.arange(-yy, yy, 0.1)

    dis_x, dis_y = np.meshgrid(dis_x, dis_y)
    z = w[0] * dis_x + w[1] * dis_y + b
    z1 = w[0] * dis_x + w[1] * dis_y + b + 1
    z2 = w[0] * dis_x + w[1] * dis_y + b - 1

    ax.contour(dis_x, dis_y, z, 0, colors='red')
    ax.contour(dis_x, dis_y, z1, 0)
    ax.contour(dis_x, dis_y, z2, 0)

    # 算分类间隔
    dis1 = 2 / math.sqrt(w[0] * w[0] + w[1] * w[1])
    #  画间隔线段

    # 在下边那条线上找一个点，令dis_x = 0
    margin1_y1 = (1 - b) / w[1]
    theta = math.atan(abs(w[1] / w[0]))
    margin1_x2 = 0 + dis1 * math.cos(theta)
    margin1_y2 = margin1_y1 + dis1 * math.sin(theta)

    ax.scatter(0, margin1_y1, marker="P", color="none", linewidths=1, s=100, edgecolors='blue')
    ax.scatter(margin1_x2, margin1_y2, marker="P", color="none", linewidths=1, s=100, edgecolors='blue')
    ax.plot([0, margin1_x2], [margin1_y1, margin1_y2], c='red')

    # 显示惩罚项系数
    ax.text(xx - 7, yy - 2, 'C= ' + str(c), fontdict={'size': '12', 'color': 'b'})

    # 显示间距
    ax.text(xx - 7, yy - 3.5, "margin=" + str('%.2f' % dis1), fontdict={'size': '12', 'color': 'b'})

    # 显示支持向量类型
    ax.scatter(xx - 6.5, yy - 6, marker='s', color="none", linewidths=1, s=100, edgecolors='red')
    ax.text(xx - 5.5, yy - 6.5, "Margin SV", fontdict={'size': '12', 'color': 'b'})

    plt.scatter(xx - 6.5, yy - 7.5, marker='^', color="none", linewidths=1, s=100, edgecolors='blue')
    ax.text(xx - 5.5, yy - 8, "Other SV", fontdict={'size': '12', 'color': 'b'})

    if os.path.exists('LinearSVM_output') == 0:
        os.mkdir('LinearSVM_output')
    plt.savefig('LinearSVM_output' + '//' + model_name + '_C' + str(n) + '.jpg', dpi=100)
    plt.show()

# 保存模型
if os.path.exists('libsvm_model') == 0:
    os.mkdir('libsvm_model')
svm_save_model('libsvm_model//' + model_name, model1)

copyright© 意疏：https://blog.csdn.net/sinat_35907936/article/details/115110787

---

参考

---

https://www.csie.ntu.edu.tw/~cjlin/libsvm/