神经网络中的矩阵求导及反向传播推导

第一部分: 一个简单的两层神经网络的反向传播

下面的代码是来自 pytorch tutorial 的一个 numpy 版本的(激活函数为relu的)两层全连接神经网络的实现, 包括网络的实现、梯度的反向传播计算和权重更新过程：

# -*- coding: utf-8 -*-import numpy as np# N is batch size; D_in is input dimension;# H is hidden dimension; D_out is output dimension.N, D_in, H, D_out = 64, 1000, 100, 10# Create random input and output datax = np.random.randn(N, D_in)y = np.random.randn(N, D_out)# Randomly initialize weightsw1 = np.random.randn(D_in, H)w2 = np.random.randn(H, D_out)learning_rate = 1e-6for t in range(500):    # Forward pass: compute predicted y    h = x.dot(w1)    h_relu = np.maximum(h, 0)    y_pred = h_relu.dot(w2)    # Compute and print loss    loss = np.square(y_pred - y).sum()    print(t, loss)    # Backprop to compute gradients of w1 and w2 with respect to loss    grad_y_pred = 2.0 * (y_pred - y)    grad_w2 = h_relu.T.dot(grad_y_pred)    grad_h_relu = grad_y_pred.dot(w2.T)    grad_h = grad_h_relu.copy()    grad_h[h < 0] = 0    grad_w1 = x.T.dot(grad_h)    # Update weights    w1 -= learning_rate * grad_w1    w2 -= learning_rate * grad_w2

这里我们主要关心其中的 反向传播 过程，核心代码如下：

h = x.dot(w1)h_relu = np.maximum(h, 0)y_pred = h_relu.dot(w2)loss = np.square(y_pred - y).sum()grad_y_pred = 2.0 * (y_pred - y)    # 64 x 10grad_w2 = h_relu.T.dot(grad_y_pred) # 100 x 10grad_h_relu = grad_y_pred.dot(w2.T) # 64 x 100grad_h = grad_h_relu.copy()         # 64 x 100grad_h[h < 0] = 0                   # 64 x 100grad_w1 = x.T.dot(grad_h)           # 1000 x 100

0. 变量关系分析

首先画出依赖图:

各个变量之间的关系如下:

h=h_relu=y_pred=loss=x⋅w1ReLU(h)h_relu⋅w2∑(y_pred−y)2=sse(y_pred)

1. grad_y_pred: 实值函数对矩阵求导

∂loss∂y_pred=∂∑ij(y_predij−yij)2∂y_pred

其中loss是实数, y_pred 是矩阵, 根据实值函数对矩阵求导规则(详见本文第二部分)有:

(∂loss∂y_pred)ij====∂loss∂y_predij∂∑ij(y_predij−yij)2∂y_predij∂(y_predij−yij)2∂y_predij2⋅(y_predij−yij)

故:

grad_y_pred=∂loss∂y_pred=2⋅(y_pred−y)

2. grad_w2: 线性变换的导数

设有 f(Y):Rm×p→R 及线性映射 X↦Y=AX+B:Rn×p→Rm×p, 其中 A∈Rm×n,B∈Rm×p.则:

∇f(AX+B)=AT∇Yf

推导: fxij=∑kl∂f∂ykl∂ykl∂xij, 其中 ∂ykl∂xij=∂∑saksxsl∂xij=akiδlj(δlj 是克罗内克函数,若 l=j 则值为1,否则为0), 将后式带入前式,得: ∂f∂xij=∑kl∂f∂yklakiδlj=∑k∂f∂ykjaki, 即矩阵 AT 的第 i 行和矩阵 ∇Yf 的第j列的内积.

类比以上公式, 则 f=sse, 线性映射为 y_pred=h_relu⋅w2, 有:

∂loss∂w2===∂sse(h_relu⋅w2)∂h_reluh_reluT⋅∂sse(y_pred)∂y_predh_reluT⋅grad_y_pred

3. grad_h_relu: 线性变换的导数

