动量法

🏷sec_momentum

在 :numref:sec_sgd一节中，我们详述了如何执行随机梯度下降，即在只有嘈杂的梯度可用的情况下执行优化时会发生什么。
对于嘈杂的梯度，我们在选择学习率需要格外谨慎。
如果衰减速度太快，收敛就会停滞。
相反，如果太宽松，我们可能无法收敛到最优解。

基础

本节将探讨更有效的优化算法，尤其是针对实验中常见的某些类型的优化问题。

泄漏平均值

上一节中我们讨论了小批量随机梯度下降作为加速计算的手段。
它也有很好的副作用，即平均梯度减小了方差。
小批量随机梯度下降可以通过以下方式计算：

为了保持记法简单，在这里我们使用\mathbf{h}_{i, t-1} = \partial_{\mathbf{w}} f(\mathbf{x}_i, \mathbf{w}_{t-1})作为样本i的随机梯度下降，使用时间t-1时更新的权重t-1。
如果我们能够从方差减少的影响中受益，甚至超过小批量上的梯度平均值，那很不错。
完成这项任务的一种选择是用泄漏平均值（leaky average）取代梯度计算：

其中\beta \in (0, 1)。
这有效地将瞬时梯度替换为多个“过去”梯度的平均值。
\mathbf{v}被称为动量（momentum），
它累加了过去的梯度。
为了更详细地解释，让我们递归地将\mathbf{v}_t扩展到

其中，较大的\beta相当于长期平均值，而较小的\beta相对于梯度法只是略有修正。
新的梯度替换不再指向特定实例下降最陡的方向，而是指向过去梯度的加权平均值的方向。
这使我们能够实现对单批量计算平均值的大部分好处，而不产生实际计算其梯度的代价。

上述推理构成了"加速"梯度方法的基础，例如具有动量的梯度。
在优化问题条件不佳的情况下（例如，有些方向的进展比其他方向慢得多，类似狭窄的峡谷），"加速"梯度还额外享受更有效的好处。
此外，它们允许我们对随后的梯度计算平均值，以获得更稳定的下降方向。
诚然，即使是对于无噪声凸问题，加速度这方面也是动量如此起效的关键原因之一。

