Implementiert in Python PRML Kapitel 7 Nichtlineare SVM

Lassen Sie uns SVM in Mustererkennung und maschinellem Lernen implementieren, Kapitel 7, "Sparse Kernel Machines". Bei der Wiedergabe von Abbildung 7.2 scheint der Gaugeian-Kernel als Kernel ausgewählt zu sein, daher folgt auch diese Implementierung.

Wie üblich überlasse ich PRML und Hajipata die Erklärung von SVM selbst, z. B. wie die Grenzfläche zu finden ist, und zeige kurz nur den Ablauf der Implementierung.

Grober Ablauf der Implementierung

(1) Lösen Sie die doppelte Darstellung (7.10) der Lagrange-Funktion mit der quadratischen Programmiermethode.

\tilde{L}(a) = \sum_{n=1}^N a_n - \frac{1}{2} \sum_{n=1}^N \sum_{m=1}^N  a_na_m t_n t_m k({\bf x_n}, {\bf x_m})    (7.10)

Es ist unvernünftig, die nichtlineare Programmiermethode selbst zu implementieren! Also habe ich den Solver leise benutzt. Die kostenlose und beliebte Bibliothek cvxopt ist gut, daher habe ich sie installiert.

Wie oben erwähnt, ist es schwer zu erkennen, ob es noch (7.10) ist, mit cvxopt gemäß dem Dokument zu modellieren, also Formeltransformation.

 - \tilde{L}(a) = \frac{1}{2} {\bf A}^{\mathrm{T}} ({\bf T}^{\mathrm{T}} k({\bf x_n}, {\bf x_m}) {\bf T}){\bf A} - {\bf A} (7.10)'

Notational ist es umständlich, aber wenn Sie es beim Vergleichen von (7.10) 'und cvxopt modellieren können, können Sie den Lagrange-Multiplikator $ a $ erhalten, indem Sie cvxopt es lösen lassen.

(2) Bei SVM ist nur der Trägervektor an der Bestimmung der Grenzfläche beteiligt, die durch die KKT-Bedingung eingegrenzt wird. Irgendwie habe ich noch keinen Hunger, aber der Vektor auf der Schnittstelle ist $ t_ny ({\ bf x}) - 1 = 0 $, also muss ab (7.16) $ a> 0 $ verwendet werden. Interpretiert als. Verwenden Sie diese Option, um den Unterstützungsvektor zu extrahieren.

Finden Sie in ③ $ b $ in (7.18) unter Verwendung des erhaltenen Unterstützungsvektors.

④ Ersetzen Sie den in (7.13) erhaltenen Wert, um die Grenzfläche zu erhalten.

Code

Dieses Mal hatte ich große Probleme, die Figur zu zeichnen. Ich konnte die Grenzfläche nicht gut darstellen, daher bezog ich mich auf [diesen Artikel] des Hilfsmittels (http://aidiary.hatenablog.com/entry/20100502/1272804952). Vielen Dank.

import math
import matplotlib.pyplot as plt
from scipy import stats
from scipy.stats.kde import gaussian_kde
from pylab import *
import numpy as np
import random
import cvxopt
from cvxopt import matrix, solvers
%matplotlib inline

def kernel(x, y):
    sigma = 5.0
    return np.exp(-norm(x-y)**2 / (2 * (sigma ** 2)))

#(7.10)' (Quadratic Programming)
def L(t, X, N):
    K = np.zeros((N, N))
    for i in range(N):
        for j in range(N):
            K[i, j] = t[i] * t[j] * kernel(X[i], X[j])
    Q = matrix(K)
    p = matrix(-1 * np.ones(N))             
    G = matrix(np.diag([-1.0]*N))       
    h = matrix(np.zeros(N))             
    A = matrix(t, (1,N))                
    b = matrix(0.0)                     
    sol = solvers.qp(Q, p, G, h, A, b)
    a = array(sol['x']).reshape(N)
    return a

#(7.13)
def y_x(a, t, X, N, b, x):
    sum = 0
    for n in range(N):
        sum += a[n] * t[n] * kernel(x, X[n])
    return sum + b

#(7.18)
def b(a, t, X, S):
    sum_A = 0
    for n in S:
        sum_B = 0
        for m in S:
            sum_B += a[m] * t[m] * kernel(X[n], X[m])
        sum_A += (t[n] - sum_B)
    return sum_A/len(S)

if __name__ == "__main__":
    N = 36
    mu_blue = [1,-1]
    cov = [[0.1,0.05], [0.05,0.1]]
    
    x_blue,y_blue = np.random.multivariate_normal(mu_blue, cov, N/2).T
    
    x_red = [0.3, 0.8, 0.9, 0.95, 1.1, 1.3, 1.6, 1.9, 1.75, 1.8, 2.0, 2.1, 2.3, 2.25, 2.4, 2.7, 3.0, 3.2]
    y_red = [-0.2, 0.1, 0.25, 0.14, -0.1, 1.6, 1.2, 0.6, 0.8, -0.6, -0.8, -0.75, 1.2, -1.15, -0.12, -0.3, -0.4, 1.4]
    
    t_blue = np.ones((1, N/2))
    t_red = -1*np.ones((1, N/2))

    blue = vstack((x_blue, y_blue))
    red = vstack((x_red, y_red))

    X = np.concatenate((blue, red), axis=1).T
    t = np.concatenate((t_blue, t_red), axis=1).T
    
    #(7.10)' (Quadratic Programming)
    a = L(t, X, N)

    #Extract Index of support vectors from (7.14) 
    S = []
    for n in range(len(a)):
        if a[n] < 0.0001: continue
        S.append(n)

    #(7.18)
    b = b(a, t, X, S)

    
    #Plot train data sets
    plt.scatter(x_blue,y_blue,color='b',marker='x')
    plt.scatter(x_red,y_red,color='r',marker='x')
    
    # Enphasize suport vectors
    for n in S:
        plt.scatter(X[n,0], X[n,1], color='g', marker='o')
    
    # Plot the decision surface
    X1, X2 = meshgrid(linspace(-10,10,100), linspace(-10,10,100))
    w, h = X1.shape
    X1.resize(X1.size)
    X2.resize(X2.size)
    Z = array([y_x(a, t, X, N, b, array([x1,x2])) for (x1, x2) in zip(X1, X2)])
    X1.resize((w, h))
    X2.resize((w, h))
    Z.resize((w, h))
    CS = contour(X1, X2, Z, [0.0], colors='k', linewidths=1, origin='lower')
    xlim(0, 4)
    ylim(-2, 2)
    title("Figure 7.2")

Ergebnis

Referenz

Nichtlineare SVM