用Python实现单层感知机

单层感知机的原理

感知机是一种二分类的线性分类模型，一个训练好的感知机会将输入的特征向量转化成+1和-1两类。在训练过程（将一个目标函数最小化）中，会根据输入特征学习到一个超平面，这个超平面由两个参数向量（$w$和$b$）决定，训练过程也就是不断改变这两个参数的过程，直到取值使得目标函数最优。

感知机模型

假设一个数据集（$x,y$），$x$是特征值，$y$是样本标签（取值为+1和-1），在数据集线性可分的情况下，一个训练好的感知机会将输入$x$利用下面的公式来得到一个正确的输出$y$。
$y=f(x)=sign(wx+b)=\{\begin{array}{lr}+1,wx+b\ge 0& \\ -1,wx+b\le 0& \end{array}$
$w$称为权值，$b$称为偏置，$wx+b=0$对应于特征空间$R^n$中的一个超平面，并将这个超平面分成两部分，以二维特征为例，$x^{(1)},x^{(2)}$是样本的两个特征，假设训练数据有22个，正例数据12个，反例数据10个，则输入数据$x$就是$22\times 2$，经过训练之后，就会得到下面这个直线（$\ge 3$维特征就是超平面），这条直线所对应的方程就是$wx+b=0$（这时候，$w$和$b$并不唯一）。

上面这张图对应的是线性可分的情况，其所对应的目标函数是误分类点到超平面的总间隔，这个模型后续的发展是往支持向量机方向上发展，比如支持向量机的软间隔、线性不可分等问题。
上面的$wx+b$还可以表示成下面这幅图，单个神经元的模型，

在这种模型下，$y=f(wx-\theta )$，$-\theta$可以被视为$b$，$f$是神经网络中的激活函数，上面模型中的$sign$函数就是激活函数的一种（实际上神经网络中不会应用这种激活函数，但此处这个函数的功能与激活函数类似）。

感知机的学习策略

对线性可分的数据集，为了找到这样的一个超平面，实际上是为了找到这样一组参数$w$和$b$，我们需要确定一个学习策略，即定义损失函数并将损失函数极小化，一个很自然的想法就是把误分类的个数最小，最好是没有误分类的个数（优化所追求的一个目标）最好，但是这样的一个损失函数并不是参数$w$和$b$的连续可导函数，不容易优化。所以另外定义的损失函数是误分类点到超平面的总距离，假设一个误分类点$x_0$，这个点到超平面$S$的距离为
$\frac{w{x}_0+b}{||w||}$，
这里，$||w||$是$w$的$L_2$范数，假设超平面$S$的误分类点集合为$M$，则目标函数就可以写成
$L(w,b)=\sum _{x_i\in M}\frac{1}{||w||}|w{x}_i+b|$，
我们希望这个总距离越短越好，最好是为0，这样就没有误分类点了，但是这种方式带着绝对值，所以在某些点是不可导的，所以就用下面这种巧妙的方式把目标函数变成连续可导的形式，而且仍然保持与上式等价的形式。
$L(w,b)=\frac{1}{||w||}\sum _{x_i\in M}-y_i(w{x}_i+b)$，
因为一个数据$(x_i,y_i)$分类错误之后，下面这两个式子是等价的（$y$取值为+1和-1的情况下）。
$|w{x}_i+b|=-y_i(w{x}_i+b)>0$，
目标函数最前面的系数对目标函数没有影响，所以最终的目标函数可以写成
$L(w,b)=\sum _{x_i\in M}-y_i(w{x}_i+b)$。

感知机的学习算法

有了这样的目标函数之后，我们来看看怎么才能找到这样的一组参数$(w,b)$来使得目标函数最小，其实这就变成了数学上的一个优化问题，我们可以给$(w,b)$随机赋初始值，然后利用各种梯度下降的算法来找到这样一组参数。下面是采用随机梯度下降算法的一个实现过程。
首先随机选取一个超平面$(w_0,b_0)$，然后用随机梯度下降法不断的极小化目标函数，这里的极小化过程不是一次使$M$中所有的误分类点梯度下降，而是一次随机选取一个误分类点使其梯度下降（其他的梯度下降还有批量梯度下降和小批量梯度下降）。损失函数的梯度由下面给出
$\frac{\partial L(w,b)}{\partial w}=-\sum _{x_i\in M}{y}_i{x}_i$
$\frac{\partial L(w,b)}{\partial b}=-\sum _{x_i\in M}{y}_i$，
随机选择一个误分类点（$x_i,y_i$），对$w,b$进行更新：
$w\leftarrow w+\eta y_i{x}_i$
$b\leftarrow b+\eta y_i$，
$\eta$是步长，在统计学中称为学习率，这样不断的迭代，就可以使得损失函数不断变小，直到为0。算法的流程如下：
（1）选取初值$w_0,b_0$
（2）在训练集中选择数据$(x_i,y_i)$
（3）如果$y_i(w{x}_i+b)\le 0$:
$w\leftarrow w+\eta y_i{x}_i$
$b\leftarrow b+\eta y_i$，
（4）转到（2），直到训练集中没有误分类点。
不同的初始值可能会使得最后得到的$w,b$不一致，但这是对结果没影响的（要说影响其实也是有的，一个超平面如果很靠近两个类别，新的数据进来可能就会误分类了，如果超平面在两类之间的话，至少对新数据的分类正确性会高一点，后面的支持向量机就是靠间隔来取胜的）。

感知机学习算法的对偶形式

至于为什么要采取这种对偶形式，知乎上面有讨论（如何理解感知机学习算法的对偶形式？），在此我简单提一下自己的想法。
对偶形式的基本想法是将$w$和$b$表示为实例$x_i$和标记$y_i$的线性组合的形式，通过求解其系数而求得$w$和$b$，不失一般性，假设初始值$w_0,b_0$均为0，对误分类点$(x_i,y_i)$通过
$w\leftarrow w+\eta y_i{x}_i$
$b\leftarrow b+\eta y_i$，
逐步修改$w,b,$设修改$n$次，则$w,b$关于$(x_i,y_i)$的增量分别是${\alpha }_i{y}_i{x}_i$和${\alpha }_i{y}_i$，这里${\alpha }_i=n_i\eta$。这样，从学习过程不难看出，最后学习到的$w,b$可以分别表示为
$w=\sum _{i=1}^N{\alpha }_i{y}_i{x}_i$，
$b=\sum _{i=1}^N{\alpha }_i{y}_i$，
这里，${\alpha }_i\ge 0,i=1,2,......,N$，当$\eta =1$时，表示第$i$个实例点由于误分而进行更新的次数。实例点更新次数越多，意味着它距离超平面越近，也就越难正确分类，换句话说，这样的实例对学习结果影响最大。我认为这里是由于步长$\eta$的选择的问题，可能是一直在最优点附近跳跃，所以这里选择的是参数${\alpha }_i$来对步长做变化。
下面的算法流程凭直觉是可以印证一部分观点的。
感知机学习算法的对偶形式
输出：$\alpha ,b$；感知机模型$f(x)=sign(\sum _{j=1}^N{\alpha }_j{y}_j{x}_j\cdot x+b)$
（1）$\alpha \leftarrow 0,b\leftarrow 0$
（2）在训练集中选取数据$(x_i,y_i)$
（3）如果$y_i(\sum _{j=1}^N{\alpha }_j{y}_j{x}_j\cdot x+b)\le 0$，
${\alpha }_i\leftarrow {\alpha }_i+\eta$，
$b\leftarrow b+\eta y_i$，
（4）转至（2）直到没有误分类数据
对偶形式中训练实例仅以内积的形式出现，为了方便，可以预先将训练集中实例间的内积形式计算出来并以矩阵的形式存储，这个矩阵就是所谓的Gram矩阵。
$G=[x_i\cdot x_j{]}_{N\times N}$。

单层感知机的Python实现

下面是代码的python实现，比较简单的写了几个数据。
首先是训练过程，会提示输入训练次数，然后利用这几个点找出超平面，最后是对新数据点的预测。
我输入的是训练30次，会得到这样的超平面。感知机学习算法实践这是一篇利用感知机进行手写数字识别的，可以继续深入的了解一下。

import numpy as np
import random
import matplotlib.pyplot as plt
#sign function
def sign(v):
    if v>0:
        return 1
    else:
        return -1
# train function to get weight and bias 
def training():
    train_data1 = [[1,3,1], [2,5,1], [3,8,1], [2,6,1]] #positive sample
    train_data2 = [[3,1,-1], [4,1,-1], [6,2,-1], [7,3,-1]] #negative sample
    train_data = train_data1 + train_data2;
       
    weight = [0,0]
    bias = 0
    learning_rate = 0.1
    
    train_num = int(input("train num:"))
    
    for i in range(train_num):
        train = random.choice(train_data)
        x1,x2,y = train;
        y_predict = sign(weight[0]*x1 + weight[1]*x2 + bias)
        print("train data:x:(%d, %d) y:%d ==>y_predict:%d" %(x1,x2,y,y_predict))
        if y*y_predict<=0:
            weight[0] = weight[0] + learning_rate*y*x1
            weight[1] = weight[1] + learning_rate*y*x2
            bias = bias + learning_rate*y
            print("update weight and bias:")
            print(weight[0], weight[1], bias)
    print("stop training :")
    print(weight[0], weight[1], bias)
    
    #plot the train data and the hyper curve
    plt.plot(np.array(train_data1)[:,0], np.array(train_data1)[:,1],'ro')
    plt.plot(np.array(train_data2)[:,0], np.array(train_data2)[:,1],'bo')
    x_1 = []
    x_2 = []
    for i in range(-10,10):
        x_1.append(i)
        x_2.append((-weight[0]*i-bias)/weight[1])
    plt.plot(x_1,x_2)
    plt.show()
    
    return weight, bias

#test function to predict
def test():
    weight, bias = training()
    while True:
        test_data = []
        data = input("enter q to quit,enter test data (x1, x2):")
        if data == 'q':
            break
        test_data += [int(n) for n in data.split(',')]
        predict = sign(weight[0]*test_data[0] + weight[1]*test_data[1] + bias)
        print("predict==>%d" %predict)

if __name__ == "__main__":
    test()