(十)再探反向传播和神经网络优化

---

欢迎访问个人网络日志🌹🌹知行空间🌹🌹

---

1.背景介绍

反向传播算法是深度学习的基石，在2017年有一段时间总会跟着CS231N的课程讲义反复推导神经网络中
反向传播，时至今日再重新回顾一下神经网络中的优化方法。

2.神经网络的模型

对于线性分类器：$s=WX$，其中，x.shape为(N,D),N表示的是样本数目，D是X的特征维度，s.shape为(N,K),K表示类别的数目，s即是每个样本对应每个类别的评分。

全连接神经网络是有很多个线性分类器组成的，如带一个隐层的2层全连接网络$s=W_{2}max(0,W_{1}X)$,W1.shape为(D,H)，其中H表示的$W_{1}X$生成的中间数据的维度，W2.shape为(H,K)。max函数是引入的作用在向量中每个元素上的非线性单元函数，如果没有max，则$s=W_2{W}_{1}X$就等同于s=WX，其中$W=W_2{W}_1$，这样若没有非线性激活函数，则无论添加多少层都依然是一个线性分类器。网络图示如下：

单个神经元来分析，类比与生物学中神经细胞，有树突，轴突等组成，神经网络神经元由参数W和非线性激活函数组成，如下图：

上图中右侧神经元表示，输入的数据为3维$X=[x_0,x_1,x_2]$,对应的权值为:$W=[W_0,W_1,W_2]$,经线性分类器后又输入到非线性激活函数$f$中得到神经元输出结果。选择合适的激活函数如sigmoid，带1个隐含层的全连接神经网络就具备了拟合任何连续函数的能力。

3.神经网络中的参数更新初探

一个机器学习分类算法的核心由三大块组成，一个是评分函数score function将原始像素图像映射成类别评分；一个是损失函数loss function，用来评估预测结果与标签的一致程度，预测结果越准确loss越低，预测结果越不准，loss越高。另外一部分就是score function的参数更新即机器学习算法的优化问题,找到一组最优的参数$W$使得loss的值尽量低。如最简单的线性分类器$Y=WX+b$,损失函数使用SOFTMAX+CrossEntropy：

$$\begin{array}{c}{p}_k=\frac{e^{f_k}}{{\sum }_j{e}^{f_j}}\\ L_i=-log(p_{y_i})\end{array}$$

理解了神经网络后，要完成神经网络的训练还需要知道神经网络的参数是如何更新的。

先生成实验样本数据：

from matplotlib import pyplot as plt
import numpy as np
N = 100 # number of points per class
D = 2 # dimensionality
K = 3 # number of classes
X = np.zeros((N*K,D)) # data matrix (each row = single example)
y = np.zeros(N*K, dtype='uint8') # class labels
for j in range(K):
  ix = range(N*j,N*(j+1))
  r = np.linspace(0.0,1,N) # radius
  t = np.linspace(j*4,(j+1)*4,N) + np.random.randn(N)*0.2 # theta
  X[ix] = np.c_[r*np.sin(t), r*np.cos(t)]
  y[ix] = j
# lets visualize the data:
plt.scatter(X[:, 0], X[:, 1], c=y, marker='o', s=20)
plt.show()

3.1随机查找取最优

一个头脑风暴想法是随机选取W/b多次，然后取结果最好的那次权重:

last_loss = float('inf')
W = None
b = None
for i in range(10000):
    Wt = 0.015 * np.random.randn(D,K)
    bt = np.random.randn(1,K)

    scores = np.dot(X, Wt) + bt
    
    num_examples = X.shape[0]
    # get unnormalized probabilities
    exp_scores = np.exp(scores)
    # normalize them for each example
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    correct_logprobs = -np.log(probs[range(num_examples),y])

    # compute the loss: average cross-entropy loss and regularization
    data_loss = np.sum(correct_logprobs)/num_examples
    reg_loss = 0.5*reg*np.sum(Wt*Wt)
    loss = data_loss + reg_loss
    if loss < last_loss:
        last_loss = loss
        W = Wt.copy()
        b = bt.copy()
    print(f"step: {i} loss: {loss} min loss: {last_loss}")

# evaluate training set accuracy
scores = np.dot(X, W) + b
predicted_class = np.argmax(scores, axis=1)
print('training accuracy: %.2f' % (np.mean(predicted_class == y)))
## training accuracy: 0.46

