In [ ]:
import numpy as np
import matplotlib.pyplot as plt
import scipy
from scipy import stats
from numpy.random import randn # Gaussian random number
import cvxopt
from cvxopt import matrix
In [ ]:
a = randn(1)
b = randn(1)
n = 100
X = randn(n, 2)
y = np.sign(a*X[:, 0] + b*X[:, 1] + 0.1*randn(n))
y = y.reshape(-1, 1) # The shape should be explicit
In [ ]:
def svm_1(X, y, C):
eps = 0.0001
n = X.shape[0]
P = np.zeros((n, n))
for i in range(n):
for j in range(n):
P[i, j] = np.dot(X[i, :], X[j, :]) * y[i] * y[j]
# Use the matrix function in the package to spacify them
P = matrix(P + np.eye(n)*eps)
A = matrix(- y.T.astype(np.float64))
b = matrix(np.array([0]).astype(np.float64))
h = matrix(np.array([C]*n + [0]*n).reshape(-1, 1).astype(np.float64))
G = matrix(np.concatenate([np.diag(np.ones(n)), np.diag(-np.ones(n))]))
q = matrix(np.array([-1]*n).astype(np.float64))
res = cvxopt.solvers.qp(P, q, A=A, b=b, G=G, h=h) # Execute the solver
alpha = np.array(res["x"]) # x is the alpha in the text
beta = ((alpha*y).T@X).reshape(2, 1)
index = (eps < alpha[:, 0]) & (alpha[:, 0] < C - eps)
beta_0 = np.mean(y[index] - X[index, :]@beta)
return {"beta":beta, "beta_0":beta_0}
In [ ]:
def f(x):
return - res["beta_0"]/res["beta"][1] - x*res["beta"][0]/res["beta"][1]
In [ ]:
a = randn(1)
b = randn(1)
n = 100
X = randn(n, 2)
y = np.sign(a*X[:, 0] + b*X[:, 1] + 0.1*randn(n))
y = y.reshape(-1, 1) # The shape should be explicit.
for i in range(n):
if y[i] == 1:
plt.scatter(X[i, 0], X[i, 1], c="red")
else :
plt.scatter(X[i, 0], X[i, 1], c="blue")
res = svm_1(X, y, C=10)
x_seq = np.arange(-3, 3, 0.5)
plt.plot(x_seq, f(x_seq))
<ipython-input-14-9ae01250edad>:7: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.) P[i, j] = np.dot(X[i, :], X[j, :]) * y[i] * y[j]
pcost dcost gap pres dres 0: -1.6325e+02 -8.9101e+03 3e+04 1e+00 1e-14 1: -3.0366e+01 -3.1994e+03 5e+03 1e-01 1e-14 2: 1.4212e+01 -4.3868e+02 6e+02 1e-02 3e-14 3: -6.3227e+01 -1.6285e+02 1e+02 1e-03 7e-15 4: -8.2897e+01 -1.3421e+02 5e+01 5e-04 5e-15 5: -9.4134e+01 -1.2559e+02 3e+01 2e-04 5e-15 6: -1.0001e+02 -1.1851e+02 2e+01 9e-05 5e-15 7: -1.0399e+02 -1.1207e+02 8e+00 3e-05 5e-15 8: -1.0717e+02 -1.0812e+02 1e+00 2e-16 5e-15 9: -1.0755e+02 -1.0763e+02 7e-02 2e-15 6e-15 10: -1.0758e+02 -1.0759e+02 2e-03 7e-15 5e-15 11: -1.0758e+02 -1.0758e+02 2e-05 2e-15 6e-15 Optimal solution found.
Out[ ]:
[<matplotlib.lines.Line2D at 0x7c2524c3e170>]
84¶
In [ ]:
def K_linear(x, y):
return x.T@y
In [ ]:
# The polynomial kernel may be constructed as folows:
def K_poly(x, y):
return (1 + x.T@y) ** 2
85¶
In [ ]:
# Execution
a = 3
b = -1
n = 200
X = randn(n, 2)
y = np.sign(a*X[:, 0] + b*X[:, 1]**2 + 0.3*randn(n))
y = y.reshape(-1, 1)
In [ ]:
def plot_kernel(K, line): # Specify the line format using ``line"
res = svm_2(X, y, 1, K)
alpha = res["alpha"][:, 0]
beta_0 = res["beta_0"]
def f(u, v):
S = beta_0
for i in range(X.shape[0]):
S = S + alpha[i]*y[i]*K(X[i, :], [u, v])
return S[0]
uu = np.arange(-2, 2, 0.1)
vv = np.arange(-2, 2, 0.1)
ww = []
for v in vv:
w = []
for u in uu:
w.append(f(u, v))
ww.append(w)
plt.contour(uu, vv, ww, levels=0, linestyles=line)
In [ ]:
def svm_2(X, y, C, K):
eps = 0.0001
n = X.shape[0]
P = np.zeros((n, n))
for i in range(n):
for j in range(n):
P[i, j] = K(X[i, :], X[j, :])*y[i]*y[j]
# Use the matrix function in the package to spacify them
P = matrix(P + np.eye(n)*eps)
A = matrix(-y.T.astype(np.float64))
b = matrix(np.array([0]).astype(np.float64))
h = matrix(np.array([C]*n + [0]*n).reshape(-1, 1).astype(np.float64))
G = matrix(np.concatenate([np.diag(np.ones(n)), np.diag(-np.ones(n))]))
q = matrix(np.array([-1]*n).astype(np.float64))
res = cvxopt.solvers.qp(P, q, A=A, b=b, G=G, h=h)
alpha = np.array(res["x"]) # x is the alpha in the text
beta = ((alpha*y).T@X).reshape(2, 1)
index = (eps < alpha[:, 0]) & (alpha[:, 0] < C - eps)
beta_0 = np.mean(y[index] - X[index, :]@beta)
return {"alpha":alpha, "beta":beta, "beta_0":beta_0}
uu = np.arange(-2, 2, 0.1)
vv = np.arange(-2, 2, 0.1)
ww = []
for v in vv:
w = []
for u in uu:
w.append(f(u, v))
ww.append(w)
plt.contour(uu, vv, ww, levels=0, linestyles=line)
In [ ]:
for i in range(n):
if y[i] == 1:
plt.scatter(X[i, 0], X[i, 1], c="red")
else:
plt.scatter(X[i, 0], X[i, 1], c="blue")
plot_kernel(K_poly, line="dashed")
plot_kernel(K_linear, line="solid")
<ipython-input-23-afc9250a7f7c>:7: DeprecationWarning: Conversion of an array with ndim > 0 to a scalar is deprecated, and will error in future. Ensure you extract a single element from your array before performing this operation. (Deprecated NumPy 1.25.) P[i, j] = K(X[i, :], X[j, :])*y[i]*y[j]
pcost dcost gap pres dres 0: -6.9533e+01 -4.9851e+02 3e+03 3e+00 2e-14 1: -4.3954e+01 -3.2535e+02 6e+02 5e-01 1e-14 2: -2.7322e+01 -1.3056e+02 2e+02 1e-01 7e-15 3: -1.9699e+01 -5.0913e+01 5e+01 4e-02 4e-15 4: -1.9472e+01 -2.7246e+01 1e+01 8e-03 4e-15 5: -2.0169e+01 -2.2575e+01 3e+00 2e-03 4e-15 6: -2.0477e+01 -2.1305e+01 1e+00 3e-04 4e-15 7: -2.0647e+01 -2.0954e+01 3e-01 6e-05 4e-15 8: -2.0766e+01 -2.0787e+01 2e-02 4e-06 4e-15 9: -2.0775e+01 -2.0775e+01 4e-04 6e-08 5e-15 10: -2.0775e+01 -2.0775e+01 8e-06 9e-10 5e-15 Optimal solution found. pcost dcost gap pres dres 0: -8.5617e+01 -4.6975e+02 2e+03 3e+00 6e-15 1: -5.7101e+01 -3.0106e+02 4e+02 3e-01 3e-15 2: -4.7762e+01 -8.7320e+01 5e+01 2e-02 7e-15 3: -5.4901e+01 -6.4934e+01 1e+01 4e-03 3e-15 4: -5.6952e+01 -6.1729e+01 5e+00 2e-03 3e-15 5: -5.7962e+01 -6.0085e+01 2e+00 5e-04 4e-15 6: -5.8414e+01 -5.9299e+01 9e-01 7e-05 3e-15 7: -5.8722e+01 -5.8943e+01 2e-01 6e-06 3e-15 8: -5.8819e+01 -5.8830e+01 1e-02 3e-07 3e-15 9: -5.8824e+01 -5.8824e+01 3e-04 1e-09 3e-15 10: -5.8824e+01 -5.8824e+01 3e-06 1e-11 3e-15 Optimal solution found.
86¶
In [ ]:
import sklearn
from sklearn import svm
In [ ]:
x = randn(200, 2)
x[0:100, ] = x[0:100, ] + 2
x[100:150, ] = x[100:150, ] - 2
y = np.concatenate(([1 for i in range(150)], [2 for i in range(50)]))
train = np.random.choice(200, 100, replace=False)
test = list(set(range(200)) - set(train))
res_svm = svm.SVC(kernel="rbf", gamma=1, C=100) # SVM without tuning
res_svm.fit(x[train, ], y[train]) # Execution
Out[ ]:
SVC(C=100, gamma=1)In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
SVC(C=100, gamma=1)
In [ ]:
res_svm.predict(x[test, ]) # Prediction of the test data
Out[ ]:
array([1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 1, 1, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 1, 1, 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 2])
In [ ]:
import mlxtend
from mlxtend.plotting import plot_decision_regions
In [ ]:
plot_decision_regions(x, y, clf=res_svm)
Out[ ]:
<Axes: >
In [ ]:
from sklearn.model_selection import GridSearchCV
In [ ]:
grid = {"C": [0.1, 1, 10, 100, 1000], "gamma" : [0.5, 1, 2, 3, 4]}
tune = GridSearchCV(svm.SVC(), grid, cv=10)
tune.fit(x[train, ], y[train])
Out[ ]:
GridSearchCV(cv=10, estimator=SVC(), param_grid={'C': [0.1, 1, 10, 100, 1000], 'gamma': [0.5, 1, 2, 3, 4]})In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
GridSearchCV(cv=10, estimator=SVC(), param_grid={'C': [0.1, 1, 10, 100, 1000], 'gamma': [0.5, 1, 2, 3, 4]})
SVC()
SVC()
87¶
In [ ]:
from sklearn.datasets import load_iris
In [ ]:
iris = load_iris()
iris.target_names
x = iris.data
y = iris.target
train = np.random.choice(150, 120, replace=False)
test = np.ones(150, dtype=bool)
test[train] = False
iris_svm = svm.SVC(kernel="rbf", gamma=1, C=10)
iris_svm.fit(x[train, ], y[train])
Out[ ]:
SVC(C=10, gamma=1)In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
SVC(C=10, gamma=1)
In [ ]:
# the function table_count in Chapter 3
def table_count(m, u, v):
n = u.shape[0]
count = np.zeros([m, m])
for i in range(n):
count[int(u[i]), int(v[i])] += 1
return (count)
In [ ]:
y_pre = iris_svm.predict(x[test, ])
table_count(3, y[test], y_pre)
Out[ ]:
array([[11., 0., 0.], [ 0., 12., 0.], [ 0., 1., 6.]])