经典机器学习方法（4）—— 感知机

• 参考：《统计学习方法（第二版）》第2章
• 之前介绍 经典机器学习方法（3）—— 多层感知机 时提到过单层感知机，考虑到单层感知机是神经网络和支持向量机的重要基础，专门写一篇文章加以介绍

---

1. 感知机

1.1 感知机模型

• 感知机（perceptron）是一种二分类的线性分类模型，其输入为样本特征向量，输出为实例类别（取 +1 和 -1 二值），有以下特点
  1. 仅适用于二分类问题，仅有一个输出
  2. 仅适用于样本线性可分的情况（因此无法处理异或问题）
  3. 属于判别模型
• 感知机：假设输入空间是 $\mathcal{X}\subseteq {\mathbb{R}}^n$，输出空间为 Y = { + 1 , − 1 } \mathcal{Y}=\{+1,-1\} Y={+1,−1}，输入 $xx\in \mathcal{X}$ 表示输入实例的特征向量，输出 $y\in \mathcal{Y}$ 表示实例的类别，由以下函数确定
  $$y=f(x)=\text{sign}(w^{\top }x{w}^{\top }x+b)$$ 其中 $ww\in {\mathbb{R}}^n,b\in \mathbb{R}$ 分别是感知机模型的权重参数和偏置参数，$\text{sign}$ 是 step 函数（符号函数），即
  $$\text{sign}(x)=\{\begin{array}{rlr}+1,& & x\ge 0\\ -1,& & x<0\end{array}$$ 感知机的结构图如下所示，可见其可以看做神经网络中一个激活函数为 step 函数的神经元
  
• 从几何角度理解，线性方程 $w^{\top }x{w}^{\top }x+b$ 对应于特征空间 ${RR}^n$ 中的一个超平面，其中 $ww$ 是超平面的法向量，$b$ 是超平面截距，它将特征空间分成两个部分，位于两边的点分属正负两类，因此它也称为分离超平面，如下所示
  
  图中绘制了二维的特征空间（这里和输入空间一致），分离超平面和特征平面正交。这里的数据集是线性可分的，即对于所有 $y_i=-1$ 的样本，有 $w^{\top }{x}_i{w}^{\top }{x}_i+b>0$；对于所有 $y_i=+1$ 的样本，有 $w^{\top }{x}_i{w}^{\top }{x}_i+b<0$。感知机模型的训练目标是找到一个能完全正确划分正负样本的超平面

1.2 感知机的学习策略

• 感知机的学习策略基于经验风险最小化，设计损失时一个自然的选择是误分类点个数，但是这样的损失不是关于 $ww,b$ 的连续可导函数，不易优化；因此我们转而使用所有误分类点到分离超平面的总距离作为损失
• 由于点到直线距离公式为
  $$\frac{1}{||ww||}|w^{\top }x{w}^{\top }x+b|$$ 且有
  $$y_i=\{\begin{array}{rlr}+1,& & w^{\top }x{w}^{\top }x+b\ge 0\\ -1,& & w^{\top }x{w}^{\top }x+b<0\end{array}$$于是误分类点有 $-y_i(w^{\top }{x}_i{w}^{\top }{x}_i+b)>0$，这样，假设超平面 $S$ 的误分类点集合为 $M$，则所有误分类点到 $S$ 的总距离为
  $$-\frac{1}{||ww||}\sum _{{\mathbf{x}}_{\mathbf{i}}\in M}{y}_i(w^{\top }{x}_i{w}^{\top }{x}_i+b)$$ 去掉系数 $\frac{1}{||ww||}$ 即得到感知机的损失函数为
  $$L(ww,b)=-\sum _{{\mathbf{x}}_{\mathbf{i}}\in M}{y}_i(w^{\top }{x}_i{w}^{\top }{x}_i+b)$$ 显然这个损失对 $ww,b$ 都是连续可导的，感知机的学习策略就是在假设空间中选取使得上述损失最小的 $ww,b$ 所对应的感知机模型

1.3 感知机的学习算法

1.3.1 原始形式

• 使用随机梯度下降法进行优化，具体而言，先任意选取一个超平面 ${ww}_0,b_0$（通常取 $0,0$），然后每次随机选取一个误分类点 $({xx}_i,y_i)$，如下执行梯度下降
  $$\begin{array}{rl}ww& \leftarrow ww+\eta y_i{x}_i{x}_i\\ b& \leftarrow b+\eta y_i\end{array}$$ 其中 $\eta \in (0,1]$ 是学习率

• 以上流程称为感知机训练的原始形式，给出其伪代码
  
  直观地看，每轮迭代中我们选择一个误分类点，调整 $ww,b$ 的值使分类平面向它靠近，从而减少该误分类点距离分离超平面的距离，直到超平面越过该点使其被正确分类

• 可以证明

  1. 对于线性可分数据集，感知机学习算法的原始形式一定收敛（一定能找到完全划分训练集的超平面）
  2. 感知机算法存在许多解，初值 ${ww}_0,b_0$ 的选择，以及随机梯度下降时选择误分类点的顺序都会影响最终得到的超平面

1.3.2 对偶形式

• 对偶形式的基本想法是：将 $ww$ 和 $b$ 表示为样本特征 ${xx}_i$ 和标记 $y_i$ 线性组合的形式，通过求解其系数而求得 $ww$ 和 $b$，对于感知机而言其原始形式和对偶形式是等价的，下面我们将原始形式转换为对偶形式

  1. 假设处理原始形式问题时设 ${ww}_0=00,b_0=0$，并且假设分类平面在对样本 ${xx}_i,y_i$ 经过 $n_i$ 次更新后移动到使其正确分类的位置，则 $ww,b$ 关于这个样本的增量分别为 $n_i\eta y_i{x}_i$ 和 $n_i\eta y_i$
  2. 设 $n_i\eta ={\alpha }_i$，原始形式最后求得的 $ww,b$ 可以表示为
     $$\begin{array}{rl}ww& =\sum _{i=1}^N{\alpha }_i{y}_i{x}_i\\ b& =\sum _{i=1}^N{\alpha }_i{y}_i\end{array}$$ 这里 ${\alpha }_i=n_i\eta \ge 0$，当 $\eta =1$ 时，它表示第 $i$ 个样本由于误分类而进行更新的次数。样本点更新次数越多，意味着其与分离超平面越近，也就越难正确分类，这样的样本对于学习结果影响最大，从公式来看就是 $\alpha$_i 越大，其对应的样本在最后得到的参数值中的占比也越大

• 下面给出对偶形式算法的伪代码
  
  注意到这里训练样本仅以内积形式出现，为了提高效率，可以先将样本间的内积都计算出来并以矩阵形式存储

  Note：这个对偶形式其实很少用，我觉得《统计学习方法（第二版）》里介绍它主要是为了引出后面的 SVM 章节。其实对比一下感知机和 SVM，会发现二者都是在以最大间隔为优化标准做分类任务，只是 SVM 中多了约束项，所以求解 SVM 时常常用拉格朗日乘数法转换为对偶问题，像感知机这种数据线性可分的情况，SVM 的原始问题和对偶问题可以同时取到最优解。事实上感知机就是 SVM 算法的重要基础，也可以类似 SVM 的处理过程那样得到它的对偶求解形式

• 和感知机学习算法的原始形式一样，对偶形式的算法也对线性可分数据集收敛，并且存在多个解

2. 实现感知机

• 用 pytorch 实现感知机训练，数据使用《统计学习方法（第二版）》例 2.1 中的数据
  

• 代码如下，可以直接复制进 vscode 运行

  import numpy as np
  import torch
  import random
  import matplotlib.pyplot as plt
  
  def data_iter(batch_size, features, labels):
      num_examples = len(features)
      indices = list(range(num_examples))
      random.shuffle(indices) # 打乱一下，样本的读取顺序是随机的
      
      # 使用 yield 关键字将此函数转为迭代器，每次访问取 batch_size 数据返回
      for i in range(0, num_examples, batch_size):
          j = torch.LongTensor(indices[i: min(i + batch_size, num_examples)]) # 最后一次可能不足一个batch
          yield  features.index_select(0, j)[0], labels.index_select(0, j)[0]
  
  # 损失函数
  def steploss(y, x, w, b): 
      return -y*(torch.dot(x, w) + b)
  
  # 模型
  def perceptron(x, w, b):  
      return 1 if (torch.dot(x, w) + b).item() >= 0 else -1
  
  # 优化方法使用小批量梯度下降
  def mbgd(params, lr, batch_size):
      for param in params:
          param.data -= lr * param.grad / batch_size # 注意这里更改param时用的param.data
  
  # 训练
  def train(net, features, labels, loss, param_w, param_b, batch_size=1, lr=1):
      w, b = param_w, param_b
      done = False
      while not done:
          for X, y in data_iter(batch_size, features, labels):
              _ = net(X, w, b)
              l = loss(y, X, w, b).sum()
              
              if w.grad != None:
                  w.grad.data.zero_()
              if b.grad != None:
                  b.grad.data.zero_()
              
              if y*(torch.dot(X, w) + b) <= 0:
                  l.backward()
                  mbgd((w, b), lr, batch_size)  
  
                  print('误分类点:', X, end='; ')
                  print('更新后参数: w={}, b={}'.format(w.tolist(), b.tolist()))
          
          # 检查是否全部分类正确了
          done = True
          for X, y in data_iter(batch_size, features, labels):
              if y*(torch.dot(X, w) + b) <= 0:
                  done = False
                  break
  
  if __name__ == '__main__':
      # 参数初始化
      w = torch.tensor([0,0], dtype=torch.float32, requires_grad=True)
      b = torch.zeros(1, dtype=torch.float32, requires_grad=True)
  
      # 训练样本 & 标记
      features = torch.tensor([[3,3],[4,3],[1,1]],dtype=torch.float32)
      labels = torch.tensor([1,1,-1],dtype=torch.float32).view(3,1)
  
      # 训练
      train(perceptron, features, labels, steploss, w, b, 1, 1)
  
      # 绘制分离超平面
      plt.figure(figsize=(6, 6))
      x1 = np.linspace(0, 5, 500)
      if w[1].item() != 0:
          x2 = (-x1*w[0].item()-b.item())/w[1].item()
      else:
          x2 = np.linspace(0, 5, 500)
          x1 = torch.full((500,), -b.item()/w[0].item())
  
      positive = torch.masked_select(features, labels==1).view(-1,2)
      negative = torch.masked_select(features, labels==-1).view(-1,2)
  
      plt.scatter(x1, x2, s=1,alpha=1,cmap="r") 
      plt.scatter(positive[:,0], positive[:,1], s=50,alpha=1,cmap="b") 
      plt.scatter(negative[:,0], negative[:,1], s=50,alpha=1,cmap="b") 
  
      plt.show()
  
  '''
  误分类点: tensor([3., 3.]); 更新后参数: w=[3.0, 3.0], b=[1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 2.0], b=[0.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[1.0, 1.0], b=[-1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[0.0, 0.0], b=[-2.0]
  误分类点: tensor([4., 3.]); 更新后参数: w=[4.0, 3.0], b=[-1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[3.0, 2.0], b=[-2.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 1.0], b=[-3.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[1.0, 0.0], b=[-4.0]
  误分类点: tensor([3., 3.]); 更新后参数: w=[4.0, 3.0], b=[-3.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[3.0, 2.0], b=[-4.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 1.0], b=[-5.0]
  '''
  '''
  误分类点: tensor([3., 3.]); 更新后参数: w=[3.0, 3.0], b=[1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 2.0], b=[0.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[1.0, 1.0], b=[-1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[0.0, 0.0], b=[-2.0]
  误分类点: tensor([3., 3.]); 更新后参数: w=[3.0, 3.0], b=[-1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 2.0], b=[-2.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[1.0, 1.0], b=[-3.0]
  '''
  '''
  误分类点: tensor([4., 3.]); 更新后参数: w=[4.0, 3.0], b=[1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[3.0, 2.0], b=[0.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[2.0, 1.0], b=[-1.0]
  误分类点: tensor([1., 1.]); 更新后参数: w=[1.0, 0.0], b=[-2.0]
  '''
  1
  2
  3
  4
  5
  6
  7
  8
  9
  10
  11
  12
  13
  14
  15
  16
  17
  18
  19
  20
  21
  22
  23
  24
  25
  26
  27
  28
  29
  30
  31
  32
  33
  34
  35
  36
  37
  38
  39
  40
  41
  42
  43
  44
  45
  46
  47
  48
  49
  50
  51
  52
  53
  54
  55
  56
  57
  58
  59
  60
  61
  62
  63
  64
  65
  66
  67
  68
  69
  70
  71
  72
  73
  74
  75
  76
  77
  78
  79
  80
  81
  82
  83
  84
  85
  86
  87
  88
  89
  90
  91
  92
  93
  94
  95
  96
  97
  98
  99
  100
  101
  102
  103
  104
  105
  106
  107
  108
  109
  110
  111
  112
  113
  114

  这里取数据是随机的，多次运行即可看到训练时使用样本顺序不同对结果造成的影响。最后的三段注释过程分别对应下面三个结果
  
  另外，这里参数都初始化为 0，修改初始参数训练结果也会有所不同