随机查找最优结果参数可以得到的分类准确率是46%，可以知道等样本3分类的随机概率是33%，因此随机查找最优得到了比盲猜要好的结果，虽然其准确率依然不高。其分类结果如下图：

3.2局部随机查找参数更新

类比于一个人在山顶迷路了要下山，其每往前走一下其是能够知道这一步是在往下走还是在往上升，只有当迈出的一步可以使得位置下降时才前进一步到达一个新位置，再重新开始寻找下一步。回到参数更新，也可以通过这种Iterative Refinement一步步找寻最优的参数：

reg = 0.001
step_size = 1.

last_loss = float('inf')
W = 0.015 * np.random.randn(D,K)
b = 0.015 * np.random.randn(1,K)
for i in range(10000):
    Wt = W + np.random.randn(D,K) * step_size
    bt = b + np.random.randn(1,K) * step_size

    scores = np.dot(X, Wt) + bt
    
    num_examples = X.shape[0]
    # get unnormalized probabilities
    exp_scores = np.exp(scores)
    # normalize them for each example
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    correct_logprobs = -np.log(probs[range(num_examples),y])

    # compute the loss: average cross-entropy loss and regularization
    data_loss = np.sum(correct_logprobs)/num_examples
    reg_loss = 0.5*reg*np.sum(Wt*Wt)
    loss = data_loss + reg_loss
    if loss < last_loss:
        last_loss = loss
        W = Wt.copy()
        b = bt.copy()
    print(f"step: {i} loss: {loss} min loss: {last_loss}")
# evaluate training set accuracy
scores = np.dot(X, W) + b
predicted_class = np.argmax(scores, axis=1)
print('training accuracy: %.2f' % (np.mean(predicted_class == y)))
# training accuracy: 0.56

可以看到使用逐步迭代更新参数的方法，分类的效果得到了一定的提升。

3.3 沿着梯度反方向更新

梯度方向是函数值增长最快的方向，那么梯度的反方向就是函数值下降最快的方向。损失函数对参数求导数即可找到损失函数值下降最快的方向，很自然的想到沿着梯度反方向更新参数。

梯度的求解有两种方法：

• 数值解：$\frac{[f(x+h)-f(x-h)]}{2h}$
• 解析解：即求得函数的导数
• 数值解求解方便，但只是近似结果，有可能有比较大的累积误差，且计算量比较大，而解析解计算较为复杂，推导容易出错。因此通常比较两种梯度计算的结果，以验证梯度计算的是否准确。

reg = 0.0001
step_size = 2

W = np.random.randn(D,K) * 0.001
b = np.random.randn(1,K) * 0.001
for i in range(1000):
    scores = np.dot(X, W) + b
    num_examples = X.shape[0]
    # get unnormalized probabilities
    exp_scores = np.exp(scores)
    # normalize them for each example
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    correct_logprobs = -np.log(probs[range(num_examples),y])

    # compute the loss: average cross-entropy loss and regularization
    data_loss = np.sum(correct_logprobs)/num_examples
    reg_loss = 0.5*reg*np.sum(Wt*Wt)
    loss = data_loss + reg_loss
    print(f"{i} loss: {loss}")
    dscores = probs
    dscores[range(num_examples),y] -= 1
    dscores /= num_examples

    dW = np.dot(X.T, dscores)
    db = np.sum(dscores, axis=0, keepdims=True)
    dW += reg*W # don't forget the regularization gradient
    
    W += -dW * step_size
    b += -db * step_size
    
# evaluate training set accuracy
scores = np.dot(X, W) + b
predicted_class = np.argmax(scores, axis=1)
print('training accuracy: %.2f' % (np.mean(predicted_class == y)))
# training accuracy: 0.57

受限于线性分类器的能力，使用梯度下降方法对准确率的提升已经不大，但模型学习速度要快的多,经过几十步迭代后，loss就逐渐稳定了，可见梯度反方向的确应该是参数更新的方向，程序中step_size是每一步更新权值大小，也被称为学习率，这也是机器学习中比较关键的参数，后续章节介绍。

