首页 > 其他分享> > 牛顿法进阶

牛顿法进阶

2020-10-08 17:01:49 作者：互联网

专栏目录

广告时间

欢迎关注《阅读经典》系列的其它专题：

编辑于 04-10 数学最优化数值最优化（书籍）

牛顿法进阶

王金戈一切很好，不缺烦恼。杨晓萌等

在之前的文章线搜索与信赖域中，我们介绍过，牛顿法是线搜索的一种，目标函数二阶近似的极值，就是每次迭代牛顿法的方向，用公式表示为

$\nabla^2 f \left( x_k \right) p_k^N = - \nabla f \left( x_k \right) \tag{6.1} \\$

其中， $p_k^N$ 是待求的牛顿方向， $\nabla^2 f \left( x_k \right)$ 是目标函数的二阶导，即海森矩阵， $\nabla f \left( x_k \right)$ 是目标函数的一阶导，即梯度。

牛顿法的显著问题在于，海森矩阵一旦不正定，求出的牛顿方向很可能无法使得目标函数下降。通过一些简单的推导可以得出这一结论。

$\begin{aligned} p_k^N &= - \nabla^2 f_k^{-1} \nabla f_k \\ \nabla f_k^T p_k^N &= - \nabla f_k^T \nabla^2 f_k^{-1} \nabla f_k \\ \end{aligned} \\$

若海森矩阵正定，那么它的逆也正定，上式第二行小于 $0$ ， $p_k^N$ 与 $\nabla f \left( x_k \right)$ 的夹角大于90°，牛顿方向是目标函数下降方向。否则，一旦海森矩阵不正定，牛顿方向就可能是目标函数上升方向，这是我们绝不愿意看到的。

--> 有没有可能，某次上升，全局下降？？

本文将介绍一些解决该问题的方法。

牛顿-共轭梯度法

前两篇文章刚刚介绍过的共轭梯度法就要展示出它的强大实力了，与牛顿法结合，可谓双剑合璧。

既然直接求解牛顿方向可能得到目标函数上升的方向，那还不如用迭代的方式，借助线性共轭梯度法求解式(6.1)。在共轭梯度法（一）：线性共轭梯度中介绍过，当海森矩阵正定时，共轭梯度的每次迭代都会在某个方向上达到极小值，但当海森矩阵不正定时，可能在某个迭代出现负曲率方向，即 $\left( p^{(i)} \right)^T A p^{(i)} \le 0$ ，沿着这一方向是无法找到极小值的，因为此时目标函数的二阶近似单调递减。所以牛顿-共轭梯度法的策略是，一旦遇到负曲率方向就提前终止。我们知道，共轭梯度的第一次迭代是梯度下降方向，所以这种策略最差的结果是在第二次迭代时就终止，牛顿-共轭梯度法退化为最速下降法。虽然收敛速度慢了点，但至少可以保证全局收敛。况且大部分情况下应该都不会那么早终止。

思路很简单，但我们需要思考另一个问题。牛顿-共轭梯度法是否具有与牛顿法相同的二阶收敛速度？

答案是肯定的，在正常情况下，即不存在提前终止的情况下，线性共轭梯度法会在N次迭代后严格收敛到线性方程的解，得到与牛顿法完全相同的结果。不过，在共轭梯度法（一）：线性共轭梯度中，我们提到，事实上不需要N次迭代就可以得到比较精确的解，所以如果每次都跑完N个迭代，会比较浪费时间。在整个牛顿-共轭梯度法的框架下，我们需要一种更方便地控制迭代次数和收敛性的方法。

考虑 $p_k$ 带来的残差，即式(6.1)左右项相减

$r_k = \nabla^2 f \left( x_k \right) p_k + \nabla f \left(x_k \right) \tag{6.2} \\$

$r_k$ 越小，表明共轭梯度法的结果越接近线性方程的解，也就越接近二次收敛。但 $r_k$ 与问题本身的规模有关，而且随着 $k$ 的增加， $r_k$ 本来就会越来越小。所以我们应当考察比值 $\| r_k \| / \| \nabla f \left( x_k \right) \|$ ，令

$\| r_k \| \le \eta_k \| \nabla f \left( x_k \right) \| \tag{6.3} \\$

选择不同的序列 $\left\{ \eta_k \right\}$ ，其中 $0 \lt \eta_k \lt 1$ ，就相当于限制了不同的收敛速度。具体而言，我们可以从式(6.3)推导出 $\eta_k$ 与收敛性的关系，过程省略，结果如下