X右乘一个矩阵时的公式:
设有 f(Y):Rm×p→R 及线性映射 X↦Y=XC+D:Rm×n→Rm×p, 其中 C∈Rn×p,D∈Rm×p.则:

∇Xf(XC+D)=(∇Yf)CT

与求grad_w2的方法类似,根据以上规则有:

∂loss∂h_relu===∂sse(h_relu⋅w2)∂h_relu∂sse(y_pred)∂y_pred⋅w2Tgrad_y_pred⋅w2T

4. grad_h: 矩阵线性变换链式求导法则

已知:

h_relu=y_pred=loss=ReLU(h)h_relu⋅w2∑(y_pred−y)2=sse(y_pred)

做变量代换有:

∂loss∂h===∂sse(ReLU(h)⋅w2)∂h∂sse(y_pred)∂y_pred⋅∂(h_relu⋅w2)∂h_relu⋅∂ReLU(h)∂hgrad_h_relu⋅∂ReLU(h)∂h

关于ReLU函数的求导问题详见本文第二部分.

5. grad_w1: 矩阵线性变换链式求导法则

∂loss∂w1====∂sse(ReLU(x⋅w1)⋅w2)∂w1∂(x⋅w1)w1⋅grad_h_relu⋅∂ReLU(h)∂h∂(x⋅w1)∂w1⋅grad_hxT⋅grad_h

注意在第二个等号中, ∂(x⋅w1)w1 放到了最左边, 而不是最右边. 这样做是为了保证求导结果与 w1 同型.

关于链式法则的顺序记忆方法:

如果待求导的变量(矩阵)在左边,那么继续在最右边写关于它的导数; 如果待求导的变量在右边, 那么在最左边写关于它的导数. 这一方法适用于以上所有求导过程!

第二部分: 神经网络中非线性部分的矩阵求导

需要强调的是, 我们认为这三种情况下导数没有定义:
- 矩阵对向量
- 向量对矩阵
- 矩阵对矩阵求导
最自然的结果当然是把结果定义成三维乃至四维张量,但是这并不好算。也有一些绕弯的解决办法(例如把矩阵抻成一个向量等),但是这些方案都不完美(例如复合函数求导的链式法则无法用矩阵乘法简洁地表达等)。凡是遇到这种情况,都通过其他手段来绕过,后面会有具体的示例。

在第一部分的神经网络中出现了两个非线性函数: sse 和 ReLU. 这一部分着重讨论非线性函数的求导问题. 通常神经网络中有三类非线性操作: loss函数(如交叉熵,IoU), 激活函数(如sigmoid,tanh), 非线性层(如池化层,卷积层,batch norm层等等).

将它们都认为是函数, 分为两类, 一类为实值函数, 即输入为矩阵,输出为实数; 另一类为矩阵/向量值函数, 它们的输入和输出都为矩阵或向量. 一般情况下,我们认为loss函数是实值函数,激活函数是向量值函数, 非线性层视情况而定.

1. 实值函数的求导: sum, mean, max/min

实值函数对矩阵求导:

• 要点: 求导结果与自变量同型,且每个元素就是 f 对自变量的相应分量求导.
• 若函数 f:Rm×n→R, 则 ∂f/∂x 也是一个 m×n 维矩阵,且 (∂f/∂x)ij=∂f/∂xij。也可使用劈形算子将导数记作 ∇Xf.

示例:

• sum
  
  ∂sum(x)∂x=(1)(same shape as x)