step 0 loss: 1.0987081944507475
step 1 loss: 1.035012790265596
step 2 loss: 0.9847131322096168
step 3 loss: 0.9446124386650137
step 4 loss: 0.9122820461040575
step 5 loss: 0.8859083863629116
step 6 loss: 0.8641456746130115
step 7 loss: 0.8459924630117203
step 8 loss: 0.8306977637849042
step 9 loss: 0.8176927238746483
step 10 loss: 0.8065417050365444
step 11 loss: 0.7969073565356761
step 12 loss: 0.7885255785711087
step 13 loss: 0.7811874419263808
step 14 loss: 0.7747260116459382
step 15 loss: 0.7690066494799522
step 16 loss: 0.7639198043202386
step 17 loss: 0.7593755983752785
step 18 loss: 0.7552997218726981
step 19 loss: 0.7516302904789174
step 20 loss: 0.74831541777118
step 21 loss: 0.7453113237504392
step 22 loss: 0.7425808488143958
step 23 loss: 0.7400922770753964
step 24 loss: 0.7378183976522287
step 25 loss: 0.7357357504859927
step 26 loss: 0.7338240163214523
step 27 loss: 0.7320655201385343
step 28 loss: 0.730444824479833
step 29 loss: 0.728948394479751
step 30 loss: 0.7275643204427976
step 31 loss: 0.7262820868888554
step 32 loss: 0.7250923793318111
step 33 loss: 0.7239869218665891
step 34 loss: 0.722958340041463
step 35 loss: 0.7220000445858086
step 36 loss: 0.7211061324211997
step 37 loss: 0.7202713020605488

这里使用的梯度下降算法直接拿所有训练样本参与训练，但对于大规模训练数据，一次性加载计算是不可能实现的，像ILSVRC有120万张图像。此时可以使用mini batch gradient descent算法，每次选择打散数据中batch_size数目的数据参与训练用于更新参数，mini batch gradient descent只所以能够工作，是因为训练数据之间是相关的，考虑极端情况，ILSVRC的120万张图像有1000个类，每个类的1200张图像都是相同的，则此时选一部分数据更新参数与利用所有数据同时参与训练，效果是等效的，实际中虽然相同类别的不同图像不可能完全相同，但其之间是有一定共同特征的，因此使用mini batch gradient descent算法训练模型依然能够收敛。

4.链式法则与反向传播

通过#3中的介绍，以线性分类器$Y=WX+b$为例，寻找最优参数$W$以得到最优分类结果的方法依赖损失对参数$W$的导数，上述例子中使用的是线性分类器，评分映射函数还比较简单梯度也比较容易求解，但对于神经网络来说，其往往包含很多层，且层与层之间还有非线性激活函数，直接对其求导几乎是不可能的，这时可以使用链式法则和反向传播算法。

以函数：$f=(x+y)z$为例，虽然这个函数还是很简单的，能够直接对$x/y/z$求偏导数，但为了说明复合函数求导和链式法则，将其写成$f=qz,q=x+y$的形式，这时：

$$\begin{array}{c}f=qz,q=x+y\\ \frac{\partial f}{\partial z}=q\\ \frac{\partial f}{\partial q}=z\\ \frac{\partial f}{\partial x}=\frac{\partial f}{\partial q}\frac{\partial q}{\partial x}=z\\ \frac{\partial f}{\partial y}=\frac{\partial f}{\partial q}\frac{\partial q}{\partial y}=z\end{array}$$

梯度求解过程：

x=-2; y=5; z=-4

q = x + y # q=3
f = q * z # f=-12

## 反向传播求梯度
dfdz = q # dfdz=3
dfdq = z # dfdq=-4
dqdx = 1
dqdy = 1

dfdx = dfdq * dqdx # dfdx = -4, 此处就是链式法则的应用
dfdy = dfdq * dqdy # dfdx = -4

链式法则的形象化图示如下：

5.梯度方向的参数更新策略

• 1.普通梯度下降策略

这里先实现如图1所示的带一层隐含层的全连接神经网络，X是二维特征的数据$X=[x_1,x_2]$,隐含层有4个神经元，则$W_1.shape=(2,4),b_1.shape=(1,4)$，共有3个类别，则$W_2.shape=(4,3),b_2.shape=(1,3)$:
$$\begin{array}{c}{Y}_1=W_{1}X+b_1\\ Y1=max(0,Y_1)\\ score=W_{2}Y1+b_2\end{array}$$

写成矩阵形式：