$\| \nabla f \left( x_{k+1} \right) \| \le \eta_k \| \nabla f \left( x_k \right) \| + O\left( \| \nabla f \left( x_k \right) \|^2 \right) \tag{6.5} \\$

这样一来，我们就可以通过设置 $\eta_k$ 来保证牛顿-共轭梯度法的收敛性。而且我们不一定非要让它二次收敛，通常，我们只需要令 $\eta_k = \min \left( 0.5, \sqrt{\| \nabla f_k \|} \right)$ ，保证其具有超线性收敛性就够了。具体操作的时候，牛顿-共轭梯度法的每次内部迭代都先判断是否满足式(6.3)，满足则提前终止，否则继续迭代。

海森校正

另一种解决海森矩阵不正定问题的方法称为海森校正（Hessian Modification）。思路是让海森矩阵加上一个微小的量，使其正定，相当于用 $B_k = \nabla^2 f \left( x_k \right) + E_k$ 替代 $\nabla^2 f \left( x_k \right)$ 。其中， $E_k$ 需要尽可能小，以避免修改后的 $p_k$ 与牛顿方向 $p_k^N$ 相差太远。

在信赖域（一）：Cauchy Point与Dogleg一文中，我们提到过令 $B_k = \nabla^2 f \left( x_k \right) + \alpha_k I$ ，当 $\alpha_k \gt | \lambda_1 |$ 时， $B_k$ 正定。其中 $\lambda_1$ 是 $\nabla^2 f \left( x_k \right)$ 最小的特征值。因此可以得到 $E_k = \tau I$ ，其中 $\tau = \max \left( 0, \delta - \lambda_1 \right)$ ， $\delta$ 是一个略大于零的常数。

但根据 $\lambda_1$ 计算 $E_k$ 的方法不止一种，对于任意对称矩阵 $A$ ，可以将其分解为 $A = Q\Lambda Q^T$ ，其中 $\Lambda$ 是由 $A$ 的特征值组成的对角矩阵。如果我们令 $A' = Q \left( \Lambda + \text{diag} \left( \tau_i \right) \right) Q^T$ ，其中 $\tau_i$ 由 $\lambda_i$ 确定，并满足 $\tau_{i}=\left\{\begin{array}{ll} {0,} & {\lambda_{i} \geq \delta} \\ {\delta-\lambda_{i},} & {\lambda_{i}<\delta} \end{array}\right.$ 。这意味着直接修改 $A$ 的特征值，将所有小于 $\delta$ 的特征值调整到 $\delta$ ，从而保证 $A'$ 的正定性。该方法相当于令 $E_k = Q \text{diag} \left( \tau_i\right) Q^T$ 。

这两种方法都是根据特征值确定 $E_k$ ，哪种更好呢？

其实说不上哪个更好，本质上，他们是在不同范数约束下的最优结果。前面提到，我们希望 $E_k$ 越小越好，因此可以用矩阵范数来衡量 $E_k$ 的大小。第一种方案其实是在欧氏范数下的最小改动，因为 $\| E_k \| = \| \tau I \| = \tau$ ，不存在使 $\| E_k \|$ 更小的其它 $\tau$ 。第二种方案则是在Frobenius范数下的最小改动。Frobenius范数定义为

$\|A\|_{F}=\left(\sum_{i=1}^{m} \sum_{j=1}^{n} a_{i j}^{2}\right)^{1 / 2} \tag{A.40} \\$

即所有矩阵元素的平方和再开根号，该范数刻画了矩阵所有元素值的大小。在该定义下， $\begin{aligned} \| E_k \|_F &= \| Q \text{diag} \left( \tau_i \right) Q^T \|_F \\ &= \| \text{diag} \left( \tau_i \right) \|_F \\ &= \sqrt{\sum_{i}^{N}{\tau_i^2}} \end{aligned} \\$

第二行去掉了矩阵 $Q$ 和 $Q^T$ ，因为实对称矩阵的特征向量构成的矩阵 $Q$ 是正交矩阵，而正交矩阵相当于旋转，不改变Frobenius范数。容易验证，第一种方案的Frobenius范数为 $\| \tau I \|_F = \sqrt{N} \tau$ ，做个简单的放缩就能发现 $\sqrt{\sum_{i}^{N}{\tau_i^2}} \le \sqrt{N} \tau$ ，所以第二种方案是Frobenius范数下的最小改动。

