机器学习笔记3-1：反向传播算法原理与实现

*注:本博客参考李宏毅老师2020年机器学习课程. 视频链接

---

1 反向传播算法数学基础

反向传播算法是用于对神经网络中的各个网络参数计算偏导值的一种算法，其核心是链式求导法则。
注：本节涉及到较复杂的数学计算，了解思想即可。

(图一)

在图一中，画出了一个神经网络的前两层结构的部分神经元。图中，$x_i$表示网络的输入。$w_{ij}$表示第i层神经元的参数，j仅用作区分不同的参数。$z_{ij}$表示第i层神经元经过激活函数之前的值，$a_{ij}$表示第i层第j个神经元的输出。$C$表示网络计算得到的损失函数值（Cost）.由此可知：
$$z_{11}=\sum _{j=1}^n{w}_{1j}{x}_j\text{\hspace{1em}(1)}$$
$$a=\sigma \left(z\right)\text{\hspace{1em}(2)}$$
其中$\sigma$表示激活函数。

如果要计算$\frac{\partial C}{\partial w_{11}}$，根据链式求导法则，有：
$$\frac{\partial C}{\partial w_{11}}=\frac{\partial z_{11}}{\partial w_{11}}\frac{\partial C}{\partial z_{11}}\text{\hspace{1em}(3)}$$
在等式3的右边，计算前一项的过程称为前向传播，计算后一项的过程称为反向传播。

1.1 前向传播

前向传播即计算$\frac{\partial z_{11}}{\partial w_{11}}$，根据式1可知，$\frac{\partial z_{11}}{\partial w_{11}}=x_1$，实际上：
$$\frac{\partial z_{2i}}{\partial w_{2i}}=a_{11},i=1,2,...m\text{\hspace{1em}(4)}$$
照此推理，只要知道上一层神经元的输出，就可以求得这一层神经元的输入z对这一层神经元的参数w的偏导数值。在顺序计算损失函数值得过程中，所有的$\frac{\partial z}{\partial w}$都可以计算完毕，因此称为前向传播。

1.2 反向传播

$\frac{\partial C}{\partial z_{11}}$的计算要稍微复杂一些，首先可以根据式2，得：
$$\frac{\partial C}{\partial z_{11}}=\frac{\partial C}{\partial a_{11}}{\sigma }^{\prime }(z_{11})\text{\hspace{1em}(5)}$$
由于$z_{11}$是一个已知数值，因此${\sigma }^{\prime }(z_{11})$是一个常数。对于上式等号右侧得第一项，从图1中可以看出，$a_{11}$得变化首先会改变$z_{2j},j=1,2,...,m$的值变化，最终导致C的值变化，那么：
$$\frac{\partial C}{\partial a_{11}}=\sum _{i=1}^m\frac{\partial C}{\partial z_{2i}}{w}_{2i}\text{\hspace{1em}(6)}$$
将式5和式6联立化简得：
$$\frac{\partial C}{\partial z_{11}}={\sigma }^{\prime }(z_{11})\sum _{i=1}^m\frac{\partial C}{\partial z_{2i}}{w}_{2i}\text{\hspace{1em}(7)}$$
上式中的${\sigma }^{\prime }(z_{11})$和$w_{2i}$的值均为已知，于是不难发现，式7的定义是递归的，也就是说，要计算损失函数值C对某一层神经元的输入$z_i$的偏导值，必须先知道C对下一层神经元的输入$z_{i+1}$的偏导值。那么自然会想到从最后一层开始，向前计算。

(图2)

在图2中给出了神经网络的最后一层的结构。
$$y_i=\frac{1}{1+e^{-z_i}},i=1,2,...,t\text{\hspace{1em}(8)}$$
当使用交叉熵作为损失函数时：
$$C=-\sum _{i=1}^t(\widehat{y_i}\ln y_i+(1-\widehat{y_i})\ln (1-y_i))\text{\hspace{1em}(9)}$$
则
$$\frac{\partial C}{\partial z_1}=\frac{\partial C}{\partial y_1}\frac{\partial y_1}{\partial z_1}=-(\frac{\widehat{y_1}}{y_1}-\frac{1-\widehat{y_1}}{1-y_1}){\sigma }^{\prime }(y_1)=y_1-\widehat{y_1}\text{\hspace{1em}(10)}$$
综上，结合式3、4、7、10，便可以求得网络中任何一个神经元对于其参数w的偏导值。对于式3，我们可以写出它的表达式：
$$\frac{\partial C}{\partial w_{11}}=x_1{\sigma }^{\prime }(z_{11})\sum _{i=1}^m\frac{\partial C}{\partial z_{2i}}{w}_{2i}\text{\hspace{1em}(11)}$$