$$\begin{array}{c}[y_{11},y_{12},y_{13},y_{14}]=[x_{11},x_{12}]\left[\begin{array}{cccc}w{1}_{11}& w{1}_{12}& w{1}_{13}& w{1}_{14}\\ w{1}_{21}& w{1}_{22}& w{1}_{23}& w{1}_{24}\end{array}\right]+[b{1}_1,b{1}_2,b{1}_3,b{1}_4]\\ [y_{11},y_{12},y_{13},y_{14}]=max\_elementwise([y_{11},y_{12},y_{13},y_{14}],0)\\ [f_{11},f_{12},f_{13}]=[y_{11},y_{12},y_{13},y_{14}]\left[\begin{array}{ccc}w{2}_{11}& w{2}_{12}& w{2}_{13}\\ w{2}_{21}& w{2}_{22}& w{2}_{23}\\ w{2}_{31}& w{2}_{32}& w{2}_{33}\\ w{2}_{41}& w{2}_{42}& w{2}_{43}\end{array}\right]+[b{2}_1,b{2}_2,b{2}_3]\end{array}$$

损失函数使用CrossEntropy,i表示第i个样本，同上述举证形式中的每行，k表示当前样本所属类别:

$$\begin{array}{c}{p}_k=\frac{e^{f_k}}{{\sum }_j{e}^{f_j}}\\ L_i=-log(p_{y_i})\end{array}$$
补充一下，多分类交叉熵损失函数的形式一般写为：

$$L_i=-\sum _{i=1}^K{y}_{true}log(y_{pre})$$
其中，$y_{pre}$是经SOFTMAX后输出的在该类别上的评分，$y_{true}$是经过One-Hot编码后的向量，如3个类时[0,0,1]。可以看到在不是当前类别的标签上都是0,因此$L_i$可以写成$L_i=-log({y_k}_{pred})$的形式。

先对上述全连接神经网络和损失函数应用链式法则求导：

$$\begin{array}{c}\frac{\partial L_i}{\partial p_{y_i}}=-\frac{1}{p_{y_i}}\\ y_i=k,p_{y_i}=p_k,p_k=\frac{e^{f_k}}{c+e^{f_k}}其中c={\sum }_j{e}^{f_j},j\ne k\\ \frac{\partial p_k}{\partial f_k}=\frac{c{e}^{f_k}}{(c+e^{f_k}{)}^2}\\ \frac{\partial L_i}{\partial f_k}=\frac{\partial L_i}{\partial p_k}\frac{\partial p_k}{\partial f_k}=-\frac{c+e^{f_k}}{e^{f_k}}\frac{c{e}^{f_k}}{(c+e^{f_k}{)}^2}=-\frac{c}{c+e^{f_k}}=\frac{e^{f_k}-(c+e^{f_k})}{c+e^{f_k}}=(p_k-1),y_i=k\\ \frac{\partial L_i}{\partial f_{y_i}}=\frac{\partial L_i}{\partial p_k}\frac{\partial p_k}{\partial f_{y_i}}=p_{y_i},y_i\ne k\\ \frac{\partial f}{\partial w}=Y1\\ \frac{\partial f}{\partial b}=1\\ \frac{\partial Y1}{\partial Y_1}=1,Y_1>0\end{array}$$

以上就根据链式法则求得了两层全连接神经网络的梯度，对于N个样本输入时，使用矩阵表示可同样推导。根据以上推导实现的两层全连接神经网络为：

H = 100
W = 0.015 * np.random.randn(D,H)
b = np.zeros((1,H))

W1 = 0.015 * np.random.randn(H,K)
b1 = np.zeros((1,K))

reg = 0.001
step_size = 1.5

def sigmoid(x):
    return 1/(1+np.exp(-x))

mu = 0.9
v1 = 0
v2 = 0
v3 = 0
v4 = 0
for i in range(5000):
    Y = np.dot(X, W) + b
    sY = np.maximum(Y, 0)
    scores = np.dot(sY, W1) + b1
    
    num_examples = X.shape[0]
    # get unnormalized probabilities
    exp_scores = np.exp(scores)
    # normalize them for each example
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    correct_logprobs = -np.log(probs[range(num_examples),y])

    # compute the loss: average cross-entropy loss and regularization
    data_loss = np.sum(correct_logprobs)/num_examples
    reg_loss = 0.5*reg*np.sum(W*W) + 0.5*reg*np.sum(W1*W1)
    loss = data_loss + reg_loss
    print(f"step: {i} loss: {loss}")
    dscores = probs
    dscores[range(num_examples),y] -= 1
    dscores /= num_examples

    print(dscores.shape)
    dW1 = np.dot(sY.T, dscores)
    db1 = np.sum(dscores, axis=0, keepdims=True)
    dW1 += reg*W1

    dsY = np.dot(dscores, W1.T)
    dsY = dsY * sigmoid(Y)*(1-sigmoid(Y))
    dsY[sY<=0] = 0

    dW = np.dot(X.T, dsY)
    db = np.sum(dsY, axis=0, keepdims=True)
    dW += reg*W # don't forget the regularization gradient
    
    # perform a parameter update
    # Vanilla update
    # step: 5403 loss: 0.30999144570577264
    W += - step_size * dW
    b += - step_size * db
    W1 += - step_size * dW1
    b1 += - step_size * db1

使用普通梯度的更新，模型的分类准确率轻松达到了97%。

• 2.动量更新策略

将损失函数类比为丘陵，则势能U=mgh,$U\propto h$,使用随机数字初始化参数，等同于将一个物体在某个位置的初速度设置为0，如此，可将优化过程类比为将物体从山顶往下滚。由势能的定义，物体所受的力可通过势能$U$对位置$x$求导得到，即$F=▽U=\frac{dU}{dh}$而$F=ma$,因此$U$对$h$的梯度与物体的加速度成正比，对于优化问题，$loss$即是$U$,$h$即是参数$W$,设初始速度为$v=0$,则参数更新不是直接作用在位置$x$上，而是作用在速度$v$上，这是与普通梯度更新的区别，其中超参数mu表示为物体的动量，通常取0.9。

# Momentum update
v = mu * v - learning_rate * dx # integrate velocity
x += v # integrate position

## 作用在前述全连接神经网络实现上
v1 = mu * v1 - step_size * dW
W += v1
v2 = mu * v2 - step_size * db
b += v2
v3 = mu * v3 - step_size * dW1
W1 += v3
v4 = mu * v4 - step_size * db1
b1 += v4

• 3.Nesterov动量更新策略

Nesterov Momentum是由俄国数学家Yurii Nesterov提出的凸优化方法，Nesterov Momentum的核心是对于位置x处的参数向量,只看动量部分，其动量对参数向量的贡献部分为mv，当更新梯度时，可以将x+mu*v当做远眺，求x+mu*v处的梯度用以更新参数x,而非使用滞后的位置x处的梯度。

x_ahead = x + mu * v
# evaluate dx_ahead (the gradient at x_ahead instead of at x)
v = mu * v - learning_rate * dx_ahead
x += v

## 通常利用`x_ahead = x + mu * v`将上述表示写为：

v_prev = v # back this up
v = mu * v - learning_rate * dx # velocity update stays the same
x += -mu * v_prev + (1 + mu) * v # position update changes form


## 应用到前述的两层全连接神经网络中
mu = 0.9
v1 = 0
v2 = 0
v3 = 0
v4 = 0

v1p = v1
v1 = mu * v1 - step_size * dW
W += -mu * v1p + (1+mu)*v1

v2p = v2
v2 = mu * v2 - step_size * db
b += -mu * v2p + (1+mu)*v2

v3p = v3
v3 = mu * v3 - step_size * dW1
W1 += -mu * v3p + (1+mu)*v3

v4p = v4
v4 = mu * v4 - step_size * db1
b1 += -mu * v4p + (1+mu)*v4

• 4.三种方法对比

将以上三种参数更新方法的优化过程中loss下降情况绘制如下图，可以看到Momentum比Vanilla收敛的要快，Nesterov Momentum收敛过程比较振荡，但收敛的loss更低，收敛的效果的更好。

6.学习率退火

6.1学习率衰减策略基础

参考如上参数更新过程，学习率过大则每一步更新的参数量过大，会导致模型振荡且收敛于loss较大的状态，学习效果较差；而学习率过小，会导致模型参数更新速度过慢，花费过多训练时间。自然而然的会想到，随着训练step的变大逐步减小学习率，这就是学习率退火的原理。通常有三种学习率的衰减策略：

• 1.按步衰减：Step decay，如每5epoches学习率衰减一半，或者20 epochs衰减0.1，一个启发式的办法是使用固定学习率，当验证集上的准确率不再上升时则减小学习率为原来的0.5倍。

• 2.指数衰减：$\alpha ={\alpha }_0{e}^{-kt}$,其中${\alpha }_0、k$是超参数，$t$是迭代$step$

• 3.$\frac{1}{t}$衰减，$\alpha ={\alpha }_0/(1+kt)$,其中${\alpha }_0、k$是超参数，$t$是迭代$step$

6.2 二阶优化方法

首先看泰勒定理：

• 牛顿法

对于函数F(x),考虑step=k，假设在$x_k$处寻找增量$△x_k$,最直观的方式是将函数在$x_k$处先进行泰勒展开：
$$F(x_k+△x_k)\approx F(x_k)+J(x_k{)}^T△x_k+\frac{1}{2}△{x_k}^{T}H(x_k)△x_k$$

$J(x)$是$F(x)$关于$x$的一阶梯度也称为雅可比矩阵(Jacobin Matrix)，$H(x)$是$F(x)$关于$x$的二阶梯度也称为海塞矩阵(Hessian Matrix)。前面讨论的优化方法使用的是一阶导数$△x^{\ast }=-\alpha J(x)$,其中$\alpha$即参数每步更新的步长，这种方法简单称为最速下降法。除此之外，还可以选择保留二阶梯度信息，此时参数更新方程可写为：

$$△x^{\ast }=\underset{△x}{\mathrm{argmin}}(F(x)+J(x{)}^T△x+\frac{1}{2}△x^{T}H(x)△x)$$

上述等式右侧是$△x$的二次多项式，求其对$△x$的导数并令其等于$0$即可求得$△x^{\ast }$
$$J(x)+H(x)△x=0\Rightarrow △x^{\ast }=-H^{-1}(x)J(x)$$

通过求解以上方程就得到了参数的更新量$△x^{\ast }$,这种方法被称为牛顿法(Newton’s method)。，可以看到这种方法中使用了二阶梯度来控制每次参数更新的步长，省去了超参数学习率$\alpha$，可以有更快的收敛速度。但是，对于神经网络动辄上百万的参数，其Hessian Matrix占据的内存太大，目前还不能够直接使用，如对于100M参数的神经网络，其Hessian Matrix维度为$[1000000,1000000]$,将需要$3725GB$RAM。实际中运用的多是牛顿法的近似，如高斯牛顿法/列文伯格马夸尔特方法。

• 高斯牛顿法

对于非线性最小二乘优化问题$F(x)=\frac{1}{2}||f(x)|{|}_2^2,\underset{x}{\min }F(x)$，除了使用上述介绍的牛顿法外，还可以使用高斯牛顿法(Gauss Newton)，避免计算Hessian Matrix。其原理是将f(x)在$x$处进行一阶泰勒展开：

$$f(x+△x)\approx f(x)+J(x{)}^T△x$$
注意这里是对$f(x)$而非目标函数$F(x)$进行泰勒展开。这里$f(x)$与$F(x)$的区别可通过一个例子来看就会更清楚，

https://zhuanlan.zhihu.com/p/42383070

根据对$f(x)$的一阶泰勒展开，优化的目标写成最小化$||f^2(x+△x)||$,目标是求$△x$能够使其最小，则：

$$△x^{\ast }=\underset{△x}{\mathrm{argmin}}\frac{1}{2}||f(x)+J(x{)}^T△x|{|}^2$$
可将上述方程写成：
$$\begin{array}{rl}\frac{1}{2}||f(x)+J(x{)}^T△x|{|}^2& =\frac{1}{2}(f(x)+J(x{)}^T△x{)}^T(f(x)+J(x{)}^T△x)\\ & =\frac{1}{2}(||f(x)|{|}_2^2+2f(x)J(x{)}^T△x+△x^{T}J(x)J^T(x)△x)\end{array}$$
根据其求极值条件，可求得$△x$为：
$$\begin{array}{c}J(x)f(x)+J(x)J^T(x)△x=0\\ \underset{H(x)}{\underbrace{J(x)J^T(x)}}△x=\underset{g(x)}{\underbrace{-J(x)f(x)}}\end{array}$$

上述方程被称为增量方程、高斯牛顿方程(Gauss Newton Equation)、正规方程(Normal Equation)。高斯牛顿方法避免了计算二阶梯度，使用$J{J}^T$近似$H$,但$J{J}^T$只有半正定，因此$J{J}^T$有可能出现奇异或病态，导致计算出来的$△x$过大，这时可以考虑使用列文伯格-马尔夸特方法，其在每一步参数更新中添加了一个参数指标$\rho =\frac{f(x+△x)-f(x)}{J^T(x)△x}$来刻画每一步参数更新的后函数$f(x)近似的好坏，并通过添加约束 $||D△x|{|}^2\le \mu$来限制每一步参数更新量的大小，以应对$H$病态的情况。

6.3自适应学习率方法

虽然前面介绍的二阶梯度方法很好，但因神经网络中参数量、训练数据量巨大，很难应用。像L-BFGS算法需要在整个数据集上同时计算二阶梯度在mini-batch上无法工作。因此，神经网络优化中使用最多的是Nestov Momentum方法加自适应学习率。

• Adagrad：该方法是John Duchi于2011年提出的方法。

  # Assume the gradient dx and parameter vector x
  cache += dx**2
  x += - learning_rate * dx / (np.sqrt(cache) + eps)
  

  其中，参数cache和dx有相同的维度，故其可以根据参数更新的大小自适应衰减学习率。缺点是研究表明神经网络过程中使用单调下降的学习率会导致参数下降过快，模型过早的停止收敛。

• RMSprop: 该方法是机器学习专家Geoff Hinton在其Cousera课程讲义中提出，目前并没有发表。

  # decay_rate是超参数，常取[0.9, 0.99, 0.999]
  cache = decay_rate * cache + (1 - decay_rate) * dx**2
  x += - learning_rate * dx / (np.sqrt(cache) + eps)    
  

• Adam：该方法是2014年12月Diederik P. Kingma提出的方法,该方法和带动量的RMSprop有些相似，

  m = beta1*m + (1-beta1)*dx
  v = beta2*v + (1-beta2)*(dx**2)
  x += - learning_rate * m / (np.sqrt(v) + eps)
  

  可以看到除了使用$m$替代了$dx$外，Adam和RMSprop是一样的，其中超参数的常用值为eps = 1e-8,beta1 = 0.9，beta2 = 0.999，一般来说Adam是目前最好的方法，比RMSprop好一些。

  #     ## raw
  #     cw, cb, cw1, cb1 = 1, 1, 1, 1
  #     ## Adagrad
  #     cw += dW**2
  #     cb += db**2
  #     cw1 += dW1**2
  #     cb1 += db1**2
      
  #     ## rmsprop
  #     cw = decay*cw + (1-decay)*dW**2
  #     cb = decay*cb + (1-decay)*db**2
  #     cw1 = decay*cw1 + (1-decay)*dW1**2
  #     cb1 = decay*cb1 + (1-decay)*db1**2
  
  ## Adam
  t = i + 1
  cw = decay*cw + (1-decay)*dW**2
  cwt = cw / (1-decay**t)
  cb = decay*cb + (1-decay)*db**2
  cbt = cb / (1-decay**t)
  cw1 = decay*cw1 + (1-decay)*dW1**2
  cw1t = cw1 / (1-decay**t)
  cb1 = decay*cb1 + (1-decay)*db1**2
  cb1t = cb1 / (1-decay**t)
  mw = dm*cw + (1-dm)*dW
  mwt = mw / (1-dm**t)
  mb = dm*cb + (1-dm)*db
  mbt = mb / (1-dm**t)
  mw1 = dm*cw1 + (1-dm)*dW1
  mw1t = mw1 / (1-dm**t)
  mb1 = dm*cb1 + (1-dm)*db1
  mb1t = mb1 / (1-dm**t)
  
  W += - step_size * mwt / (np.sqrt(cwt) + eps)
  b += - step_size * mbt / (np.sqrt(cbt) + eps)
  W1 += - step_size * mw1t / (np.sqrt(cw1t) + eps)
  b1 += - step_size * mb1t / (np.sqrt(cb1t) + eps)
  

如上实验中，RMSprop和Adam优化效果较好，其中RMSprop比Adam的优化效果还要好些。

---

欢迎访问个人网络日志🌹🌹知行空间🌹🌹

---

参考资料

• 1.https://cs231n.github.io/neural-networks-3/
• 2.视觉SLAM十四讲第六讲
• 3.https://zhuanlan.zhihu.com/p/42383070