• max
  
  ∂max(x)∂x={1,xij=max(x)0,otherwise

2. 向量值函数的求导: ReLU, Sigmoid

由于激活函数函数实际上是分别对每个分量进行计算,因此实际上激活函数对向量的求导结果的形式与其一维形式相同.

• ReLU:
  
  ∂ReLU(x)∂x={1,xij>00,otherwise
• sigmoid:
  
  ∂sigmoid(x)∂x=(sigmoid′(xij))

关于ReLU函数的求导问题:

Relu是一个非常优秀的激活函数，相比较于传统的Sigmoid函数，有三个优点：

1. 防止梯度弥散
   sigmoid的导数只有在0附近的时候有比较好的激活性，在正负饱和区的梯度都接近于0，所以这会造成梯度弥散，而relu函数在大于0的部分梯度为常数，所以不会产生梯度弥散现象。
2. 稀疏激活性
   
   • relu函数在负半区的导数为0 ，所以一旦神经元激活值进入负半区，那么梯度就会为0，也就是说这个神经元不会经历训练，即所谓的稀疏性。
   • 这也可以认为是它的一个缺点, 理论上不能用Gradient-Based方法。同时如果de-active了，容易无法再次active。
   • 如果你使用 ReLU，那么一定要小心设置learning rate， 而且要注意不要让你的网络出现很多“dead”神经元，如果这个问题不好解决，那么可以试试 Leaky ReLU、PReLU 或者Maxout.
3. 加快计算
   由于relu函数分段线性性质，导致其前传，后传，求导都是分段线性, 只需要一个if-else语句即可实现; 而sigmoid函数要进行浮点四则运算，并且由于两端饱和，在传播过程中容易丢弃信息.

虽然理论上不能用ReLU来进行梯度计算,但在实际中是大量运用的,因为ReLU虽然不是严格处处有导数的, 但是它确是处处可计算次梯度(subgradient)的, 这样已经足够进行反向传播计算了.

3. pooling 的求导

pooling 操作本来是不可导的, pooling操作使得feature map的尺寸发生变化, 假如做2×2的池化，假设那么第l+1层的feature map有16个梯度，那么第l层就会有64个梯度，这使得梯度无法对位的进行传播下去。其实解决这个问题的思想也很简单，就是把1个像素的梯度传递给4个像素，但是需要保证传递的梯度总和不变。根据这条原则，mean pooling和max pooling的反向传播也是不同的。

mean pooling

mean pooling 的前向传播就是把一个patch中的值取平均来做pooling, 那么反向传播的过程也就是把某个元素的梯度等分为n份分配给前一层，这样就保证池化前后的梯度（残差）之和保持不变，还是比较理解的，图示如下:

mean pooling比较容易让人理解错的地方就是会简单的认为直接把梯度复制N遍之后直接反向传播回去，但是这样会造成loss之和变为原来的N倍，网络是会产生梯度爆炸的。

max pooling

max pooling也要满足梯度之和不变的原则，max pooling的前向传播是把patch中最大的值传递给后一层，而其他像素的值直接被舍弃掉。那么反向传播也就是把梯度直接传给前一层某一个像素，而其他像素不接受梯度，也就是为0。所以max pooling操作和mean pooling操作不同点在于需要记录下池化操作时到底哪个像素的值是最大，也就是max id.

---

补充知识: 矩阵求导

1. 矩阵范数的导数

F范数：

∥A∥F=∑ij|Aij|2−−−−−−−−√=Tr(AAH)−−−−−−−−√

F范数平方的导数：

∂∂A∥A∥2F=∂∂ATr(AAH)=2A

2. 矩阵、向量、标量的导数

(黑体表示向量或矩阵，正常体为单个变量)

∂xTa∂x∂aTXb∂X∂aTXTb∂X∂X∂Xij=∂aTx∂x=a=abT=baT=Jij

3. 梯度与Hessian矩阵

fΔxf∂∂x(AB)∂2f∂x∂xT=xTAx+bTx=∂f∂x=(A+AT)x+b=∂A∂xB+A∂B∂x=A+AT

4. 链式法则

求解矩阵的函数的导数，例如 U=f(X)，要求 g(U) 关于 X 的导数：

∂g(U)∂X=∂g(f(X))∂X

那么链式法则如下：

∂g(U)∂X=∂g(U)∂Xij=∑k=1M∑l=1N∂g(U)∂ukl∂ukl∂Xij

写成矩阵形式如下：

∂g(U)Xij=Tr[(∂g(U)∂U)T∂U∂Xij]

参考资料:

1. 机器学习中的矩阵、向量求导
2. Matrix Cookbook
3. pooling的反向传播