了解上面两种方法，其实并不能在实际中帮上忙。因为海森矩阵的特征值我们根本无从得知，那是一个需要花大代价才能知道的数字。在实际中，我们一般采取 $E_k = \tau I$ 这种方案，但不通过特征向量确定 $\tau$ ，而是以迭代的方式把可行的 $\tau$ 试出来。

具体怎么做呢？我们可以先从一个稍小的数开始尝试，比如令 $\tau_0 = \| A \|_F / 2$ ，在每次迭代中对校正后的海森矩阵进行Cholesky分解，

$LL^T = A + \tau_k I \\$

如果分解成功，那么分解的结果 $L$ 就可以用来快速求解式(6.1)。如果分解失败，那么将 $\tau_k$ 扩大一倍，进入下一次尝试，直到分解成功。

Modified Cholesky分解

本文到目前为止已经介绍了两类方案，一是用迭代的方式求解式(6.1)，二是对海森矩阵进行校正，即叠加一个小矩阵，使其正定。但这些方法似乎都不够完美，迭代方法依赖于线性共轭梯度，也就带来了共轭梯度本身的一些缺陷，比如需要预处理、受计算机数值精度影响。海森校正苦于无法得知特征值的大小，只能多试几次，试的次数越多，浪费的时间就越多。

但问题难不倒聪明的研究者，有人提出了一种新的海森校正方法，称为modified Cholesky分解，可以在做Cholesky分解的过程中调整原矩阵的值，不需要迭代，不需要试错，分解完毕后直接得到新的正定海森矩阵。让我们看看这种方法是如何实现的。

对于任意一个正定对称矩阵，一定可以分解为如下形式

$A = LDL^T \tag{6.15} \\$

其中， $L$ 是对角线元素全为1的下三角矩阵， $D$ 是对角矩阵。为了直观地解释分解的过程，我们以3×3矩阵为例，式(6.15)展开后得到

$\begin{aligned} \left[\begin{array}{ccc} {a_{11}} & {a_{21}} & {a_{31}} \\ {a_{21}} & {a_{22}} & {a_{32}} \\ {a_{31}} & {a_{32}} & {a_{33}} \end{array}\right] &= \left[\begin{array}{ccc} {1} & {0} & {0} \\ {l_{21}} & {1} & {0} \\ {l_{31}} & {l_{32}} & {1} \end{array}\right]\left[\begin{array}{ccc} {d_{1}} & {0} & {0} \\ {0} & {d_{2}} & {0} \\ {0} & {0} & {d_{3}} \end{array}\right]\left[\begin{array}{ccc} {1} & {l_{21}} & {l_{31}} \\ {0} & {1} & {l_{32}} \\ {0} & {0} & {1} \end{array}\right] \\ &= \left[\begin{array}{ccc} {d_1} & {0} & {0} \\ {l_{21} d_1} & {l_{21}^2 d_1 + d_2} & {0} \\ {l_{31} d_1} & {l_{31} l_{21} d_1 + l_{32} d_2} & {l_{31}^2 d_1 + l_{32}^2 d_2 + d_3} \end{array}\right] \end{aligned} \\$

令等号左右第一列相等，可以得到

$\begin{aligned} d_1 &= a_{11} \\ l_{21} &= a_{21} / d_1 \\ l_{31} &= a_{31} / d_1 \end{aligned} \\$

再令第二列、第三列相等，可以得到

$\begin{aligned} d_2 &= a_{22} - l_{21}^2 d_1 \\ l_{32} &= \left( a_{32} - l_{31} l_{21} d_1 \right) / d_2 \\ d_3 &= a_{33} - l_{31}^2 d_1 - l_{32}^2 d_2 \end{aligned} \\$

如果矩阵 $A$ 不正定，那么计算出的 $d_i$ 很可能是负的，同时导致 $L$ 矩阵的元素值不稳定，很可能过大。也就是说，分解成 $LDL^T$ 形式后再修改 $d_i$ 的值已经来不及了，我们需要在分解的过程中随时修改 $d_i$ 。以上面的例子为例，如果 $d_1$ 和 $d_2$ 都是正数，而 $d_3$ 是个负数，那么我们就直接把 $d_3$ 替换为一个正数。由于 $d_3$ 是最后一步算出来的，所以它的改变不影响其它分解结果。但假设另一种情况， $d_1$ 是正数， $d_2$ 就已经是负数了，这时我们仍然直接把 $d_2$ 替换为正数，不过这会影响下一步 $l_{32}$ 的计算，所以这种情况下 $L$ 矩阵会随之改变。

