本系列文章面向程序员,希望通过Image Caption Generation,一个有意思的具体任务,深入浅出地介绍深度学习的知识,涉及到很多深度学习流行的模型,如CNN,RNN/LSTM,Attention等。本文为第三篇。

作者:李理,目前就职于环信,即时通讯云平台和全媒体智能客服平台,在环信从事智能客服和智能机器人相关工作,致力于用深度学习来提高智能机器人的性能。 
相关文章: 
从Image Caption Generation理解深度学习(part I) 
从Image Caption Generation理解深度学习(part II)

2.2.5 反向传播算法的推导

前面我们用很简单的几十行python代码基本上完成了一个多层神经网络。但是还差最重要的部分,那就是计算loss function对参数的偏导数,也就是反向传播算法。下面我们来仔细的完成公式的推导,以及接下来会讲怎么用代码来实现。这一部分数学公式多一些,可能很多读者会希望跳过去,不过我还是建议大家仔细的阅读,其实神经网络用到的数学相比svm,bayes network等机器学习算法,已经非常简单了。请读者阅读的时候最好准备一支笔和几张白纸,每一个公式都能推导一下。如果坚持下来,你会觉得其实挺简单的。

(1) feedforward阶段的矩阵参数表示和计算

之前我们讨论的是一个神经元的计算,而在代码里用到的却是矩阵向量乘法。而且细心的读者会发现我们在构造参数矩阵weights的时候,行数和列数分别是后一层的节点数和前一层的节点数。这似乎有点不自然,为什么不反过来呢?看过下面这一部分就会明白了。

图片描述

首先我们熟悉一下第L(因为小写的L和1太像,所以我用大写的L)层的参数w_jk。它表示第L-1层的第k个神经元到第L层的第j个神经元的权重。比如第3层的w_24,参考上面的图,它表示的是第2层的第4个神经元到第3层的第二个神经元。

对bias和激活函数后的结果a也采用类似的记号,如下图所示。

图片描述

b_32表示第2层的第3个神经元的bias,而a_13第3层的第1个神经元的激活。

使用上面的记号,我们就可以计算第L层的第j个神经元的输出a_jl

图片描述

第L层的第j个神经元的输入是L-1层的a_1,a_2,...;对应的权值是w_j1,w_j2,...;bias是b_jL。所以a_jL就是上面的公式,k的范围是从1到第L-1层的神经元的个数。 
为了用矩阵向量乘法来一次计算第L层的所有神经元的输出,我们需要定义第L层的参数矩阵w_l,它的大小是m*n,其中m是第L层的神经元个数;而n则是第L-1层的个数。它的第i行第j列就是我们上面定义的w_jk。此外我们还要定义向量b_l,它的大小是m(也就是第L层神经元的个数),它的第j个元素就是我们上面定义的b_j。 
最后,我们定义element-wise的函数,比如f(x) = x^2,如果输入是一个向量,那么结果是和输入一样大小的向量,它的每个元素是对输入向量的每一个元素应用这个函数的结果。

图片描述

有了上面的定义,我们就可以一次计算出第L层的输出(一个长度为m的向量)

图片描述

下面是对上面这个公式的详细证明(说明): 
我们需要证明的是向量aL的第j个元素就是前面的a_jL

图片描述

此外,为了方便后面的求解,我们把加权累加和也用一个符号z_l来表示。

图片描述

其中,它的第j个元素就是第L层的第j个神经元的加权累加和:

图片描述

这样a_l就可以简单的对z_l的每个元素计算激活函数

图片描述

现在我们再回顾一下feedforward的代码就非常直观了:

def feedforward(self, a):"""Return the output of the network if a is input."""for b, w in zip(self.biases, self.weights):
a = sigmoid(np.dot(w, a)+b)return a

传给函数feedforward的参数a就是输入向量x,第一层就是x,第二层就是第一个隐层,每一层的计算就是非常简单的参数矩阵w_l乘以上一层的激活a_l-1在加上b_l,然后用激活函数计算。

初始化的时候w的大小是 (后一层的神经元个数) * (前一层的神经元个数),再回顾一下初始化参数的代码:

