In [1]:

    

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline


    

In [2]:

    

def GenerateDummyData(nPos=100, nNeg=100):
    N = nPos + nNeg
    Xpos = np.random.normal(+6, 1, (nPos, 2))
    Xneg = np.random.normal(+3, 1, (nNeg, 2))
    Ypos =  np.ones(nPos)
    Yneg = -np.ones(nNeg)
    X = np.r_[Xpos, Xneg]
    Y = np.r_[Ypos, Yneg]
    perm = np.random.permutation(N)
    X = X[perm,:]
    Y = Y[perm]
    return X, Y


    

In [3]:

    

def Primal(X, Y, lam, w):
    r = lam * w.dot(w) / 2
    l = 1 - Y * X.dot(w)
    l = ( l * (l > 0) ).sum() / l.size
    return r + l

def Dual(X, Y, lam, w, alp):
    r = -lam * w.dot(w) / 2
    l = (Y * alp).sum() / alp.size
    return r + l

def SolveSVM(X, Y, lam, T=10):
    N, d = X.shape
    
    alpha = np.zeros(N)
    w = np.zeros(d) # in general cace s.t. alpha != 0 , w = alpha.dot(X) / (lam*N)
    p, d, sv = [ np.zeros(T) for _ in range(3) ] # to save learning history.
    
    for t in range(T):
        perm = np.random.permutation(N)
        for i in perm:
            x, y, alp_old = X[i], Y[i], alpha[i]
            alp_new = y * np.clip( y*alp_old + lam*N*( 1-y*x.dot(w) ) / x.dot(x), 0, 1 )
            alpha[i] = alp_new
            w += (alp_new - alp_old) * x / (lam*N)

        p[t] = Primal(X,Y,lam,w)
        d[t] = Dual(X,Y,lam,w,alpha)
        sv[t] = (alpha != 0).sum()
    
    return w, alpha, p, d, sv


    

In [4]:

    

%%time
np.random.seed(2018)
lam = 1e-3
I = int(3e3) # the number of iteration.

X, Y = GenerateDummyData()
N = Y.size
X = np.c_[X, np.ones(X.shape[0]) ] # Add bias

w,alpha,p,d,sv = SolveSVM(X, Y, lam, I)


    

    




CPU times: user 9.82 s, sys: 4 ms, total: 9.82 s
Wall time: 9.82 s

In [5]:

    

# Compute separating line.
x = X[:,0]
x = np.linspace(x.min(), x.max())
y = -( w[0]*x + w[2] ) / w[1]

# Plot
fig = plt.figure(figsize=(10, 4))
ax1 = plt.subplot2grid((2,3), (0,0), rowspan=2)
ax2 = plt.subplot2grid((2,3), (0,1))
ax3 = plt.subplot2grid((2,3), (1,1))
ax4 = plt.subplot2grid((2,3), (0,2), rowspan=2)

# Plot
ax = ax1
ax.set_title('Data distribution and SVs.')
ax.plot(x, y, '--')
ax.plot(*X[alpha != 0,:2].T, 'go', alpha=0.3)
ax.plot(*X[Y==+1,:2].T, 'r+')
ax.plot(*X[Y==-1,:2].T, 'b+')

ax = ax2
ax.set_title('Objective value')
ax.plot(p, label='primal')
ax.plot(d, label='dual')
ax.set_xscale('log')
ax.set_yscale('log')
ax.legend()

ax = ax3
ax.set_title('Duality gap')
ax.plot(p-d)
ax.set_xscale('log')
ax.set_yscale('log')

ax = ax4
ax.set_title('# of Supporting vectors')
ax.plot(sv)

for ax in [ax1,ax2,ax3,ax4] : ax.grid()
fig.tight_layout()