正如人们所期望的，由于其功效，动量是深度学习及其后优化中一个深入研究的主题。
例如，请参阅[文章](https://distill.pub/2017/momentum/）（作者是 :cite:Goh.2017)，观看深入分析和互动动画。
动量是由 :cite:Polyak.1964提出的。
:cite:Nesterov.2018在凸优化的背景下进行了详细的理论讨论。
长期以来，深度学习的动量一直被认为是有益的。
有关实例的详细信息，请参阅 :cite:Sutskever.Martens.Dahl.ea.2013的讨论。

条件不佳的问题

为了更好地了解动量法的几何属性，我们复习一下梯度下降，尽管它的目标函数明显不那么令人愉快。
回想我们在 :numref:sec_gd中使用了f(\mathbf{x}) = x_1^2 + 2 x_2^2，即中度扭曲的椭球目标。
我们通过向x_1方向伸展它来进一步扭曲这个函数

与之前一样，f在(0, 0)有最小值，
该函数在x_1的方向上非常平坦。
让我们看看在这个新函数上执行梯度下降时会发生什么。

%matplotlib inline
from d2l import mxnet as d2l
from mxnet import np, npx
npx.set_np()

eta = 0.4
def f_2d(x1, x2):
    return 0.1 * x1 ** 2 + 2 * x2 ** 2
def gd_2d(x1, x2, s1, s2):
    return (x1 - eta * 0.2 * x1, x2 - eta * 4 * x2, 0, 0)

d2l.show_trace_2d(f_2d, d2l.train_2d(gd_2d))
​

#@tab pytorch
%matplotlib inline
from d2l import torch as d2l
import torch

eta = 0.4
def f_2d(x1, x2):
    return 0.1 * x1 ** 2 + 2 * x2 ** 2
def gd_2d(x1, x2, s1, s2):
    return (x1 - eta * 0.2 * x1, x2 - eta * 4 * x2, 0, 0)

d2l.show_trace_2d(f_2d, d2l.train_2d(gd_2d))
​

#@tab tensorflow
%matplotlib inline
from d2l import tensorflow as d2l
import tensorflow as tf

eta = 0.4
def f_2d(x1, x2):
    return 0.1 * x1 ** 2 + 2 * x2 ** 2
def gd_2d(x1, x2, s1, s2):
    return (x1 - eta * 0.2 * x1, x2 - eta * 4 * x2, 0, 0)

d2l.show_trace_2d(f_2d, d2l.train_2d(gd_2d))
​

#@tab paddle
%matplotlib inline
from d2l import paddle as d2l
import warnings
warnings.filterwarnings("ignore")
import paddle

eta = 0.4
def f_2d(x1, x2):
    return 0.1 * x1 ** 2 + 2 * x2 ** 2
def gd_2d(x1, x2, s1, s2):
    return (x1 - eta * 0.2 * x1, x2 - eta * 4 * x2, 0, 0)

d2l.show_trace_2d(f_2d, d2l.train_2d(gd_2d))
​

从构造来看，x_2方向的梯度比水平x_1方向的梯度大得多，变化也快得多。
因此，我们陷入两难：如果选择较小的学习率，我们会确保解不会在x_2方向发散，但要承受在x_1方向的缓慢收敛。相反，如果学习率较高，我们在x_1方向上进展很快，但在x_2方向将会发散。
下面的例子说明了即使学习率从0.4略微提高到0.6，也会发生变化。
x_1方向上的收敛有所改善，但整体来看解的质量更差了。

#@tab all
eta = 0.6
d2l.show_trace_2d(f_2d, d2l.train_2d(gd_2d))
​

动量法

动量法（momentum）使我们能够解决上面描述的梯度下降问题。
观察上面的优化轨迹，我们可能会直觉到计算过去的平均梯度效果会很好。
毕竟，在x_1方向上，这将聚合非常对齐的梯度，从而增加我们在每一步中覆盖的距离。
相反，在梯度振荡的x_2方向，由于相互抵消了对方的振荡，聚合梯度将减小步长大小。
使用\mathbf{v}_t而不是梯度\mathbf{g}_t可以生成以下更新等式：

请注意，对于\beta = 0，我们恢复常规的梯度下降。
在深入研究它的数学属性之前，让我们快速看一下算法在实验中的表现如何。

#@tab all
def momentum_2d(x1, x2, v1, v2):
    v1 = beta * v1 + 0.2 * x1
    v2 = beta * v2 + 4 * x2
    return x1 - eta * v1, x2 - eta * v2, v1, v2

eta, beta = 0.6, 0.5
d2l.show_trace_2d(f_2d, d2l.train_2d(momentum_2d))
​

正如所见，尽管学习率与我们以前使用的相同，动量法仍然很好地收敛了。
让我们看看当降低动量参数时会发生什么。
将其减半至\beta = 0.25会导致一条几乎没有收敛的轨迹。
尽管如此，它比没有动量时解将会发散要好得多。

#@tab all
eta, beta = 0.6, 0.25
d2l.show_trace_2d(f_2d, d2l.train_2d(momentum_2d))
​

请注意，我们可以将动量法与随机梯度下降，特别是小批量随机梯度下降结合起来。
唯一的变化是，在这种情况下，我们将梯度\mathbf{g}_{t, t-1}替换为\mathbf{g}_t。
为了方便起见，我们在时间t=0初始化\mathbf{v}_0 = 0。

有效样本权重

回想一下\mathbf{v}_t = \sum_{\tau = 0}^{t-1} \beta^{\tau} \mathbf{g}_{t-\tau, t-\tau-1}。
极限条件下，\sum_{\tau=0}^\infty \beta^\tau = \frac{1}{1-\beta}。
换句话说，不同于在梯度下降或者随机梯度下降中取步长\eta，我们选取步长\frac{\eta}{1-\beta}，同时处理潜在表现可能会更好的下降方向。
这是集两种好处于一身的做法。
为了说明\beta的不同选择的权重效果如何，请参考下面的图表。

#@tab all
d2l.set_figsize()
betas = [0.95, 0.9, 0.6, 0]
for beta in betas:
    x = d2l.numpy(d2l.arange(40))
    d2l.plt.plot(x, beta ** x, label=f'beta = {beta:.2f}')
d2l.plt.xlabel('time')
d2l.plt.legend();
​

实际实验

让我们来看看动量法在实验中是如何运作的。
为此，我们需要一个更加可扩展的实现。

从零开始实现

相比于小批量随机梯度下降，动量方法需要维护一组辅助变量，即速度。
它与梯度以及优化问题的变量具有相同的形状。
在下面的实现中，我们称这些变量为states。

#@tab mxnet, pytorch
def init_momentum_states(feature_dim):
    v_w = d2l.zeros((feature_dim, 1))
    v_b = d2l.zeros(1)
    return (v_w, v_b)
​

#@tab tensorflow
def init_momentum_states(features_dim):
    v_w = tf.Variable(d2l.zeros((features_dim, 1)))
    v_b = tf.Variable(d2l.zeros(1))
    return (v_w, v_b)
​

#@tab paddle
def init_momentum_states(feature_dim):
    v_w = d2l.zeros((feature_dim, 1))
    v_b = d2l.zeros([1])
    return (v_w, v_b)
​

def sgd_momentum(params, states, hyperparams):
    for p, v in zip(params, states):
        v[:] = hyperparams['momentum'] * v + p.grad
        p[:] -= hyperparams['lr'] * v
​

#@tab pytorch
def sgd_momentum(params, states, hyperparams):
    for p, v in zip(params, states):
        with torch.no_grad():
            v[:] = hyperparams['momentum'] * v + p.grad
            p[:] -= hyperparams['lr'] * v
        p.grad.data.zero_()
​

#@tab tensorflow
def sgd_momentum(params, grads, states, hyperparams):
    for p, v, g in zip(params, states, grads):
            v[:].assign(hyperparams['momentum'] * v + g)
            p[:].assign(p - hyperparams['lr'] * v)
​

#@tab paddle
def sgd_momentum(params, states, hyperparams):
    a = []
    for p, v in zip(params, states):
        with paddle.no_grad():
            v[:] = hyperparams['momentum'] * v + p.grad
            p[:] -= hyperparams['lr'] * v
        p.grad.zero_()
        a.append(p)
    return a
​

让我们看看它在实验中是如何运作的。

#@tab all
def train_momentum(lr, momentum, num_epochs=2):
    d2l.train_ch11(sgd_momentum, init_momentum_states(feature_dim),
                   {'lr': lr, 'momentum': momentum}, data_iter,
                   feature_dim, num_epochs)

data_iter, feature_dim = d2l.get_data_ch11(batch_size=10)
train_momentum(0.02, 0.5)
​

当我们将动量超参数momentum增加到0.9时，它相当于有效样本数量增加到\frac{1}{1 - 0.9} = 10。
我们将学习率略微降至0.01，以确保可控。

#@tab all
train_momentum(0.01, 0.9)
​

降低学习率进一步解决了任何非平滑优化问题的困难，将其设置为0.005会产生良好的收敛性能。

#@tab all
train_momentum(0.005, 0.9)
​

简洁实现

由于深度学习框架中的优化求解器早已构建了动量法，设置匹配参数会产生非常类似的轨迹。

d2l.train_concise_ch11('sgd', {'learning_rate': 0.005, 'momentum': 0.9},
                       data_iter)
​

#@tab pytorch
trainer = torch.optim.SGD
d2l.train_concise_ch11(trainer, {'lr': 0.005, 'momentum': 0.9}, data_iter)
​

#@tab tensorflow
trainer = tf.keras.optimizers.SGD
d2l.train_concise_ch11(trainer, {'learning_rate': 0.005, 'momentum': 0.9},
                       data_iter)
​

#@tab paddle
trainer = paddle.optimizer.Momentum
d2l.train_concise_ch11(trainer, {'learning_rate': 0.005, 'momentum': 0.9}, data_iter)
​

理论分析

f(x) = 0.1 x_1^2 + 2 x_2^2的2D示例似乎相当牵强。
下面我们将看到，它在实际生活中非常具有代表性，至少最小化凸二次目标函数的情况下是如此。

二次凸函数

考虑这个函数

这是一个普通的二次函数。
对于正定矩阵\mathbf{Q} \succ 0，即对于具有正特征值的矩阵，有最小化器为\mathbf{x}^* = -\mathbf{Q}^{-1} \mathbf{c}，最小值为b - \frac{1}{2} \mathbf{c}^\top \mathbf{Q}^{-1} \mathbf{c}。
因此我们可以将h重写为

梯度由\partial_{\mathbf{x}} f(\mathbf{x}) = \mathbf{Q} (\mathbf{x} - \mathbf{Q}^{-1} \mathbf{c})给出。
也就是说，它是由\mathbf{x}和最小化器之间的距离乘以\mathbf{Q}所得出的。
因此，动量法还是\mathbf{Q} (\mathbf{x}_t - \mathbf{Q}^{-1} \mathbf{c})的线性组合。

由于\mathbf{Q}是正定的，因此可以通过\mathbf{Q} = \mathbf{O}^\top \boldsymbol{\Lambda} \mathbf{O}分解为正交（旋转）矩阵\mathbf{O}和正特征值的对角矩阵\boldsymbol{\Lambda}。
这使我们能够将变量从\mathbf{x}更改为\mathbf{z} := \mathbf{O} (\mathbf{x} - \mathbf{Q}^{-1} \mathbf{c})，以获得一个非常简化的表达式：

这里b' = b - \frac{1}{2} \mathbf{c}^\top \mathbf{Q}^{-1} \mathbf{c}。
由于\mathbf{O}只是一个正交矩阵，因此不会真正意义上扰动梯度。
以\mathbf{z}表示的梯度下降变成

这个表达式中的重要事实是梯度下降在不同的特征空间之间不会混合。
也就是说，如果用\mathbf{Q}的特征系统来表示，优化问题是以逐坐标顺序的方式进行的。
这在动量法中也适用。

在这样做的过程中，我们只是证明了以下定理：带有和带有不凸二次函数动量的梯度下降，可以分解为朝二次矩阵特征向量方向坐标顺序的优化。

标量函数

鉴于上述结果，让我们看看当我们最小化函数f(x) = \frac{\lambda}{2} x^2时会发生什么。
对于梯度下降我们有

每|1 - \eta \lambda| < 1时，这种优化以指数速度收敛，因为在t步之后我们可以得到x_t = (1 - \eta \lambda)^t x_0。
这显示了在我们将学习率\eta提高到\eta \lambda = 1之前，收敛率最初是如何提高的。
超过该数值之后，梯度开始发散，对于\eta \lambda > 2而言，优化问题将会发散。

#@tab all
lambdas = [0.1, 1, 10, 19]
eta = 0.1
d2l.set_figsize((6, 4))
for lam in lambdas:
    t = d2l.numpy(d2l.arange(20))
    d2l.plt.plot(t, (1 - eta * lam) ** t, label=f'lambda = {lam:.2f}')
d2l.plt.xlabel('time')
d2l.plt.legend();
​

为了分析动量的收敛情况，我们首先用两个标量重写更新方程：一个用于x，另一个用于动量v。这产生了：

我们用\mathbf{R}来表示2 \times 2管理的收敛表现。
在t步之后，最初的值[v_0, x_0]变为\mathbf{R}(\beta, \eta, \lambda)^t [v_0, x_0]。
因此，收敛速度是由\mathbf{R}的特征值决定的。
请参阅文章 :cite:Goh.2017了解精彩动画。
请参阅 :cite:Flammarion.Bach.2015了解详细分析。
简而言之，当0 < \eta \lambda < 2 + 2 \beta时动量收敛。
与梯度下降的0 < \eta \lambda < 2相比，这是更大范围的可行参数。
另外，一般而言较大值的\beta是可取的。

小结

• 动量法用过去梯度的平均值来替换梯度，这大大加快了收敛速度。
• 对于无噪声梯度下降和嘈杂随机梯度下降，动量法都是可取的。
• 动量法可以防止在随机梯度下降的优化过程停滞的问题。
• 由于对过去的数据进行了指数降权，有效梯度数为\frac{1}{1-\beta}。
• 在凸二次问题中，可以对动量法进行明确而详细的分析。
• 动量法的实现非常简单，但它需要我们存储额外的状态向量（动量\mathbf{v}）。

练习

1. 使用动量超参数和学习率的其他组合，观察和分析不同的实验结果。
2. 试试梯度下降和动量法来解决一个二次问题，其中有多个特征值，即f(x) = \frac{1}{2} \sum_i \lambda_i x_i^2，例如\lambda_i = 2^{-i}。绘制出x的值在初始化x_i = 1时如何下降。
3. 推导h(\mathbf{x}) = \frac{1}{2} \mathbf{x}^\top \mathbf{Q} \mathbf{x} + \mathbf{x}^\top \mathbf{c} + b的最小值和最小化器。
4. 当我们执行带动量法的随机梯度下降时会有什么变化？当我们使用带动量法的小批量随机梯度下降时会发生什么？试验参数如何？

:begin_tab:mxnet
Discussions
:end_tab:

:begin_tab:pytorch
Discussions
:end_tab:

:begin_tab:tensorflow
Discussions
:end_tab:

:begin_tab:paddle
Discussions
:end_tab: