In [1]:
# To support both python 2 and python 3
from __future__ import division, print_function, unicode_literals
# Common imports
import numpy as np
import os
# to make this notebook's output stable across runs
np.random.seed(42)
# To plot pretty figures
%matplotlib inline
import matplotlib
import matplotlib.pyplot as plt
plt.rcParams['axes.labelsize'] = 14
plt.rcParams['xtick.labelsize'] = 12
plt.rcParams['ytick.labelsize'] = 12
SVM 是完成多功能的机器学习模型，可以完成线性或者非线性的分类，回归，甚至异常值检测。且 SVM 特别适合应用在复杂但中小规模的数据集的分类问题上。
In [2]:
from sklearn.svm import SVC
from sklearn import datasets
iris = datasets.load_iris()
# petal length, width
X = iris['data'][:, (2, 3)]
y = iris['target']
# setosa 或者 versicolor
setosa_or_versicolor = (y == 0) | (y == 1)
X = X[setosa_or_versicolor]
y = y[setosa_or_versicolor]
# 初始化 svm
svm_clf = SVC(kernel='linear', C=float('inf'))
svm_clf.fit(X, y)
Out[2]:
In [3]:
# bad model
x0 = np.linspace(0, 5.5, 200)
pred_1 = 5 * x0 - 20
pred_2 = x0 - 1.8
pred_3 = 0.1 * x0 + 0.5
# 绘制 svm 的决策边界
def plot_svc_decision_boundary(svm_clf, xmin, xmax):
# 特征权重 w
w = svm_clf.coef_[0]
# 决策函数常量 b
b = svm_clf.intercept_[0]
x0 = np.linspace(xmin, xmax, 200)
# 决策边界， w0*x0 + w1*x1 + b = 0
# => x1 = -w0/w1 * x0 - b/w1
decision_boundary = -w[0]/w[1] * x0 - b/w[1]
margin = 1/w[1]
gutter_up = decision_boundary + margin
gutter_down = decision_boundary - margin
# 触碰点
svs = svm_clf.support_vectors_
# 绘制触碰点
plt.scatter(svs[:, 0], svs[:, 1], s=180, facecolors='#FFAAAA')
plt.plot(x0, decision_boundary, "k-", linewidth=2)
# 上界
plt.plot(x0, gutter_up, "k--", linewidth=2)
# 下界
plt.plot(x0, gutter_down, "k--", linewidth=2)
plt.figure(figsize=(12,2.7))
plt.subplot(121)
# 三种可能的线性分类边界
plt.plot(x0, pred_1, "g--", linewidth=2)
plt.plot(x0, pred_2, "m-", linewidth=2)
plt.plot(x0, pred_3, "r-", linewidth=2)
# versicolor
plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs", label="Iris-Versicolor")
# setosa
plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo", label="Iris-Setosa")
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.legend(loc="upper left", fontsize=14)
plt.axis([0, 5.5, 0, 2])
plt.subplot(122)
plot_svc_decision_boundary(svm_clf, 0, 5.5)
plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs")
plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo")
plt.xlabel("Petal length", fontsize=14)
plt.axis([0, 5.5, 0, 2])
plt.show()
SVM 分类器在两种类别之间保持一条尽可能宽敞的街道，其被称为最大间隔分类。
我们添加更多的样本点到街道中，也不会影响判定边界，因为判定边界由位于街道边缘的样本点确定，这些样本点被称作 支持向量 。
svm 对特征缩放敏感¶
In [4]:
# 初始化 svm
Xs = np.array([[1, 50], [5, 20], [3, 80], [5, 60]]).astype(np.float64)
ys = np.array([0, 0, 1, 1])
svm_clf = SVC(kernel='linear', C=100)
svm_clf.fit(Xs, ys)
# 画图
plt.figure(figsize=(12, 3.2))
plt.subplot(121)
# 取出 ys 为 1 的 Xs 的坐标
plt.plot(Xs[:, 0][ys==1], Xs[:, 1][ys==1], 'bo')
# 同理，取出 ys 为 0 的 Xs 的坐标
plt.plot(Xs[:, 0][ys==0], Xs[:, 1][ys==0], "ms")
# 现在情况下，支持向量的位置进行打印
plot_svc_decision_boundary(svm_clf, 0, 6)
plt.xlabel("$x_0$", fontsize=20)
plt.ylabel("$x_1$ ", fontsize=20, rotation=0)
plt.title("Unscaled", fontsize=16)
plt.axis([0, 6, 0, 90])
# 下面进行特征缩放 (标准化)
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
# 特征缩放
X_scaled = scaler.fit_transform(Xs)
svm_clf.fit(X_scaled, ys)
plt.subplot(122)
# 同理，取出 ys 为 1 的 Xs 的坐标
plt.plot(X_scaled[:, 0][ys==1], X_scaled[:, 1][ys==1], "bo")
# 同理，取出 ys 为 0 的 Xs 的坐标
plt.plot(X_scaled[:, 0][ys==0], X_scaled[:, 1][ys==0], "ms")
plot_svc_decision_boundary(svm_clf, -2, 2)
plt.xlabel("$x_0$", fontsize=20)
plt.title("Scaled", fontsize=16)
plt.axis([-2, 2, -2, 2])
plt.show()
虽然说 svm 对特征缩放敏感，但是可以发现，使用特征缩放以后，判定边界看起来好的多了。
软间隔分类¶
如果我们严格要求所有数据都不在「街道」上，都在正确的两边，这称为 硬分类 ，但是存在两个问题：
- 只对线性可分的数据起作用
- 对异常点敏感
In [5]:
X_outliers = np.array([[3.4, 1.3], [3.2, 0.8]])
y_outliers = np.array([0, 0])
# 添加第一个点和分类结果
Xo1 = np.concatenate([X, X_outliers[:1]], axis=0)
yo1 = np.concatenate([y, y_outliers[:1]], axis=0)
# 添加第二个点和分类结果
Xo2 = np.concatenate([X, X_outliers[1:]], axis=0)
yo2 = np.concatenate([y, y_outliers[1:]], axis=0)
plt.figure(figsize=(12,2.7))
plt.subplot(121)
plt.plot(Xo1[:, 0][yo1==1], Xo1[:, 1][yo1==1], "bs")
plt.plot(Xo1[:, 0][yo1==0], Xo1[:, 1][yo1==0], "yo")
plt.text(0.3, 1.0, "Impossible!", fontsize=24, color="red")
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
# 指明离群点
plt.annotate("Outlier",
xy=(X_outliers[0][0], X_outliers[0][1]), # 点位置
xytext=(2.5, 1.7), # 文字位置
ha="center",
arrowprops=dict(facecolor='black', shrink=0.1),
fontsize=16,
)
plt.axis([0, 5.5, 0, 2])
# 初始化 svm
svm_clf2 = SVC(kernel='linear', C=10**9)
# 拟合数据到 svm
svm_clf2.fit(Xo2, yo2)
plt.subplot(122)
plt.plot(Xo2[:, 0][yo2==1], Xo2[:, 1][yo2==1], "bs")
plt.plot(Xo2[:, 0][yo2==0], Xo2[:, 1][yo2==0], "yo")
# svm 边界
plot_svc_decision_boundary(svm_clf2, 0, 5.5)
plt.xlabel("Petal length", fontsize=14)
# 指明离群点
plt.annotate("Outlier",
xy=(X_outliers[1][0], X_outliers[1][1]), # 点位置
xytext=(3.2, 0.08), # 文字的位置
ha="center",
arrowprops=dict(facecolor='black', shrink=0.1),
fontsize=16,
)
plt.axis([0, 5.5, 0, 2])
plt.show()
上面的图中，左图很难找到硬间隔，右图很难一般化。
为了避免上面的问题，可以使用更加软性的模型，目的在保持街道尽可能大和避免间隔违规（出错的分类）之间找到一个良好的平衡，这就是软间隔分类。
下面的 sklearn 代码加载了内置的 Iris 数据集，进行特征缩放，并训练一个线性 SVM 模型，来检测 Virginica 鸢尾花。
In [6]:
import numpy as np
from sklearn import datasets
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.svm import LinearSVC
iris = datasets.load_iris()
X = iris['data'][:, (2, 3)]
y = (iris["target"] == 2).astype(np.float64) # Iris-Virginica
svm_clf = Pipeline([
('scalar', StandardScaler()),
('linear_svc', LinearSVC(C=1, loss='hinge', random_state=42)),
])
svm_clf.fit(X, y)
svm_clf.predict([[5.5, 1.7]])
Out[6]:
作为一种选择，我们可以使用 SVC(kernal='linear', C=1) ，但是它比较慢，尤其是在较大的训练集上，所以一般不被推荐。
另一种选择是 SGDClassifier ,使用 SGDClassifier(loss='hinge', alpha=1/(m*C))。它还是没有 LinearSVC 那样快速收敛，但是对于那些不适合放在内存的大数据集是非常有用的。
在 svm 类中，可以使用 c 超参数来控制街道尽可能大和避免间隔违规这种平衡：较小的 c 会导致更宽的街道，但更多间隔违规。
In [7]:
scaler = StandardScaler()
svm_clf1 = LinearSVC(C=1, loss='hinge', random_state=42)
svm_clf2 = LinearSVC(C=100, loss="hinge", random_state=42)
# 初始化 svm
scaled_svm_clf1 = Pipeline([
("scaler", scaler),
("linear_svc", svm_clf1),
])
scaled_svm_clf2 = Pipeline([
("scaler", scaler),
("linear_svc", svm_clf2),
])
# 拟合数据
scaled_svm_clf1.fit(X, y)
scaled_svm_clf2.fit(X, y)
Out[7]:
$$DecisionFunction=w1\cdot x1+w2\cdot x2+...+wn\cdot xn+b$$
如果 df 结果大于等于0，那么就是积极的，反之消极，这些事缩放以后的参数，我们希望使用未缩放的参数，x' 与 w' 来进行表示。
$$DecisionFunction=w'1\cdot x'1+w'2\cdot x'2+...+w'n\cdot x'n+b'$$
我们知道在 StandardScaler 中：
$$scaler:xi=(x'i-mean)/scale=>x'i=xi*scale+mean$$
因此：
$$DecisionFunction=w'1(x1*scale+mean)+...+w'n(xn*scale+mean)+b'$$ $$=w'1*scale*x1 + ... + w'n*scale*xn + w'1*mean + ... + w'n*mean + b'$$
我们可以发现：
$$wi=w'i*scale=>w'i=wi/scale$$
多余项全部归属于 b
$$b'= - w'1*mean - ... - w'n*mean + b$$ $$= - w1*mean/scale - ... - wn*mean/scale + b$$ $$= DecisionFunction(-mean/scale)$$
现在我们获得了 w' 和 b' ，我们要知道一点：支持向量是那些不是在正确一段的数据。
就好像是明明是积极的，得分却 < 1；或者明明是消极的，得分却 > 1。
这里的分很好算，就是决策分。我们可以简写成：
$$X.dot(w)+b$$
y 里存放着，正 1 和 误 0 ，但是我们可以把它转换成 正 1 和 误 -1 。用 $t = 2y-1$ 来转换。
所以说我们的目标是：
$$X.dot(w) + b < -1 是消极$$ $$X.dot(w) + b > +1 是积极$$
换句话说，就是：
$$-1 * X.dot(w) + b > 1 是消极$$ $$+1 * X.dot(w) + b > 1 是积极$$
因此下面的公式得出了支持向量：
$$t * X.dot(w) + b<1$$
In [8]:
# 这里的 b 和 w 求得是上面的 b' 和 w'
b1 = svm_clf1.decision_function([-scaler.mean_ / scaler.scale_])
b2 = svm_clf2.decision_function([-scaler.mean_ / scaler.scale_])
w1 = svm_clf1.coef_[0] / scaler.scale_
w2 = svm_clf2.coef_[0] / scaler.scale_
svm_clf1.intercept_ = np.array([b1])
svm_clf2.intercept_ = np.array([b2])
svm_clf1.coef_ = np.array([w1])
svm_clf2.coef_ = np.array([w2])
# y 里设置了 0 是消极， 1 是积极，如果 2y-1 ，那么 -1 是消极， 1 是积极
t = y * 2 - 1
support_vectors_idx1 = (t * (X.dot(w1) + b1) < 1).ravel()
support_vectors_idx2 = (t * (X.dot(w2) + b2) < 1).ravel()
svm_clf1.support_vectors_ = X[support_vectors_idx1]
svm_clf2.support_vectors_ = X[support_vectors_idx2]
In [9]:
plt.figure(figsize=(12,3.2))
plt.subplot(121)
plt.plot(X[:, 0][y==1], X[:, 1][y==1], "g^", label="Iris-Virginica")
plt.plot(X[:, 0][y==0], X[:, 1][y==0], "bs", label="Iris-Versicolor")
plot_svc_decision_boundary(svm_clf1, 4, 6)
plot_svc_decision_boundary(svm_clf1, 4, 6)
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.legend(loc="upper left", fontsize=14)
plt.title("$C = {}$".format(svm_clf1.C), fontsize=16)
plt.axis([4, 6, 0.8, 2.8])
plt.subplot(122)
plt.plot(X[:, 0][y==1], X[:, 1][y==1], "g^")
plt.plot(X[:, 0][y==0], X[:, 1][y==0], "bs")
plot_svc_decision_boundary(svm_clf2, 4, 6)
plt.xlabel("Petal length", fontsize=14)
plt.title("$C = {}$".format(svm_clf2.C), fontsize=16)
plt.axis([4, 6, 0.8, 2.8])
Out[9]:
上面的右图，使用了较大的 C ，导致更少的间隔违规，但是间隔变小了；左图反之。
非线性 svm 分类¶
很多数据不是线性可分的。我们在某些情况下可以将多项式变成线性可分的数据，下面的例子就是线性可分和不可分，如果我们只有一个特征 x1 那就是不能分的，但是如果是 x2 = (x1)^2 就变得可分了：
In [10]:
X1D = np.linspace(-4, 4, 9).reshape(-1, 1)
X2D = np.c_[X1D, X1D**2]
y = np.array([0, 0, 1, 1, 1, 1, 1, 0, 0])
plt.figure(figsize=(11, 4))
plt.subplot(121)
plt.grid(True, which='both')
plt.axhline(y=0, color='k')
plt.plot(X1D[:, 0][y==0], np.zeros(4), "bs")
plt.plot(X1D[:, 0][y==1], np.zeros(5), "g^")
plt.gca().get_yaxis().set_ticks([])
plt.xlabel(r"$x_1$", fontsize=20)
plt.axis([-4.5, 4.5, -0.2, 0.2])
plt.subplot(122)
plt.grid(True, which='both')
plt.axhline(y=0, color='k')
plt.axvline(x=0, color='k')
plt.plot(X2D[:, 0][y==0], X2D[:, 1][y==0], "bs")
plt.plot(X2D[:, 0][y==1], X2D[:, 1][y==1], "g^")
plt.xlabel(r"$x_1$", fontsize=20)
plt.ylabel(r"$x_2$", fontsize=20, rotation=0)
plt.gca().get_yaxis().set_ticks([0, 4, 8, 12, 16])
plt.plot([-4.5, 4.5], [6.5, 6.5], "r--", linewidth=3)
plt.axis([-4.5, 4.5, -1, 17])
# 对齐
plt.subplots_adjust(right=1)
plt.show()
In [11]:
# 我们在卫星数据集上可以验证这一想法。
from sklearn.datasets import make_moons
X, y = make_moons(n_samples=100, noise=0.15, random_state=42)
def plot_dataset(X, y, axes):
plt.plot(X[:, 0][y==0], X[:, 1][y==0], 'bs')
plt.plot(X[:, 0][y==1], X[:, 1][y==1], "g^")
# 限制范围
plt.axis(axes)
plt.grid(True, which='both')
plt.xlabel(r"$x_1$", fontsize=20)
plt.ylabel(r"$x_2$", fontsize=20, rotation=0)
plot_dataset(X, y, [-1.5, 2.5, -1, 1.5])
plt.show()
In [12]:
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures
polynomial_svm_clf = Pipeline([
('poly_features', PolynomialFeatures(degree=3)),
('scaler', StandardScaler()),
('svm_clf', LinearSVC(C=10, loss='hinge', random_state=42)),
])
polynomial_svm_clf.fit(X, y)
Out[12]:
In [13]:
def plot_predictions(clf, axes):
x0s = np.linspace(axes[0], axes[1], 100)
x1s = np.linspace(axes[2], axes[3], 100)
x0, x1 = np.meshgrid(x0s, x1s)
# 网格坐标
X = np.c_[x0.ravel(), x1.ravel()]
# 预测标签
y_pred = clf.predict(X).reshape(x0.shape)
# 计算置信度
y_decision = clf.decision_function(X).reshape(x0.shape)
# 边界
plt.contourf(x0, x1, y_pred, cmap=plt.cm.brg, alpha=0.2)
# 内部轮廓
plt.contourf(x0, x1, y_decision, cmap=plt.cm.brg, alpha=0.1)
plot_predictions(polynomial_svm_clf, [-1.5, 2.5, -1, 1.5])
plot_dataset(X, y, [-1.5, 2.5, -1, 1.5])
plt.show()
多项式核¶
添加多项式特征很容易实现，不仅仅在 svm ，在各种机器学习算法中都有不错的表现，但是低次数的多项式不能处理非常复杂的数据集，而高次数的多项式却产生了大量的特征，会使得模型变得慢。幸运的是，在使用 SVM 中，我们可以运用一种被称为「核技巧」的神奇数学技巧，它可以添加很多多项式，甚至有高次数的多项式的时候，但是一样产生好的结果。
不会有大量特征导数的组合爆炸，因为你并没有增加任何特征，这个技巧可以用 SVC 来实现，我们再次用卫星数据集来测试一下效果。
In [14]:
# SVC 中运用了核技巧，在使用高次数的多项式的时候，也能产生好的结果，因为他没有添加任何特征。
from sklearn.svm import SVC
poly_kernel_svm_clf = Pipeline([
('scaler', StandardScaler()),
('svm_clf', SVC(kernel='poly', degree=3, coef0=1, C=5)),
])
poly_kernel_svm_clf.fit(X, y)
Out[14]:
In [15]:
poly100_kernel_svm_clf = Pipeline([
("scaler", StandardScaler()),
("svm_clf", SVC(kernel="poly", degree=10, coef0=100, C=5))
])
poly100_kernel_svm_clf.fit(X, y)
Out[15]:
In [16]:
plt.figure(figsize=(11, 4))
plt.subplot(121)
plot_predictions(poly_kernel_svm_clf, [-1.5, 2.5, -1, 1.5])
plot_dataset(X, y, [-1.5, 2.5, -1, 1.5])
plt.title(r"$d=3, r=1, C=5$", fontsize=18)
plt.subplot(122)
plot_predictions(poly100_kernel_svm_clf, [-1.5, 2.5, -1, 1.5])
plot_dataset(X, y, [-1.5, 2.5, -1, 1.5])
plt.title(r"$d=10, r=100, C=5$", fontsize=18)
plt.show()
上左图中用了3阶的多项式核来训练一个 SVM 分类器，右边使用了10阶的多项式核 SVM 分类器，显然，我们的模型过拟合了，超参数 cof0 控制了高阶多项式和低阶多项式对模型的影响。
回顾： c 是超参数，用来控制街道尽可能大和避免间隔违规的平衡。又或者说升降拟合度，C 增大，增加拟合； C 减小，降低拟合。
上面我们对于 svm 解决非线性问题的解答是通过 PolynomialFeatures 将特征提高维度，然后 StandardScaler 对特征进行放缩，最后用 LinearSVC 来拟合预测。如果需要高次数的多项式，可以用 SVC 的核技巧，比如 SVC(kernel='poly', degree=10) ，这是多项式核。
增加相似特征¶
除了上面的方法解决非线性问题，我们还能使用相似函数来计算每个样本与特定 landmark 的相似度。
例如，我们前面的一维数据集，并在 x1=-2 和 x2=1 之间添加两个 landmark ，接下来我们定义一个相似函数，即高斯径向基函数 RBF ，设置 gamma=0.3 ， RBF:
$$\varphi _{\gamma }(x,\iota )=exp(-\gamma \left \| x-\iota \right \|^{2})$$
In [17]:
# 计算 RBF
def gaussian_rbf(x, landmark, gamma):
return np.exp(-gamma * np.linalg.norm(x - landmark, axis=1)**2)
gamma = 0.3
x1s = np.linspace(-4.5, 4.5, 200).reshape(-1, 1)
# 增加的两个 landmark
x2s = gaussian_rbf(x1s, -2, gamma)
x3s = gaussian_rbf(x1s, 1, gamma)
XK = np.c_[gaussian_rbf(X1D, -2, gamma), gaussian_rbf(X1D, 1, gamma)]
yk = np.array([0, 0, 1, 1, 1, 1, 1, 0, 0])
plt.figure(figsize=(11, 4))
plt.subplot(121)
plt.grid(True, which='both')
plt.axhline(y=0, color='k')
plt.scatter(x=[-2, 1], y=[0, 0], s=150, alpha=0.5, c="red") # 画点
plt.plot(X1D[:, 0][yk==0], np.zeros(4), "bs")
plt.plot(X1D[:, 0][yk==1], np.zeros(5), "g^")
plt.plot(x1s, x2s, "g--")
plt.plot(x1s, x3s, "b:")
plt.gca().get_yaxis().set_ticks([0, 0.25, 0.5, 0.75, 1])
plt.xlabel(r"$x_1$", fontsize=20)
plt.ylabel(r"Similarity", fontsize=14)
plt.annotate(r'$\mathbf{x}$',
xy=(X1D[3, 0], 0),
xytext=(-0.5, 0.20),
ha="center",
arrowprops=dict(facecolor='black', shrink=0.1),
fontsize=18,
)
plt.text(-2, 0.9, "$x_2$", ha="center", fontsize=20)
plt.text(1, 0.9, "$x_3$", ha="center", fontsize=20)
plt.axis([-4.5, 4.5, -0.1, 1.1])
plt.subplot(122)
plt.grid(True, which='both')
plt.axhline(y=0, color='k')
plt.axvline(x=0, color='k')
plt.plot(XK[:, 0][yk==0], XK[:, 1][yk==0], "bs")
plt.plot(XK[:, 0][yk==1], XK[:, 1][yk==1], "g^")
plt.xlabel(r"$x_2$", fontsize=20)
plt.ylabel(r"$x_3$ ", fontsize=20, rotation=0)
plt.annotate(r'$\phi\left(\mathbf{x}\right)$',
xy=(XK[3, 0], XK[3, 1]),
xytext=(0.65, 0.50),
ha="center",
arrowprops=dict(facecolor='black', shrink=0.1),
fontsize=18,
)
plt.plot([-0.1, 1.1], [0.57, -0.1], "r--", linewidth=3)
plt.axis([-0.1, 1.1, -0.1, 1.1])
plt.subplots_adjust(right=1)
plt.show()
上面创建了两个 landmark ，可以看到，它是一个从 0 到 1 形函数，值为 0 的里 landmark 远， 1 就在 landmark 上。现在我们准备计算新特征，例如，我们看一下样本 x1 = -1: 它距离第一个 landmark 距离是 1，距离第二个 landmark 是 2，因此，它的新特征是：
x2=exp(-0.3 × (1^2))≈0.74
和 x3=exp(-0.3 × (2^2))≈0.30
右图展示特征转换以后的数据集，可以发现，它已经变得线性可分了。
ok ，我们如果选择 landmark 呢？最简单的方法是在数据集的每一个样本的位置到创建一个 landmark ，这将产生更多的维度从而增加了转换后数据集的线性可分的可能性，但缺点是， m 个样本，n 个特征的数据集，被转换成了 m 个实例， m 个特征的数据集，因为每个点到某一点成为了一个特征。这样一来，训练集非常大，也会得到同样大小的特征。
In [18]:
x1_example = X1D[3, 0]
for landmark in (-2, 1):
k = gaussian_rbf(np.array([[x1_example]]), np.array([[landmark]]), gamma)
print("Phi({}, {}) = {}".format(x1_example, landmark, k))
高斯 RBF 核¶
就像多项式特征法一样，相似特征法对各种机器学习算法同样也有不错的表现，但是在所有额外特征上的计算成本可能很高，特别是在大规模的训练集上。然而，「核」技巧再次显现了它在 SVM 上的神奇之处，高斯核让你可以获得同样好的结果成为可能，就像你在相似特征法添加许多相似特征那样，事实上，你并不需要 RBF 来添加他们，使用 SVC 类的高斯 RBF 核来检验一下。
In [19]:
X, y = make_moons(n_samples=100, noise=0.15, random_state=42)
In [20]:
gamma1, gamma2 = 0.1, 5
C1, C2 = 0.001, 1000
hyperparams = (gamma1, C1), (gamma1, C2), (gamma2, C1), (gamma2, C2)
svm_clfs = []
for gamma, C in hyperparams:
rbf_kernel_svm_clf = Pipeline([
('scaler', StandardScaler()),
('svm_clf', SVC(kernel='rbf', gamma=gamma, C=C)),
])
rbf_kernel_svm_clf.fit(X, y)
svm_clfs.append(rbf_kernel_svm_clf)
plt.figure(figsize=(11, 7))
for i, svm_clf in enumerate(svm_clfs):
plt.subplot(221 + i)
plot_predictions(svm_clf, [-1.5, 2.5, -1, 1.5])
plot_dataset(X, y, [-1.5, 2.5, -1, 1.5])
gamma, C = hyperparams[i]
plt.title(r"$\gamma = {}, C = {}$".format(gamma, C), fontsize=16)
plt.show()
从上面的图中，我们可以发现超参数 gamma 和 C 对模型的影响。
增大 gamma 将导致钟形变窄，导致每个样本的影响范围变得更小：判定边界最终变得更不规则，在单个样本周围环绕。相反，较小的 gamma 使得曲线变得更宽，样本有更大的影响范围，判定边界最终更加平滑。所以说 gamma 是可调整的超参数：如果模型过拟合，应该减小 gamma ；反之，增大。
除了这些，还有很多可供选择的超参数，比如：对文本文档和 DNA 序列进行分类的时候，可以使用字符串核 (String kernels) ，或者 ssk (string subsequence kernel)或者 基于编辑距离 levenshtein distance 的核函数。
我们可以先使用线性核函数， LinearSVC 比 SVC(kernel='linear') 快得多。
计算复杂度¶
LinearSVC 基于 liblinear 库，它实现了线性 SVM 的优化算法，它并不支持核技巧，但是它样本和特征的数量几乎是线性的，训练复杂度大约是 O(m*n)
如果需要非常高的精度，这个算法需要花费更多的时间，这是由容差值超参数 tol 控制的，大多数分类问题中，默认就够了。
SVC 类基于 libsvm 库，它实现了支持核技巧的算法。训练时间复杂度通常介于 O(m^2 n) 和 O(m^3 n)之间，不幸的是，这意味着当训练样本变大的时候，它将变得极其慢，这个算法对于复杂但小型或中等数量的训练集是完美的。它能对特征数量进行很好地缩放，尤其对系数矩阵来说。在这个情况下，算法对每个样本的非零特征的平均数量进行大概的缩放。
SVM 回归¶
SVM 算法应用广泛，不仅仅支持线性和非线性的分类任务，还支持线性和非线性的回归任务。
技巧在于逆转我们的目标：限制间隔违规的情况下，不是试图在两个类别之间找到尽可能大的街道（间隔）。
SVM 回归任务是限制间隔违规情况下，尽量放置更多样本在「街道」上。
街道的宽度由 tol 超参数控制，下图显示了一些随机生成的线性数据上，两个线性 SVM 回归模型的训练情况。一个有较大的间隔 tol=1.5 ，另一个间隔较小 tol=0.5
In [21]:
np.random.seed(42)
m = 50
X = 2 * np.random.rand(m, 1)
y = (4 + 3*X + np.random.randn(m, 1)).ravel()
In [22]:
from sklearn.svm import LinearSVR
svm_reg1 = LinearSVR(epsilon=1.5, random_state=42)
svm_reg2 = LinearSVR(epsilon=0.5, random_state=42)
svm_reg1.fit(X, y)
svm_reg2.fit(X, y)
def find_support_vectors(svm_reg, X, y):
y_pred = svm_reg.predict(X)
off_margin = (np.abs(y - y_pred) >= svm_reg.epsilon)
return np.argwhere(off_margin)
svm_reg1.support_ = find_support_vectors(svm_reg1, X, y)
svm_reg2.support_ = find_support_vectors(svm_reg2, X, y)
eps_x1 = 1
eps_y_pred = svm_reg1.predict([[eps_x1]])
eps_y_pred
Out[22]:
In [23]:
def plot_svm_regression(svm_reg, X, y, axes):
x1s = np.linspace(axes[0], axes[1], 100).reshape(100, 1)
y_pred = svm_reg.predict(x1s)
plt.plot(x1s, y_pred, 'k-', linewidth=2, label=r"$\hat{y}$")
plt.plot(x1s, y_pred + svm_reg.epsilon, "k--")
plt.plot(x1s, y_pred - svm_reg.epsilon, "k--")
# 画出支持向量
plt.scatter(X[svm_reg.support_], y[svm_reg.support_], s=180, facecolors='#FFAAAA')
plt.plot(X, y, "bo")
plt.xlabel(r"$x_1$", fontsize=18)
plt.legend(loc="upper left", fontsize=18)
plt.axis(axes)
plt.figure(figsize=(9, 4))
plt.subplot(121)
plot_svm_regression(svm_reg1, X, y, [0, 2, 3, 11])
plt.title(r"$\epsilon = {}$".format(svm_reg1.epsilon), fontsize=18)
plt.ylabel(r"$y$", fontsize=18, rotation=0)
plt.annotate(
'', xy=(eps_x1, eps_y_pred), xycoords='data',
xytext=(eps_x1, eps_y_pred - svm_reg1.epsilon),
textcoords='data', arrowprops={'arrowstyle': '<->', 'linewidth': 1.5}
)
plt.text(0.91, 5.6, r"$\epsilon$", fontsize=20)
plt.subplot(122)
plot_svm_regression(svm_reg2, X, y, [0, 2, 3, 11])
plt.title(r"$\epsilon = {}$".format(svm_reg2.epsilon), fontsize=18)
plt.show()
所以，从上面的例子中，我们也从侧面发现一点，那就是添加更多数据样本在间隔之内并不会影响模型的预测。因此也认为该模型是不敏感的。
我们可以使用 sklearn 的 LinearSVR 类来实现线性回归。
同时 SVR 类也支持核技巧，在回归任务上， SVR 类和 SVC 类是一样的，并且 LinearSVR 和 LinearSVC 是等价的。
In [25]:
np.random.seed(42)
m = 100
X = 2 * np.random.rand(m, 1) - 1
y = (0.2 + 0.1 * X + 0.5 * X**2 + np.random.randn(m, 1)/10).ravel()
In [27]:
from sklearn.svm import SVR
svm_poly_reg1 = SVR(kernel='poly', degree=2, C=100, epsilon=0.1)
svm_poly_reg2 = SVR(kernel='poly', degree=2, C=0.01, epsilon=0.1)
svm_poly_reg1.fit(X, y)
svm_poly_reg2.fit(X, y)
Out[27]:
In [31]:
plt.figure(figsize=(9, 4))
plt.subplot(121)
plot_svm_regression(svm_poly_reg1, X, y, [-1, 1, 0, 1])
plt.title(r"$degree={}, C={}, \epsilon = {}$".format(svm_poly_reg1.degree, svm_poly_reg1.C, svm_poly_reg1.epsilon), fontsize=18)
plt.ylabel(r"$y$", fontsize=18, rotation=0)
plt.subplot(122)
plot_svm_regression(svm_poly_reg2, X, y, [-1, 1, 0, 1])
plt.title(r"$degree={}, C={}, \epsilon = {}$".format(svm_poly_reg2.degree, svm_poly_reg2.C, svm_poly_reg2.epsilon), fontsize=18)
plt.show()
为了明白 C 的实际作用，我现在从吴恩达 过拟合问题 出发，寻找答案。
我从 吴恩达2014 的 第三周第七讲 到 第7七周第12讲 中找到了现在我所需要的东西。
我们发现上面的 C 在变大的时候拟合的更好，这是因为正则化不明显，或者说，这是约束力更差的表现，更接近无约束的情况，因为 C 作为正则化参数的倒数存在，和强度成反比。
背后机制¶
现在开始，我们让偏置项叫做 b ，将权重称作 w 。
决策函数和预测¶
线性 svm 分类器简单地使用决策 w*b 来预测新样本的类别：如果结果为正，就是正类，为1；反之为负类，为0。
In [32]:
iris = datasets.load_iris()
X = iris['data'][:, (2, 3)]
y = (iris["target"] == 2).astype(np.float64) # Iris-Virginica
In [54]:
from mpl_toolkits.mplot3d import Axes3D
def plot_3D_decision_function(ax, w, b, x1_lim=[4, 6], x2_lim=[0.8, 2.8]):
x1_in_bounds = (X[:, 0] > x1_lim[0]) & (X[:, 0] < x1_lim[1])
X_crop = X[x1_in_bounds]
y_crop = y[x1_in_bounds]
x1s = np.linspace(x1_lim[0], x1_lim[1], 20)
x2s = np.linspace(x2_lim[0], x2_lim[1], 20)
x1, x2 = np.meshgrid(x1s, x2s)
xs = np.c_[x1.ravel(), x2.ravel()]
df = (xs.dot(w) + b).reshape(x1.shape)
m = 1 / np.linalg.norm(w)
boundary_x2s = -x1s*(w[0]/w[1])-b/w[1]
margin_x2s_1 = -x1s*(w[0]/w[1])-(b-1)/w[1]
margin_x2s_2 = -x1s*(w[0]/w[1])-(b+1)/w[1]
ax.plot_surface(x1s, x2, np.zeros_like(x1),
color="b", alpha=0.2, cstride=100, rstride=100)
# 边界和街道
ax.plot(x1s, boundary_x2s, 0, "k-", linewidth=2, label=r"$h=0$")
ax.plot(x1s, margin_x2s_1, 0, "k--", linewidth=2, label=r"$h=\pm 1$")
ax.plot(x1s, margin_x2s_2, 0, "k--", linewidth=2)
ax.plot(X_crop[:, 0][y_crop==1], X_crop[:, 1][y_crop==1], 0, "g^")
ax.plot(X_crop[:, 0][y_crop==0], X_crop[:, 1][y_crop==0], 0, "bs")
# 决策平面
ax.plot_wireframe(x1, x2, df, alpha=0.3, color="k")
ax.axis(x1_lim + x2_lim)
ax.text(4.5, 2.5, 3.8, "Decision function $h$", fontsize=14)
ax.set_xlabel(r"Petal length", fontsize=14)
ax.set_ylabel(r"Petal width", fontsize=14)
ax.set_zlabel(r"$h = \mathbf{w}^T \mathbf{x} + b$", fontsize=14)
ax.legend(loc="upper left", fontsize=14)
fig = plt.figure(figsize=(11, 6))
ax1 = fig.add_subplot(111, projection='3d')
plot_3D_decision_function(ax1, w=svm_clf2.coef_[0], b=svm_clf2.intercept_[0])
plt.show()
从图中可以发现，现在我们的决策函数的平面 h 的位置。
训练目标¶
在决策函数下的斜率：它等于权重向量的范数 ||w|| 。如果我们把这个斜率除以2，决策函数等于 +-1 的点将会离决策边界到原来的两倍大。
In [55]:
def plot_2D_decision_function(w, b, ylabel=True, x1_lim=[-3, 3]):
x1 = np.linspace(x1_lim[0], x1_lim[1], 200)
y = w * x1 + b
m = 1 / w
plt.plot(x1, y)
plt.plot(x1_lim, [1, 1], "k:")
plt.plot(x1_lim, [-1, -1], "k:")
plt.axhline(y=0, color='k')
plt.axvline(x=0, color='k')
plt.plot([m, m], [0, 1], "k--")
plt.plot([-m, -m], [0, -1], "k--")
plt.plot([-m, m], [0, 0], "k-o", linewidth=3)
plt.axis(x1_lim + [-2, 2])
plt.xlabel(r"$x_1$", fontsize=16)
if ylabel:
plt.ylabel(r"$w_1 x_1$ ", rotation=0, fontsize=16)
plt.title(r"$w_1 = {}$".format(w), fontsize=16)
plt.figure(figsize=(12, 3.2))
plt.subplot(121)
plot_2D_decision_function(1, 0)
plt.subplot(122)
plot_2D_decision_function(0.5, 0, ylabel=False)
plt.show()
Hinge loss¶
In [58]:
t = np.linspace(-2, 4, 200)
h = np.where(1 - t < 0, 0, 1 - t) # max(0, 1-t)
plt.figure(figsize=(5,2.8))
plt.plot(t, h, "b-", linewidth=2, label="$max(0, 1 - t)$")
plt.grid(True, which='both')
plt.axhline(y=0, color='k')
plt.axvline(x=0, color='k')
plt.yticks(np.arange(-1, 2.5, 1))
plt.xlabel("$t$", fontsize=16)
plt.axis([-2, 4, -1, 2.5])
plt.legend(loc="upper right", fontsize=16)
plt.show()