# sizes = [784, 30, 10]def __init__(self, sizes):self.num_layers = len(sizes)
self.sizes = sizes
self.biases = [np.random.randn(y, 1) for y in sizes[1:]]
self.weights = [np.random.randn(y, x)for x, y in zip(sizes[:-1], sizes[1:])]

x, y in zip(sizes[:-1], sizes[1:]) x是第一层到最后倒数第二层,y是第二层到最后一层,比如上面的sizes=[784, 30, 10] 
x是[784, 30], y是[30, 10],注意随机的矩阵是(y,x),所以self.weights是两个矩阵,大小分别是30*78410*30

(2) 关于损失函数C的两个假设 
1. 损失函数是每个训练数据的损失的平均 
也就是C是这样的形式:

图片描述

对于之前我们使用的MSE损失函数,这是满足的。我们使用batch的梯度下降的时候需要求C对参数w的偏导数,因为损失函数是每个训练数据的损失的平均,所以我们只需要求每个数据的偏导数,然后加起来平均就行。这个假设几乎所有的损失函数都是满足的【我是没见过损失函数不满足这个条件】

  1. 损失函数是最后一层输出的函数

图片描述

这个条件几乎常见的损失函数都是这样的,我们之前时候的MSE就是计算最后一层的输出aL和正确的y(one-hot)的均方误差,显然是满足的。

(3) Hadamard product 
这个名字看起来很复杂,其实很简单,就是两个向量elementwise的乘法。看一个例子就清楚了:

图片描述

(4) 反向传播算法(back propagation)的4个公式

回顾一下,我们之前说了,梯度下降其实最核心的问题就是求损失函数对每一个参数的偏导数。那我们就直接一个一个求好了,为什么又要搞出一个反向传播算法呢?其实这个算法在不同的领域被不同的人重复“发现”过很多次,有过很多不同的名字,最本质的应该就是逆向求导(reverse-mode differentiation)或者叫做自动求导(automatic differentiation)。自动求导(AD)是非常通用的一种求偏导数的方法,很早就在流体力学和大气物理等领域使用,反向传播算法可以认为是AD在神经网络中的应用。不过最早发现这个算法的人(是谁最早好像还有点争议)并不是先知道AD可以直接用于神经网络,他发现这个算法是基于错误的反向传播而得到的,所有命名为(错误的)反向传播算法。后面我们会讲到AD,这是一个强大的算法,任何一个函数,你能把它分解成有向无环图的计算图【函数一般都能分解成一些无依赖的最基础的变量的复合函数,因此肯定可以表示成这样一个有向无环图】,然后每个节点都表示一个函数。只要你能求出这个函数在特定点的梯度【也就是这个函数对所以自变量的偏导数】(不需要求解析的偏导数,当然很多情况,这些函数都是能直接求出解析解,然后代入这个特定点就行,但理论上我们是可以用其他方法,比如数值梯度近似来求的),就能自动的计算损失函数对每一个参数的偏导数(也是在这个点的),而且只要反向根据拓扑排序遍历这个图一次就行,非常高效和简单。后面我们会详细的介绍AD。这个方法非常通用,TensorFlow的核心就是AD。使用AD的框架就比较灵活,我想“创造”一种新的网络结构,我又不想【其实更可能是不会】推导出梯度的公式,那么我只需要把我的网络能用这样一个有向无环图表示就行。当然节点必须要能够求出梯度来,一般我们的函数比如矩阵的运算,卷积等等TensorFlow都封装好了——它把它叫做一个op。我们只需要搭积木一样把这个计算图定义出来,TensorFlow就自动的能根据AD计算出损失函数对所有参数的梯度来了。当然如果你要用到一个TensorFlow没有的op,那你就需要根据它的规范实现这个op,一个op最核心的接口就是两个,一个是输入x,求f(x);另一个就是求f在某个x0点的梯度。

不过这里,我们还是沿着神经网络的发展历史,从错误的反向传播角度来理解和推导这个算法。

首先,我们会对每一个神经元比如第L层的第j个,都定义一个错误δ_jL

图片描述

也就是损失函数对z也就是线性累加和的偏导数。为什么定义这样一个东西呢?我们假设在第L层的第j个神经元上有一个精灵(Daemon)

图片描述