2 矩阵运算

下面我们将上述计算过程写成矩阵运算的形式，这种形式更适合在实际的编码过程中使用各种矩阵运算库进行加速运算（如numpy、tensorflow、PyTorch等）。

• 首先回到图1，我们将神经网络的输入写成$X_1\in {\mathbb{R}}^{n\times 1}$，表示它是一个n行1列的矩阵，对应输入具有n个维度。
• 对于神经网络的第一层，假设第一层有k个神经元，则每一个神经元都应该具有n个参数w。因此将第一层神经元的所有参数w写成矩阵$W_1\in {\mathbb{R}}^{k\times n}$，表示它是一个k行n列的矩阵，其中每一行表示同一个神经元对输入$X$的不同维度的对应参数，每一列表示不同神经元对输入$X$的相同维度的对应参数，下标1代表它是第一个隐藏层的参数。
• 同理，第一层的每一个神经元还有一个参数b，将该曾所有的b写成矩阵为$B_1\in {\mathbb{R}}^{k\times 1}$。
• 第一个隐藏层的激活函数前的输出$Z$也是一个矩阵，表示为$Z_1\in {\mathbb{R}}^{k\times 1}$，且$Z_1=W_1{X}_1+B_1$。该层的输出表示为$A_1\in {\mathbb{R}}^{k\times 1}$，且$A_1=\sigma (Z_1)$，同时该层的输出也是第二层的输入，因此$A_1=X_2$
• 其余各层的运算与上述过程类似，不再赘述。
• 对于神经网络的最后一层，如图2，将得到整个网络的输出$Y$，与之前不同的是，将激活函数$\sigma$换成了Softmax函数，假设网络输出包含t维，那么$Y\in {\mathbb{R}}^{t\times 1}$，$Y$的各个维度的值的计算遵循公式8。

因此，矩阵运算计算第二层输入的过程可以表示如下：

(图3)

其中，dot表示矩阵乘法。在前向传播的过程中，不断重复以上过程，就可以得到所有的层的输入$X$。该过程如下图4，前向传播的任务就是完成矩阵$X$的计算：

(图4)

将损失函数$C$对第i层(非最后一层)的参数$W_i$求导的计算式11加以扩展，改写为矩阵运算形式如下：
$$\frac{\partial C}{\partial W_i}=[{\sigma }^{\prime }(Z_i)\ast (W_{i+1}^T\frac{\partial C}{\partial Z_{i+1}})]X_i^T\text{\hspace{1em}(12)}$$
上式中，$W_i$为第i层的网络权重，$X_i$为第i层神经元的输入，$Z_i=W_i{X}_i+B_i$为第i层神经元激活前的输出值，$\ast$表示元素乘法（即对应位置元素相乘），$T$表示矩阵转置。
假设该网络的第i层有m个神经元，第i+1层有n个神经元，第i-1层有p个神经元,那么式11中各变量的维度分别如下表：

变量	维度	变量	维度
$X_i$	(p,1)	$Z_i$	(m,1)
$Z_{i+1}$	(n,1)	$W_i$	(m,p)
$\frac{\partial C}{\partial W_i}$	(m,p)	$\frac{\partial C}{\partial Z_{i+1}}$	(n,1)
$W_{i+1}^T\frac{\partial C}{\partial Z_{i+1}}$	(m,1)	$W_{i+1}^T$	(m,n)

上述讨论中我们没有考虑对B的偏导数，实际上，如果将式3中等号右边第一项替换为$\frac{\partial z_{11}}{\partial b_{11}}$即可得出对B的偏导，并且这一项始终为1，因此：
$$\frac{\partial C}{\partial B_i}={\sigma }^{\prime }(Z_i)\ast (W_{i+1}^T\frac{\partial C}{\partial Z_{i+1}})\text{\hspace{1em}(13)}$$

3 模型实现

根据上述数学推理，设计并实现一个神经网络已经变得可能。下面我们开始着手设计实现这个数据结构。

3.1 模型结构设计

首先明确我们要设计的神经网络的结构，假设该模型有m层，m最小为3，因为除了输入层与输出层之外，网络应该包含至少一个隐藏层。
每个隐藏层都有自己的W和B。若第$i\in [1,m-1]$层有$n_i$个神经元，那么其参数满足：$W_i\in {\mathbb{R}}^{n_i\times n_{i-1}},B_i\in {\mathbb{R}}^{n_i\times 1}$，整个网络的参数一共2(m-2)个矩阵（其中m-2个矩阵为权重W的矩阵，m-2个矩阵为偏置B的矩阵）。网络结构定义如下：