为什么要这样做呢？回到我们最初的目标，如果 $A$ 不正定，我们想要找到一个正定矩阵 $A' = A + E$ 来替代 $A$ ，这样一来，虽然 $A = LDL^T$ 可能不存在，但我们可以找到使 $A' = L'D'L'^T$ 成立的 $A'$ 。上面的分解步骤相当于找到了一个 $A' = L' \left( D + \text{diag} \left( \tau_i \right) \right) L'^T$ ，其中 $\tau_{i}=\left\{\begin{array}{ll} {0,} & {d_i \geq \delta} \\ {\delta-d_{i},} & {d_{i}<\delta} \end{array}\right.$ 。

到这里还没完，这一方法有一个神奇的性质， $L'D'L'^T = LDL^T + \text{diag} \left( \tau_i \right)$ 恒成立，即 $A' = A + \text{diag} \left( \tau_i \right)$ 。我们在 $D$ 的对角线元素上的改动竟然直接相当于修改 $A$ 的对角线元素，同时保持非对角线元素不变。这个性质非常完美，可以最大限度地保留 $A$ 本身的值，只对某些对角线元素做微小的修改。与上一节介绍的两种方案相比，它的改动最小，也最高效。至于为什么会存在这样的性质，观察上面的推导，可以发现 $A$ 的非对角线元素都包含 $l_{ij} d_j$ 项，而这一项在整个推导过程中其实是不变的。

由于篇幅有限，modified Cholesky的实际操作上的细节我们就不展开讲了，感兴趣的读者可以参考原书。

除了modified Cholesky方法以外，书中还介绍了其它一些在线修改海森矩阵的方法（所谓在线，是指在分解的过程中修改），我们也不再赘述。

信赖域牛顿法

前面介绍的各种方法都是线搜索牛顿法的变形，而且我们关注的是线搜索中如何确定方向，步长则仍然用线搜索（一）：步长的选取中的方法来确定。

既然有线搜索牛顿法，自然也有信赖域牛顿法。其实所谓牛顿法，无非就是取 $B_k = \nabla^2 f_k$ 而已，在线搜索中如此，在信赖域中也如此。曾经介绍过的Dogleg、二维子空间最小化、CG-Steihaug、确切方法等信赖域方法都可以把牛顿法套进去。所以我们不必一一介绍这些方法，只谈谈这些方法需要注意的地方。

对于Dogleg和二维子空间最小化这两种方法来说，它们的近似程度很依赖于无约束极小值点的位置。如果海森矩阵正定，显然不会出问题。但如果海森矩阵不正定，如果用前面介绍的海森校正方法将其变为正定，再计算无约束极小值点，就可能产生较大的偏差。所以这两种方法只适用于海森矩阵正定的情况。

CG-Steihaug方法结合牛顿法后称为信赖域牛顿共轭梯度法，它与本文第一部分介绍的牛顿-共轭梯度法相比，只是多了个信赖域半径。但恰恰就是这个信赖域，会使得该方法优于牛顿-共轭梯度法。回顾前面提到的内容，牛顿-共轭梯度法为了解决海森矩阵不正定的问题，采取提前结束的策略，一旦发现负曲率方向，就结束共轭梯度的迭代，进入下一次线搜索的迭代。这并非什么好策略，顶多算是个不得已而为之的办法。但在信赖域中，我们有更好的方案。由于信赖域半径的存在，即使遇到了负曲率方向，我们可以在终止迭代的同时，沿该方向前进到信赖域边沿，因为我们认为信赖域半径内的区域都是可信区域，而目标函数沿着该负曲率方向又是单调下降的。这种方案在实践中明显优于牛顿-共轭梯度法。

不过需要注意的是，只要使用共轭梯度法，就一定需要预处理，否则收敛性会受到严重影响。

总结

本文介绍了牛顿法应用在线搜索和信赖域中时遇到的实际问题和解决方案。我们越来越发现，各种方法交织在一起，它们可能只有细微的差别，这需要大家对前面介绍的信赖域、线搜索、共轭梯度的基础知识足够熟悉，才能理清这些层出不穷的新方法。

上一篇：共轭梯度法（二）：非线性共轭梯度
下一篇：数值求导和自动求导

编辑于 04-11

标签：进阶,迭代,牛顿,矩阵,海森,梯度,共轭
来源： https://www.cnblogs.com/cx2016/p/13781901.html