当这个神经元得到来自上一次的输入累加计算出z_jL的时候,它会恶作剧的给一点很小的干扰Δz_jL。原来它应该输出的是σ(z_jL),现在变成了σ(z_jL +Δz_jL)。这个微小的变化逐层传播,最终导致损失函数C也发生如下的变化:

图片描述

这个其实就是导数的直觉定义:微小的Δx引起微小的Δy,Δy/Δx约等于导数。 
不过这个精灵是个好精灵,它想帮助我们减少损失。 当

图片描述

大于0的时候,它让Δz_jL小于0,反之当它小于0的时候它让Δz_jL大于0。这样

图片描述

总是小于0 
因此我们的loss就会变小。而其绝对值越大,我们的损失减少的越多。 
当然你会说为什么不能让Δz_jL非常大,这样我们的损失总是减少很多?可惜这个精灵是个数学家,它说如果Δx太大,那么Δy=df/dx *Δx就不准确了。

所以我们可以这样认为:它就是第L层的第j个神经元“引起”的“错误”。如果绝对值大,则它的“责任”也大,它就得多做出一些调整;反之如果它趋近于0,说明它没有什么“责任”,也就不需要做出什么改变。 
因此通过上面的启发,我们定义出δ_jL来。

图片描述

接下来我们逐个介绍反向传播算法的4个公式。

公式1. 第L层(最后一层) 的错误

图片描述

这个公式的第一项,就是损失C对a_jL的导数,它越大,说明C受a_jL的影响也就越大,如果有了错误,第a_jL的“责任”也就越大,错误也就越大。第二项是a_jLz_jL的影响。两者乘起来就是z_jL对最终损失的影响,也就是它的“责任”的大小。

这个公式很好计算,首先第二项就是把z_jL的值(这个在feedforward节点就算出来并存储下来了)代入σ'(x)。如果σ是sigmoid函数,我们前面也推导过它的导数:σ’(x)=σ(x)*(1-σ(x))。第一项当然依赖于损失函数的定义,一般也很好求。比如我们的MSE损失:

图片描述
图片描述

具体的推导我在纸上写了一下,虽然很简单,我们也可以练练手,尤其是对于求和公式的展开,希望大家能熟悉它,以后的推导我可能就不展开求和公式了,你需要知道求和公式里哪些项是和外面的自变量无关的。

图片描述

公式BP1是elementwise的,我们需要变量j来计算每一个δ_jL。我们也可以把它写成向量的形式,以方便利用线性代数库,它们可以一次计算向量或者矩阵,可以用很多技术利用硬件特性来优化(包括GPU,SSE等)速度。

图片描述

右边δ'(z_L)很容易理解,左边的记号可能有些费解,其实我们把∇aC当成一个整体就好了,它是一个向量,第一个元素是∂C/∂a_1L,第二个就是∂C/∂a_2L,… 
如果算上函数C是MSE的话,上面的公式就可以简化成:

图片描述

公式2. 第l层(非最后一层) 的错误

图片描述

等下我们会证明这个公式,不过首先我们来熟悉一下公式。如果我们想“背”下这个公式的话,似乎看起来比第一个BP1要复杂很多 。我们先检查一下矩阵和向量的维度,假设l+1层有m个元素,l层n个。则w_l+1的大小是m*n,转置之后是n*m, δ_l+1的大小是n*1,所以矩阵相乘后是m*1,这和δ_l是一样的,没有问题。

接下来我们仔细观察一下BP2这个公式,首先第二项σ'(z_l)和前面的含义一样,代表a_l对于z_l的变化率。 
而第一项复杂一点,我们知道第l层的第j个神经元会影响第l+1层的所有神经元,从而也影响最终的损失C。这个公式直接给了一个矩阵向量的形式,看起来不清楚,所以我在草稿纸上展开了:

图片描述

最终第L层的第j个神经元的损失就是如下公式:

图片描述

这下应该就比较清楚了,第l层的第j个神经元的损失,就是把l+1层的损失“反向传播”回来,当然要带上权重,权重越大,“责任”也就越大。

如果要“背”出这个公式也没有那么复杂了,先不看σ'(z_l),第一项应该是矩阵w_l+1乘以δ_l+1。由于矩阵是m*n,而 
向量δ_l+1是m*1,为了能让矩阵乘法成立,那么就只能把w转置一下,变成n*m,然后就很容易记住这个公式了。