import numpy as np
class Network():
    def __init__(self, *layers):
        # 至少包含输入层，输出层，和一个隐藏层，所以至少有三层
        assert(len(layers) > 2)
        # W的维度为(当前层神经元数，上一层神经元数)
        self.W = [np.random.randn(y, x)
                  for x, y in zip(layers[:-1], layers[1:])]
        # B的维度为(当前层神经元数，1)
        self.B = [np.random.randn(y, 1)
                  for x, y in zip(layers[:-1], layers[1:])]

根据式12可知，在前向传播的过程中，需要记录$X_i$和$Z_i$的值，以方便计算偏导值。

def forward(network, X):
    Zs = []
    Xs = []
    for W, B in zip(network.W, network.B):
        Xs.append(X)
        Z = W.dot(X)+B
        X = sigmoid(Z)
        Zs.append(Z)
    Y = sigmoid(Z)
    return Xs, Zs, Y
# sigmoid函数
def sigmoid(Z):
    A = 1/(1+np.exp(-Z))
    return A

反向传播的计算$\frac{\partial C}{\partial Z_i}$的表达式就是式12等式右边第一项，计算过程中需要保存所有的$\frac{\partial C}{\partial Z_i}$。

def backward(network, pcpz, Zs):
    P_z = [pcpz]
    for W, Z in zip(network.W[::-1][:-1], Zs[::-1][1:]):
        pcpz = sigma_(Z)*(W.T.dot(pcpz))
        P_z.append(pcpz)
    P_z.reverse()
    return P_z
# sigmoid函数的导数
def sigma_(Z):
    return sigmoid(Z)*(1-sigmoid(Z))

在开始反向传播之前，首先需要按照式10计算$\frac{\partial C}{\partial Z_i}$的最后一项：

def get_last_pcpz(Y0, Y1):
    return Y1-Y0

下面给出训练过程的代码：

def start_train(train_x, train_y,epoch,lr,*hidden_layers):
    nn = Network(train_x[0].shape[0], *hidden_layers, train_y[0].shape[0])
    loss = [] #记录损失函数值的变化
    for i in range(epoch):
        # 初始化W和B在当轮迭代中的变化量
        delta_W = [np.zeros_like(w) for w in nn.W]
        delta_B = [np.zeros_like(b) for b in nn.B]
        l = 0
        # 对训练数据中的每一个样本进行计算
        for X, Y in zip(train_x, train_y):
            Xs, Zs, Y1 = forward(nn, X) #前向传播
            pcpz = get_last_pcpz(Y, Y1) #计算最后一项
            P_z = backward(nn, pcpz, Zs) #反向传播
            l += np.abs(Y-Y1).sum() #计算当前样本的损失函数值
            # 更新当前轮次的W和B由于该样本所带来的变化量
            delta_B = [db+pz for db, pz in zip(delta_B, P_z)] 
            delta_W = [dw+pz.dot(x.T) for dw, pz, x in zip(delta_B, P_z, Xs)]
        loss.append(l)
        # 更新模型参数
        nn.W = [w-dw*lr/len(train_x) for w, dw in zip(nn.W, delta_W)]
        nn.B = [b-db*lr/len(train_x) for b, db in zip(nn.B, delta_B)]
    return nn, loss

3.2 模型测试

使用一组简单的数据对所实现的模型进行测试。该组数据由2个样本构成，其值分别为(1,1)和(0,0)，我们将前者归为类别1，后者归为类别2，对于类别1，我们用(1,0)表示，类别2，我们用(0,1)表示，下面对着一组简单的数据进行训练，并给出训练过程中loss的变化情况。

import matplotlib.pyplot as plt
train_x = [np.ones((2, 1)), np.zeros((2, 1))]
train_y = [np.asarray([[1], [0]]), np.asarray([[0], [1]])]
nn, loss = start_train(train_x, train_y,3000,1,3,4)
for X, Y in zip(train_x, train_y):
    Xs, Zs, Y1 = forward(nn, X)
    print("input:{} output: {} class:{}".format(X.tolist(),Y1.tolist(),np.argmax(Y1)+1))
plt.plot(loss)
plt.show()

input:[[1.0], [1.0]] output: [[0.9989013915763175], [0.00127671782739844]] class:1
input:[[0.0], [0.0]] output: [[0.0008225495145913559], [0.9990508703257248]] class:2