import numpy as np from utils import common import matplotlib.pyplot as plt from tqdm import tqdm import scipy.spatial.distance as ssd

导入所需的库

def RBF_kernel(sigma): def kernel(x, y): ''' :param x: N x d :param y: M x d :return: N x M ''' dist = ssd.cdist(x, y, 'euclidean') # 计算欧氏距离 return np.exp(-np.square(dist) * sigma) # 计算RBF核函数的值

return kernel

定义RBF核函数

def auto_scale(X): X_var = np.var(X,ddof=1) # 计算X的方差 sigma = 1.0 / (X.shape[1] * X_var) if X_var != 0 else 1.0 # 计算自动缩放的参数sigma return sigma

自动缩放数据

def linear_kernel(): def kernel(x, y): ''' :param x: N x d :param y: M x d :return: N x M ''' return x @ y.T # 计算线性核函数的值

return kernel

定义线性核函数

class SVM: def init(self, C, kernel, max_iter=500,show_fitting_bar=True): self.kernel = kernel # 核函数 self.C = C # 惩罚参数 self.max_iter = max_iter # 最大迭代次数 self.biary_encoder = common.BinaryEncoder() # 二进制编码器 self.tqdm = tqdm if show_fitting_bar else lambda x:x # 进度条显示

def error(self, index):
    pred = (self.alpha * self._label).T @ self.Q[:, index] + self.b  # 计算预测值
    return pred - self._label[index]  # 计算误差

def fit(self, X: np.ndarray, y: np.ndarray, tol: float = 1e-3):
    '''
    :param X: N x d
    :param y: N
    :param tol: scalar
    :return:
    '''
    N, d = X.shape
    self._X = X  # 训练数据
    self._label = self.biary_encoder.fit_transform(y)  # 标签编码
    self.Q = self.kernel(X, X)  # 计算核矩阵
    self.alpha = np.zeros(N)  # 初始化拉格朗日乘子
    self.b = 0.0  # 初始化截距
    satisfied_count = 0

    for iter_count in self.tqdm(range(self.max_iter)):  # 迭代训练
        alpha_changed_count = 0
        for i in range(N):
            Ei = self.error(i)  # 计算第i个样本的误差
            alphai_old = self.alpha[i].copy()  # 保存旧的拉格朗日乘子
            yi = self._label[i]  # 第i个样本的标签

            # 选择违反KKT条件的alpha i
            if (yi * Ei < -tol and alphai_old < self.C) or (yi * Ei > tol and alphai_old > 0):
                # 随机选择alpha j
                j = np.random.randint(0, N)
                while j == i:
                    j = np.random.randint(0, N)
                yj = self._label[j]  # 第j个样本的标签
                Ej = self.error(j)  # 计算第j个样本的误差
                alphaj_old = self.alpha[j].copy()  # 保存旧的拉格朗日乘子

                # 计算边界L和H
                if yi == yj:
                    L = max(0.0, alphai_old + alphaj_old - self.C)
                    H = min(self.C, alphai_old + alphaj_old)
                else:
                    L = max(0.0, alphaj_old - alphai_old)
                    H = min(self.C, self.C + alphaj_old - alphai_old)

                # 计算eta
                eta = 2.0 * self.Q[i, j] - self.Q[i, i] - self.Q[j, j]
                if eta >= 0:
                    continue

                # 更新alpha j
                self.alpha[j] =alphaj_old - (yj * (Ei - Ej) / eta)

                # 限制alpha j的取值范围
                self.alpha[j] = np.clip(self.alpha[j], L, H)

                # 更新alpha i
                delta_j = alphaj_old - self.alpha[j]
                if np.abs(delta_j) < tol:
                    continue
                self.alpha[i] =alphai_old + yi * yj * delta_j

                # 更新截距b
                bi = self.b - Ei \
                     - yi * (self.alpha[i] - alphai_old) * self.Q[i, i] \
                     - yj * (self.alpha[j] - alphaj_old) * self.Q[i, j]
                bj = self.b - Ej \
                     - yi * (self.alpha[i] - alphai_old) * self.Q[i, j] \
                     - yj * (self.alpha[j] - alphaj_old) * self.Q[j, j]
                if 0.0+tol < self.alpha[i] < self.C-tol:
                    self.b = bi
                elif 0.0+tol < self.alpha[j] < self.C-tol:
                    self.b = bj
                else:
                    self.b = (bi + bj) / 2.0

                alpha_changed_count += 1
        if alpha_changed_count == 0:
            satisfied_count += 1
        else:
            satisfied_count = 0
        if satisfied_count >= 5:
            break

    mask = self.alpha > tol  # 选取大于阈值的拉格朗日乘子
    self.support_vectors = mask  # 支持向量
    self.X_sup = X[mask]  # 支持向量的特征
    self.Y_sup = self._label[mask]  # 支持向量的标签
    self.alpha_sup = self.alpha[mask]  # 支持向量的拉格朗日乘子

    del self._X
    del self._label
    del self.alpha
    del self.Q

def predict_proba(self, X):
    sigma = self.kernel(X, self.X_sup)  # 计算测试样本与支持向量的核矩阵
    sigma = self.Y_sup[np.newaxis, :] * self.alpha_sup[np.newaxis, :] * sigma  # 计算预测值
    y_hat = np.sum(sigma, axis=1) + self.b  # 计算预测概率
    return y_hat

def predict(self, X):
    r = self.predict_proba(X)  # 预测概率
    r[r>0]=1
    r[r<=0]=-1
    r = self.biary_encoder.inverse_transform(r)  # 将预测概率转换为标签
    return r

定义支持向量机类

def plot_boundary():

def generate_data(n_sample=50):
    x1 = np.random.rand(n_sample, 1) * 6 - 2
    y1 = x1 + np.random.normal(0, 0.7, (n_sample, 1))
    x2 = np.random.rand(n_sample, 1) * 6 - 3
    y2 = x2 + np.random.normal(0, 0.7, (n_sample, 1)) + 3
    c1 = np.hstack([x1, y1, np.ones_like(x1) * -1])
    c2 = np.hstack([x2, y2, np.ones_like(x2)])
    data_with_label = np.vstack([c1, c2])
    data, label = data_with_label[:, :-1], data_with_label[:, -1]

    return data, label

import matplotlib.pyplot as plt
np.random.seed(1)
X,y = generate_data(50)
kernels = [RBF_kernel(1.0),linear_kernel()]
plt.figure(figsize=(10,5))
for i,kernel in enumerate(kernels):
    plt.subplot(1,2,i+1)
    model = SVM(C=1, kernel=kernel,max_iter=1000)
    model.fit(X, y, tol=1e-5)
    sv = model.support_vectors

    plt.scatter(X[sv, 0], X[sv, 1], c=y[sv], marker='*',
                linewidths=0.5, edgecolors=(0, 0, 0, 1))
    plt.scatter(X[~sv, 0], X[~sv, 1], c=y[~sv],
                linewidths=0.5, edgecolors=(0, 0, 0, 1))
    xvals = np.linspace(-3, 6, 200)
    yvals = np.linspace(-3, 6, 200)
    xx, yy = np.meshgrid(xvals, yvals)

    pred = model.predict_proba(np.c_[xx.ravel(), yy.ravel()])
    zz = np.reshape(pred, xx.shape)
    plt.pcolormesh(xx, yy, zz, zorder=0)
    plt.contour(xx, yy, zz, levels=(-1, 0, 1), colors='w', linewidths=1.5, zorder=1, linestyles='solid')

    plt.xlim([