注意,BP2的计算是从后往前的,首先根据BP1,最后一层的δ_L我们已经算出来了,因此可以向前计算L-1层的δ_L-1, 
有了δ_L-1就能计算δ_L-2,…,最终能算出第一个隐层(也就是第2层)δ_1来。

公式3. 损失函数对偏置b的梯度

这前面费了大力气求δ_l,不要忘了我们的最终目标是求损失函数对参数w和b的偏导数,而不是求对中间变量z的偏导数。

因此这个公式就是对b的偏导数。

图片描述

或者写成向量的形式:

图片描述

∂C/∂b就是δ!

公式4. 损失函数对w的梯度

图片描述

或者参考下图写成好记的形式:

图片描述
图片描述

也就是说对于一条边w_jkL∂C/∂w_ij就是这条边射出的点的错误δ乘以进入点的激活。非常好记。

我们把这四个公式再总结一下:

图片描述

(5) 这四个公式的证明

首先是BP1,请参考下图:

图片描述

然后是BP2: 
这里用到了chain rule,其实也非常简单和直观,就是复合函数层层组合。最简单的方法就是用图画出来,比如y最终 
是x的函数,我们要求∂y/∂x,如果y是u,v的函数,然后u,v才是x的函数,那么我们把变量x,y,u,v都画成图上的点,y是u,v的函数,那么我们画上从u和v到y的边,同样,我们画上从x到u和v的边,然后从y到x的每条路径,我们经过的边都是一个偏导数,我们把它累加起来就行【这其实就是后面我们会讲的AD】。因此∂y/∂x=∂y/∂u * ∂u/∂x +∂y/∂v * ∂v/∂x。

图片描述

剩下的BP3和BP4也非常类似,我就不证明了。

反向传播算法 
1. a_1 = 输入向量x 
2. Feedforward 根据公式

图片描述

图片描述

计算z_la_l并存储下来(反向传播时要用的) 
3. 计算最后一层的错误

图片描述

  1. 向前计算所有层的错误

图片描述

  1. 计算损失对所有参数的偏导数

图片描述

图片描述

2.2.6 代码实现反向传播算法

我们已经把公式推导出来了,那怎么用代码实现呢?我们先把代码复制一下,然后说明部分都是作为代码的注释了, 
请仔细阅读。

class Network(object):
  def update_mini_batch(self, mini_batch, eta):
    # mini_batch是batch大小,eta是learning rate

    nabla_b = [np.zeros(b.shape) for b in self.biases]    # 构造和self.biases一样大小的向量,比如前面的例子 sizes=[784,30,10],则
    # nabla_b是两个向量,大小分别是30和10

    nabla_w = [np.zeros(w.shape) for w in self.weights]    # 构造和self.weights一样大小的矩阵,比如前面的例子 sizes=[784,30,10],则
    # nabla_w是两个矩阵,大小分别是30*784和10*30

    for x, y in mini_batch: #对于每个训练样本x和y
 delta_nabla_b, delta_nabla_w = self.backprop(x, y) # 用backprop函数计算损失函数对每一个参数的偏导数。
 # backprop函数下面会详细讲解

 nabla_b = [nb+dnb for nb, dnb in zip(nabla_b, delta_nabla_b)] # 把返回的对b偏导数累加到nabla_b中

 nabla_w = [nw+dnw for nw, dnw in zip(nabla_w, delta_nabla_w)] # 把返回的对w的偏导数累加到nabla_w中

    self.weights = [w-(eta/len(mini_batch))*nw  for w, nw in zip(self.weights, nabla_w)]    # 计算完一个batch后更新参数w

    self.biases = [b-(eta/len(mini_batch))*nb for b, nb in zip(self.biases, nabla_b)]    # 更新b... def backprop(self, x, y):
    # 输入是x和y,返回损失函数C对每个参数w和b的偏导数
    # 返回的格式是两个元组,第一个是b的偏导数,第二个是w的。
    nabla_b = [np.zeros(b.shape) for b in self.biases]    # 构造和self.biases一样大小的向量,比如前面的例子 sizes=[784,30,10],则
    # nabla_b是两个向量,大小分别是30和10

    nabla_w = [np.zeros(w.shape) for w in self.weights]    # 构造和self.weights一样大小的矩阵,比如前面的例子 sizes=[784,30,10],则
    # nabla_w是两个矩阵,大小分别是30*784和10*30

    # feedforward
    activation = x
    activations = [x] # 用一个list保存所有层的激活,下面backward会有用的
    zs = [] # 同样的用一个list保存所有层的加权累加和z,下面也会用到。

    #下面这段代码在feedward也有,不过那里是用来predict用的不需要保存zs和activations
    for b, w in zip(self.biases, self.weights):
 z = np.dot(w, activation)+b
 zs.append(z)
 activation = sigmoid(z)
 activations.append(activation)   # backward pass
  #1. 首先计算最后一层的错误delta,根据公式BP1,它是损失函数对a_L的梯度乘以σ'(z_L)
    #  sigmoid_prime就是σ'(z_L),而∂C/∂a_L就是函数cost_derivative,对于MSE的损失函数,
    #  它就是最后一层的激活activations[-1] - y
    delta = self.cost_derivative(activations[-1], y) * \
 sigmoid_prime(zs[-1])    # 2. 根据公式BP3,损失对b的偏导数就是delta
    nabla_b[-1] = delta    # 3. 根据公式BP4,损失对w的偏导数时delta_out * activation_in
    #  注意,我们的公式BP4是elementwise的,我们需要写成矩阵向量的形式
    #  那怎么写呢?我们只需要关心矩阵的大小就行了。
    #  假设最后一层有m(10)个神经元,前一层有n(30)个,
    #  则delta是10*1, 倒数第二层的激活activations[-2]是30*1
    #  我们想求的最后一层的参数nabla_w[-1]是10*30,那么为了能够正确的矩阵乘法,
    #  只要一种可能就是 delta 乘以 activations[-2]的转置,其实也就是向量delta和activations[-2]的外积
    nabla_w[-1] = np.dot(delta, activations[-2].transpose())    # 接下来从倒数第二层一直往前计算delta,同时也把对w和b的偏导数求出来。
    # 这里用到一个比较小的trick就是python的下标是支持负数的,-1表示最后一个元素,-2是倒数第二个
    # l表示倒数第l层,2就表示倒数第2层,num_layers - 1就表示顺数第2层(也就是第1个隐层)
    # 比如我们的例子:sizes=[784, 30, 10],那么l就是从2到3(不包含3),l就只能是2,页就是第1个(也是唯一的一
    # 个)隐层   
    for l in xrange(2, self.num_layers): # 倒数第l层的z
 z = zs[-l] # 计算σ'(z_l)
 sp = sigmoid_prime(z) # 根据BP2,计算delta_l,注意weights[-l+1]表示倒数第l层的下一层
 delta = np.dot(self.weights[-l+1].transpose(), delta) * sp # 同上,根据BP3
 nabla_b[-l] = delta # BP4,矩阵乘法参考前面的说明
 nabla_w[-l] = np.dot(delta, activations[-l-1].transpose())    return (nabla_b, nabla_w)

2.2.7 为什么反向传播算法是一个高效的算法?

分析完代码,我们发现一次backprop函数调用需要feedforward一次,网络有多少边,就有多少次乘法,有多少个点就有多少次加分和激活函数计算(不算第一层输入层)。反向计算也是一样,不过是从后往前。也就是说这是时间复杂度为O(n)的算法。 
如果我们不用反向传播算法,假设我们用梯度的定义计算数值梯度。对于每一个参数wj,

我们都用公式 limit (f(w1, w2, …, wj+Δ wj, …) - f(w1, w2, …, wj, …)/Δwj

f(w1, w2, wj, …)只需要feedforward一次,但是对于每个参数wj,都需要feedforward一层来计算f(w1, w2, …, wj+Δ wj, …),它的时间复杂度是O(n),那么对所有的参数的计算需要O(n^2)的时间复杂度。

假设神经网络有1百万个参数,那么每次需要10^12这个数量级的运算,而反向传播算法只需要10^6,因此这个方法比反向传播算法要慢1百万倍。