Numerical Optimization Ch3. Line Search Methods

第三章: 线搜索算法

线搜索的迭代格式, 如前所述 $$x_{k+1}=x_k+{\alpha }_k{p}_k.$$在 第二章中, 我们提到了常用的搜索方向: 最速下降方向、牛顿方向、拟牛顿方向、共轭方向. 这章我们补齐 关于搜索步长${\alpha }_k$的讨论, 并探讨 ${\alpha }_k,p_k$的选取与 全局收敛性的关系, 证明最速下降法、拟牛顿法和牛顿法的 (局部) 收敛速度. 至于共轭梯度法, 我们放到后面.

1. 搜索步长的选取

在选取步长时, 我们将面临取舍: 既想给予$f$充分的下降, 又不想在做决策时消耗过多. 最理想的步长自当是$$\arg \underset{\alpha }{\min }\varphi (\alpha )=f(x_k+\alpha p_k),\alpha >0.$$但一般来说, 找到$\varphi (\alpha )$的全局极小甚至是局部极小都是非常昂贵的 (需多次获取$f,\nabla f$). 因此, 不精确线搜索往往更加实用.
典型的步长搜索算法通过尝试一系列$\alpha$值, 直至满足停机准则. 一般包括两个阶段:

1. 确定区间 (bracketing phase). 即找到包含满足条件$\alpha$的区间.
2. 二分或插值 (bisection or interpolation phase) . 用于在区间内计算满足条件的$\alpha$.

我们在细化这两个阶段前先来谈谈停机准则. 直观上, ${\alpha }_k$需要能带来$f$的下降, 即$$f(x_k+{\alpha }_k{p}_k)<f(x_k).$$ 但我们有反例说明, 光下降是不够的, 比如

图中函数的最小值为$f^{\ast }=-1$, 但一系列迭代点 { x k } \{x_k\} {xk​}上的函数值为$f(x_k)=5/k,k=0,1,\dots .$ 此例中, 尽管每一步迭代都带来了函数值下降, 但最后却没能收敛到极小点. 为此, 我们需要选出的$\alpha$能带来"充分(sufficient)"下降.

1.1 Wolfe条件

由之前讨论, 不精确线搜索应当保证$f$充分下降, 即$$f(x_k+{\alpha }_k{p}_k)\le f(x_k)+c_1{\alpha }_k\nabla f_k^T{p}_k,c_i\in (0,1).$$也就是说$x_k+\alpha p_k$处的函数值应当位于直线$l(\alpha )=f(x_k)+c_1\alpha \nabla f_k^T{p}_k$下面. 这点用图像说明更加直观. 充分下降条件也称为Armijo条件.

图中$l(\alpha )$以下的部分均是$\alpha$的"合理"采集点. 一般$c_1$取$10^{-4}$. 但充分下降条件似乎不能保证算法较好地收敛, 譬如从上图就可看出, 任一充分小的$\alpha$都能带来所谓的充分下降. 为了避免选中特别小的$\alpha$, 我们再附加一个曲率条件 (curvature condition) :$$\nabla f(x_k+{\alpha }_k{p}_k{)}^T{p}_k\ge c_2\nabla f_k^T{p}_k,c_2\in (c_1,1).$$注意上面不等式的左边就是${\varphi }^{\prime }({\alpha }_k)$, 右边是${\varphi }^{\prime }(0)$. 这一点可以解释成, 若当前的${\varphi }^{\prime }(\alpha )$为"相当负"的数, 那就说明$f$沿着这个方向走更远能够得到更多的下降; 若${\varphi }^{\prime }({\alpha }_k)$不是"那么负"或者干脆是正数, 那也就是说我们不能期望继续沿着当前方向搜索会得到更好的函数值下降, 从而终止搜索. $c_2$一般取0.9或0.1.¹ 为什么称此条件为曲率条件? 个人认为, 它间接地表现了函数在局部极小点附近的凸性.² 关于曲率条件的图示如下.

将充分下降条件和曲率条件合起来就是Wolfe条件:$$\begin{array}{rl}f(x_k+{\alpha }_k{p}_k)& \le f(x_k)+c_1{\alpha }_k\nabla f_k^T{p}_k,\\ \nabla f(x_k+{\alpha }_k{p}_k{)}^T{p}_k& \ge c_2\nabla f_k^T{p}_k,\end{array}$$其中$0<c_1<c_2<1$.
Wolfe条件还可加强为强Wolfe条件, 而这区别在于曲率条件:$$\begin{array}{rl}f(x_k+{\alpha }_k{p}_k)& \le f(x_k)+c_1{\alpha }_k\nabla f_k^T{p}_k,\\ |\nabla f(x_k+{\alpha }_k{p}_k{)}^T{p}_k|& \le c_2|\nabla f_k^T{p}_k|.\end{array}$$这一条件强制${\alpha }_k$位于局部极小点的一个邻域中.³ 此外, 我们不再允许${\varphi }^{\prime }({\alpha }_k)$取正值, 从而排除了距稳定点较远的点.

利用微分中值定理不难证明, 对于光滑有下界的函数$f$, 满足(强)Wolfe条件的${\alpha }_k$总是存在的.
定理1 设$f:{\mathbb{R}}^n\to \mathbb{R}$连续可微, $p_k$为在$x_k$处的一个下降方向. 假设$f$沿直线 { x k + α p k ∣ α &gt; 0 } \{x_k+\alpha p_k|\alpha&gt;0\} {xk​+αpk​∣α>0}有界. 则对于$0<c_1<c_2<1$, 存在满足(强)Wolfe条件的步长 (区间) .
证明: 由条件得$\varphi (\alpha )=f(x_k+\alpha p_k),\alpha >0$下有界. 由于$0<c_1<1$, 且$l(\alpha )=f(x_k)+\alpha c_1\nabla f_k^T{p}_k$下无界, 因此$l(\alpha )$必与$\varphi (\alpha )$交于一点. 设${\alpha }^{\prime }$为第一个交点, 即有$$f(x_k+{\alpha }^{\prime }{p}_k)=f(x_k)+{\alpha }^{\prime }{c}_1\nabla f_k^T{p}_k.$$对于小于${\alpha }^{\prime }$的$\alpha$, 充分下降条件显然成立. 由微分中值定理, 存在${\alpha }^{\prime \prime }\in (0,{\alpha }^{\prime })$使得$$f(x_k+{\alpha }^{\prime }{p}_k)-f(x_k)={\alpha }^{\prime }\nabla f(x_k+{\alpha }^{\prime \prime }{p}_k{)}^T{p}_k.$$组合上面两式, 就有$$\nabla f(x_k+{\alpha }^{\prime \prime }{p}_k{)}^T{p}_k=c_1\nabla f_k^T{p}_k\ge c_2\nabla f_k^T{p}_k.$$这是由于$c_1<c_2$且$\nabla f_k^T{p}_k<0$. 因此, 满足Wolfe条件的步长存在 (由$f$得光滑性, 也必定存在满足Wolfe条件的步长区间). 另外, 由上式左端是负数, 因此同时证明了强Wolfe条件的情形.

1.2 Goldstein条件

如同Wolfe条件, Goldstein条件也保证充分下降以及不会取到过小的步长. 数学表达如下:$$f(x_k)+(1-c){\alpha }_k\nabla f_k^T{p}_k\le f(x_k+{\alpha }_k{p}_k)\le f(x_k)+c{\alpha }_k\nabla f_k^T{p}_k,0<c<1/2.$$图示如下.

如图, 我们发现Goldstein条件可能会排除$\varphi$的所有极小点. 不过两条件还是有很多相似的地方, 包括它们的收敛理论. Goldstein条件常用于牛顿类算法, 但却不适用于拟牛顿算法. 这是因为后者需要保持Hessian逼近的正定性.

1.3 充分下降与回溯

回溯可以在一定程度上替代曲率条件. 它运作的机制是, 从较大的$\bar{\alpha }$开始, 不断乘以缩减因子$\rho$, 直到缩减后的$\alpha$满足充分下降条件. 其最基本的算法形式如下:

算法1 回溯线搜索
Choose $\bar{\alpha }>0,\rho \in (0,1),c\in (0,1)$; Set $\alpha \leftarrow \bar{\alpha }$;
repeat until $f(x_k+\alpha p_k)\le f(x_k)+c\alpha \nabla f_k^T{p}_k$
$\alpha \leftarrow \rho \alpha$
end (repeat)
Terminate with ${\alpha }_k=\alpha$.

对于牛顿或拟牛顿法, 初始$\bar{\alpha }$选取为1, 而对于其他诸如最速下降法的算法, 选取各有不同标准. 缩减因子$\rho$在每步迭代也可改变, 事实上只需保证$\rho \in [{\rho }_{lo},{\rho }_{hi}]$, 其中${\rho }_{lo},{\rho }_{hi}$分别为区间下界和上界.
这一简单的搜索策略适用于牛顿法, 但不适用于拟牛顿法和共轭梯度法. 这在后面会进一步讨论.

2. 线搜索的收敛性

本节讨论线搜索算法的全局收敛性. 这里要用到在第二章中定义的, $p_k$与负梯度$-\nabla f_k$的夹角${\theta }_k$, 有$$\cos {\theta }_k=\frac{-\nabla f_k^T{p}_k}{\Vert \nabla f_k\Vert \Vert p_k\Vert }.$$下面叙述并证明著名的Zoutendijk定理, 其中提出Zoutendijk条件. 该条件可用于证明广泛算法的全局收敛性. 从定理的叙述和证明的过程中, 我们可以体会到${\alpha }_k$和$p_k$选取的重要性.

定理2 (Zoutendijk定理) 考虑任一一种以$$x_{k+1}=x_k+{\alpha }_k{p}_k$$迭代的算法, 其中$p_k$为下降方向, ${\alpha }_k$满足Wolfe条件. 设$f$在${\mathbb{R}}^n$下有界, 且$f$在包含水平集 L = { x : f ( x ) ≤ f ( x 0 ) } \mathcal{L}=\{x:f(x)\le f(x_0)\} L={x:f(x)≤f(x0​)}的一个开集$\mathcal{N}$中连续可微, 其中$x_0$为迭代的初始点. 设$\nabla f$在$\mathcal{N}$上Lipschitz连续, 即存在$L>0$使得$$\Vert \nabla f(x)-\nabla f(\tilde{x})\Vert \le L\Vert x-\tilde{x}\Vert ,\forall x,\tilde{x}\in \mathcal{N}.$$则该算法满足Zoutendijk条件$$\sum _{k\ge 0}{\cos }^2{\theta }_k\Vert \nabla f_k{\Vert }^2<\infty .$$证明: 由曲率条件得$$(\nabla f_{k+1}-\nabla f_k{)}^T{p}_k\ge (c_2-1)\nabla f_k^T{p}_k,$$而由Lipschitz条件可得$$(\nabla f_{k+1}-\nabla f_k{)}^T{p}_k\le {\alpha }_{k}L\Vert p_k{\Vert }^2.$$结合上面两式, 得到$${\alpha }_k\ge \frac{c_2-1}{L}\frac{\nabla f_k^T{p}_k}{\Vert p_k{\Vert }^2}.$$将此不等式代入充分下降条件, 推出$$f_{k+1}\le f_k+c_1{\alpha }_k\nabla f_k^T{p}_k\le f_k-c_1\frac{1-c_2}{L}\frac{(\nabla f_k^T{p}_k{)}^2}{\Vert p_k{\Vert }^2}=f_k-c{\cos }^2{\theta }_k\Vert \nabla f_k{\Vert }^2.$$其中$c=c_1(1-c_2)/L$. 上面的不等式对$k$做累和, 得到$$f_{k+1}\le f_0-c\sum _{j=0}^k{\cos }^2{\theta }_j\Vert \nabla f_j{\Vert }^2.$$因为$f$下有界, 所以$f_0-f_{k+1}$对于$k$有一致上界 (小于正无穷) . 因此令$k\to \infty$, Zoutendijk条件成立:$$\sum _{k=0}^{\infty }{\cos }^2{\theta }_k\Vert \nabla f_k{\Vert }^2<\infty .$$
对于Goldstein条件和强Wolfe条件也有类似的定理, 且最终都会满足Zoutendijk条件. 注意Zoutendijk定理的条件并不苛刻: $f$下有界可以说是优化的必要条件; 而$\nabla f$的Lipschitz连续性则对体现了函数的光滑性, 实际应用时一般也是满足的.
特别地, 从Zoutendijk条件我们可以得到对于线搜索算法的一类全局收敛判定方法. 我们先定义什么叫做全局收敛.

定义1 (全局收敛) 若迭代算法满足$$\underset{k\to \infty }{\lim }\Vert \nabla f_k\Vert =0,$$则称算法是全局收敛的.

算法的全局收敛并非指字面义上的收敛到全局极小点. 从定义上看, 事实上指的是收敛到稳定点. 只有在$p_k$上加上二阶信息, 我们才有可能使得算法收敛到局部极小点. 回过头来看Zoutendijk条件:$$\sum _{k=0}^{\infty }{\cos }^2{\theta }_k\Vert \nabla f_k{\Vert }^2<\infty .$$若$\cos {\theta }_k\ge \delta >0,\forall k$, 则直接由Zoutendijk条件得到全局收敛性. 这就是说, 只要我们能够保证

1. 算法产生的搜索方向$p_k$与$-\nabla f_k$的正交方向有一致的距离;
2. 找寻的${\alpha }_k$满足Wolfe条件,

那么算法就是全局收敛的. 注意我们这只用到了级数收敛的必要条件, 还未将有限和条件纳入. 我们将利用它来证明序列 { cos ⁡ 2 θ k ∥ ∇ f k ∥ 2 } \{\cos^2\theta_k\Vert\nabla f_k\Vert^2\} {cos2θk​∥∇fk​∥2}的快速收敛. 下面我们来讨论之前叙述的算法的全局收敛性.

1. 最速下降法显然是全局收敛的, 因为夹角${\theta }_k\equiv 0$.
2. 牛顿类算法. 假设模型中的$B_k$ (牛顿法中是Hessian矩阵, 拟牛顿法中是Hessian矩阵的近似)正定且条件数一致有界, 即存在$M>0$, $$\Vert B_k\Vert \Vert B_k^{-1}\Vert \le M,\forall k.$$从这一点可以证得在2-范数下 (限于本人水平, 未证出其他范数下是否成立) , $\cos {\theta }_k\ge \frac{1}{M},\forall k$. 事实上, $$\begin{array}{rl}\cos {\theta }_k=\frac{-\nabla f_k^T{p}_k}{\Vert \nabla f_k\Vert \Vert p_k\Vert }& =\frac{\nabla f_k^T{B}_k^{-1}\nabla f_k}{\Vert B_k^{-1}\nabla f_k\Vert \Vert \nabla f_k\Vert }\\ & \ge \frac{\nabla f_k^T{B}_k^{-1}\nabla f_k}{\Vert B_k^{-1}\nabla f_k\Vert \Vert B_k\Vert \Vert B_k^{-1}\nabla f_k\Vert }(\because \Vert \nabla f_k\Vert \le \Vert B_k\Vert \Vert B_k^{-1}\nabla f_k\Vert ,范数相容性)\\ & \ge \frac{(\nabla f_k^T{B}_k^{-1}\nabla f_k)\Vert B_k^{-1}\Vert }{\nabla f_k^T{B}_k^{-1}{B}_k^{-1}\nabla f_k}\cdot \frac{1}{M}.(\because 条件数的一致有界性)\end{array}$$设$B_k$的特征值分解为$B_k=U_k{\Sigma }_k{U}_k^{-1}$, 其中$U$为正交矩阵, $\Sigma$为对角矩阵, 对角元为${\lambda }_1^{(k)}\ge \cdots \ge {\lambda }_n^{(k)}>0$. 记$$\nabla f_k^T{U}_k=(z_1,\dots ,z_n)̸=0$$且根据2-范数的性质, $\Vert B_k^{-1}\Vert$就是$B_k^{-1}$的最大特征(奇异)值, 即$\Vert B_k^{-1}\Vert =1/{\lambda }_n^{(k)}$. 所以上面不等式的最后右端项为$$\frac{(\nabla f_k^T{B}_k^{-1}\nabla f_k)\Vert B_k^{-1}\Vert }{\nabla f_k^T{B}_k^{-1}{B}_k^{-1}\nabla f_k}\cdot \frac{1}{M}=\frac{\left(\sum _{i=0}^n\frac{1}{{\lambda }_i^{(k)}}{z}_i^2\right)\frac{1}{{\lambda }_n^{(k)}}}{\sum _{i=0}^n{\left(\frac{1}{{\lambda }_i^{(k)}}\right)}^2{z}_i^2}\cdot \frac{1}{M}\ge \frac{1}{M},\forall k.$$根据之前的推导就可得全局收敛性.
3. 共轭梯度法. 此时我们只能证明弱全局收敛性 , 即$$\underset{k\to \infty }{\mathrm{liminf}}\Vert \nabla f_k\Vert =0.$$这也就是说, 我们只能证明有梯度范数序列的一个子列 { ∥ ∇ f k j ∥ } \{\Vert\nabla f_{k_j}\Vert\} {∥∇fkj​​∥}收敛. 与之前的正面证明不同, 我们将在后面用反证法证之. 那么所谓的"弱全局收敛性"有什么用呢?事实上, 我们可以从定理的证明中得到一些规律: 简洁的定理难证, 复杂的定理好证. 所谓的"弱全局收敛性"一般要比全局收敛性要好证得多, 而在实际应用中, 我们以梯度范数序列的收敛作为算法收敛的判定条件. 若我们证明了"弱全局收敛性", 那就说明存在梯度范数序列的子列收敛, 算法最终必定是收敛的; 另一方面, "弱全局收敛"与全局收敛的应用界限并没有那么明确. 我们往往可以加上较为实际的条件达成全局收敛.

至此, 我们已经证明了一些通用算法的全局收敛性. 事实上, 为了保证算法的弱全局收敛, 我们可以这样构造线搜索算法:

1. 每步迭代都产生函数值的下降;
2. 每$m$步迭代就采用一次负梯度作为搜索方向, 并在满足Wolfe条件或Goldstein条件时选取步长${\alpha }_k$.

这样尽管偶尔的最速下降不一定能产生多大的下降, 但却可以保证"总体"上的收敛性 (因为有子列收敛) . 当然我们可以在非最速下降的迭代步上做一些功夫, 加速收敛.

3. 收敛速度

从之前的讨论我们可以得出结论: 设计全局收敛的算法似乎是很简单的. 我们仅仅需要搜索方向$p_k$不会与负梯度$-\nabla f_k$垂直就行了, 而这点完全可以在每次生成搜索方向时附加角度检验 (angle test) 达到. 先给定一个阈值 (threshold) $\delta >0$, 如果$\cos {\theta }_k\ge \delta$就接受这一搜索方向; 如果$\cos {\theta }_k<\delta$就不妨直接采用最速下降方向. 乍一看这样的设计天衣无缝. 事实上, 我们有两点原因使我们不去这么做:

1. 能得到较快收敛速度的$p_k$往往就与负梯度$-\nabla f_k$差不多垂直, 这点在病态问题上尤其突出. 并且$\delta$的不恰当选取也会强行"过滤"掉一些好的搜索方向;
2. 角度检验可能会影响拟牛顿算法的不变性质.

全局收敛性与局部的收敛速度是矛盾的, 我们需要根据问题找平衡 (tradeoffs) . 一方面, 具有全局收敛性的算法可能收敛速度不快, 比如下面要证实的, 最速下降法保证全局收敛却只有线性收敛速度; 另一方面, 具有较快收敛速度的算法不一定具有全局收敛性. 譬如上一节提到的牛顿法的全局收敛性, 是基于Hessian矩阵正定且条件数一致有界得到的. 换句话说若$x^{\ast }$是局部极小点, $x^{\ast }$处的Hessian正定且 (比如说) Lipschitz连续, 则在离$x^{\ast }$充分近的$x_k$的确可以保证这些条件成立, 也即其收敛是局部性质. 若$x_k$离$x^{\ast }$较远, 我们无法证实其全局收敛, 甚至可能产生的$p_k$都不是下降方向 (当然这可以通过做矩阵修正避免) . 因此真正的挑战在于: 设计算法, 使得它既有全局收敛性, 也有较好的收敛速度.

我们先来谈谈之前提到的算法的收敛速度如何.

3.1 最速下降法的收敛速度

先考虑目标函数为二次函数的特殊情形, 即$$f(x)=\frac{1}{2}{x}^{T}Qx-b^{T}x,$$其中$Q$是对称正定矩阵. 则全局极小点$x^{\ast }$就是线性系统$Ax=b$的唯一解. 下面计算步长${\alpha }_k$. 函数$$f(x_k-\alpha \nabla f_k)=\frac{1}{2}(x_k-\alpha \nabla f_k{)}^{T}Q(x_k-\alpha \nabla f_k)-b^T(x_k-\alpha \nabla f_k)$$对$\alpha$求导置零可得$${\alpha }_k=\frac{\nabla f_k^T\nabla f_k}{\nabla f_k^{T}Q\nabla f_k}.$$如果我们使用${\alpha }_k$的精确表达式 (即精确线搜索步长) , 我们有迭代式$$x_{k+1}=x_k-\left(\frac{\nabla f_k^T\nabla f_k}{\nabla f_k^{T}Q\nabla f_k}\right)\nabla f_k.$$由多元微积分的知识, 我们知道最速下降法的搜索方向与函数曲面等高线垂直, 因此在搜索时会出现锯齿现象 (zigzag) . 直观上, 这种现象影响了它的收敛速度, 尤其是到了迭代后期接近$x^{\ast }$时.
为量化其收敛速度, 设加权范数$\Vert x{\Vert }_Q^2=x^{T}Qx$. 所以$$\begin{array}{rl}\frac{1}{2}\Vert x-x^{\ast }{\Vert }_Q^2& =\frac{1}{2}{x}^{T}Qx+\frac{1}{2}(x^{\ast }{)}^{T}Q{x}^{\ast }-x^{T}Q{x}^{\ast }\\ & =\frac{1}{2}{x}^{T}Qx+\frac{1}{2}{b}^T{x}^{\ast }-b^{T}x(\because Q{x}^{\ast }=b)\\ & =\frac{1}{2}{x}^{T}Qx-b^{T}x-\left(\frac{1}{2}(x^{\ast }{)}^{T}Q{x}^{\ast }-b^T{x}^{\ast }\right)\\ & =f(x)-f(x^{\ast }).\end{array}$$因此这一范数度量了函数值的差. 代入迭代公式, 并注意$\nabla f_k=Q(x_k-x^{\ast })$我们有$$\Vert x_{k+1}-x^{\ast }{\Vert }_Q^2=\left\{1-\frac{(\nabla f_k^T\nabla f_k{)}^2}{(\nabla f_k^{T}Q\nabla f_k)(\nabla f_k^T{Q}^{-1}\nabla f_k)}\right\}\Vert x_k-x^{\ast }{\Vert }_Q^2.$$上面式子进一步可以得到$$\Vert x_{k+1}-x^{\ast }{\Vert }_Q^2\le {\left(\frac{{\lambda }_n-{\lambda }_1}{{\lambda }_n+{\lambda }_1}\right)}^2\Vert x_k-x^{\ast }{\Vert }_Q^2,$$其中$0<{\lambda }_1\le \cdots \le {\lambda }_n$是$Q$的特征值. 关于上面两个式子的证明可见David Luenberger所著Introduction to Linear and Nonlinear Programming. 下面我们从另一个角度证明. 先证明一个引理.

引理1 设$P(t)$是$t$的一个多项式, 则$$\Vert P(Q)x{\Vert }_Q\le \underset{1\le i\le n}{\max }|P({\lambda }_i)|\Vert x{\Vert }_Q,x\in {\mathbb{R}}^n.$$证明: 设$u_1,\dots ,u_n$是$Q$对应于${\lambda }_1,\dots ,{\lambda }_n$的特征向量, 它们构成${\mathbb{R}}^n$的一组标准正交基, 则对任一$x\in {\mathbb{R}}^n$, 有$x={\sum }_{i=1}^n{\alpha }_i{u}_i$, 从而$$\begin{array}{rl}\Vert P(Q)x{\Vert }_Q^2=x^{T}P(Q)QP(Q)x& ={\left(\sum _{i=1}^n{\alpha }_{i}P({\lambda }_i)u_i\right)}^{T}Q\left(\sum _{i=1}^n{\alpha }_{i}P({\lambda }_i)u_i\right)\\ & =\sum _{i=1}^n{\lambda }_i{\alpha }_i^2{P}^2({\lambda }_i)\le \underset{1\le i\le n}{\max }{P}^2({\lambda }_i)\sum _{i=1}^n{\lambda }_i{\alpha }_i^2\\ & =\underset{1\le i\le n}{\max }{P}^2({\lambda }_i)\Vert x{\Vert }_Q^2.\end{array}$$因此得证. 下证误差估计式. $$\begin{array}{rl}\Vert x_{k+1}-x^{\ast }{\Vert }_Q^2& =(x_{k+1}-x^{\ast })Q(x_{k+1}-x^{\ast })\\ & \le (x_k-\alpha \nabla f_k-x^{\ast }{)}^{T}Q(x_k-\alpha \nabla f_k-x^{\ast })\\ & =[(I-\alpha Q)(x_k-x^{\ast }){]}^{T}Q[(I-\alpha Q)(x_k-x^{\ast })]\\ & =\Vert (I-\alpha Q)(x_k-x^{\ast }){\Vert }_Q^2.\end{array}$$记$P_{\alpha }(t)=1-\alpha t$, 利用引理可得$$\begin{array}{rl}\Vert x_{k+1}-x^{\ast }{\Vert }_Q^2& \le \Vert P_{\alpha }(Q)(x_k-x^{\ast }){\Vert }_Q^2\\ & \le \underset{1\le i\le n}{\max }{P}_{\alpha }^2({\lambda }_i)\Vert x_k-x^{\ast }{\Vert }_Q^2\end{array}$$对$\forall \alpha$成立. 利用Chebyshev多项式的性质, 可得$$\underset{\alpha }{\min }\underset{1\le i\le n}{\max }|1-\alpha {\lambda }_i|=\frac{{\lambda }_n-{\lambda }_1}{{\lambda }_n+{\lambda }_1}.$$代入即得$$\Vert x_{k+1}-x^{\ast }{\Vert }_Q^2\le {\left(\frac{{\lambda }_n-{\lambda }_1}{{\lambda }_n+{\lambda }_1}\right)}^2\Vert x_k-x^{\ast }{\Vert }_Q^2.$$从此误差估计可得:

1. 最速下降法线性收敛. 特别地, 当所有的特征值相同时, 收敛一步完成, 此时$Q$为纯量阵. 等高线呈圆形, 而负梯度则直接指向圆心;
2. 由于$\frac{{\lambda }_n-{\lambda }_1}{{\lambda }_n+{\lambda }_1}=\frac{{\lambda }_n/{\lambda }_1-1}{{\lambda }_n/{\lambda }_1+1}=\frac{{\kappa }_2(Q)-1}{{\kappa }_2(Q)+1},$ 其中${\kappa }_2(\cdot )$指代2-范数下的条件数, 因此当$Q$的最大特征值与最小特征值差得越大, 速度就越慢, 极端情况下会出现严重的锯齿现象. 因此最速下降法常作为头几步的运行算法.

对于一般的非线性函数, 我们在假设${\alpha }_k$为精确线搜索步长时, 有
定理3 假设$f:{\mathbb{R}}^n\to \mathbb{R}$二次连续可微, 最速下降法中的采用精确步长线搜索, 迭代项收敛到$x^{\ast }$, 且Hessian矩阵${\nabla }^{2}f(x^{\ast })$正定. 令$r$为任一满足$$r\in \left(\frac{{\lambda }_n-{\lambda }_1}{{\lambda }_n+{\lambda }_1},1\right)$$的标量, 其中${\lambda }_1\le \cdots \le {\lambda }_n$为${\nabla }^{2}f(x^{\ast })$的特征值. 则对于充分大的$k$, 我们有$$f(x_{k+1}-f(x^{\ast })\le r^2[f(x_k)-f(x^{\ast })].$$
一般说来, 我们无法保证不精确线搜索步长下的收敛速度.

3.2 牛顿法

3.2.1 牛顿法的收敛速度

基于之前的讨论, 由于Hessian矩阵${\nabla }^2{f}_k$不一定正定, 因此产生的搜索方向$p_k$不一定是下降方向. 这里我们先假设${\nabla }^{2}f(x)$在$x^{\ast }$附近连续 (或Lipschitz连续) , ${\nabla }^{2}f(x^{\ast })$正定从而$x^{\ast }$的某个邻域内恒有${\nabla }^{2}f(x)$正定. 我们将证明牛顿法具有局部二次收敛性, 步长${\alpha }_k\equiv 1$.

定理4 假设$f$二次可微, Hessian矩阵${\nabla }^{2}f(x)$在$x^{\ast }$的一个邻域内Lipschitz连续, 在$x^{\ast }$满足二阶充分条件. 考虑迭代式$x_{k+1}=x_k+p_k$, 其中$p_k$为牛顿步. 则

1. 若初始点$x_0$充分靠近$x^{\ast }$, 则迭代序列收敛到$x^{\ast }$;
2. 序列 { x k } \{x_k\} {xk​}二次收敛;
3. 范数序列 { ∥ ∇ f k ∥ } \{\Vert\nabla f_k\Vert\} {∥∇fk​∥}二次收敛到0.

证明: 证明的过程中将大量应用Taylor定理. 由最优性条件$\nabla f_{\ast }=0$有$$\begin{array}{rl}{x}_k+p_k-x^{\ast }& =x_k-x^{\ast }-{\nabla }^2{f}_k^{-1}\nabla f_k\\ & ={\nabla }^2{f}_k^{-1}[{\nabla }^2{f}_k(x_k-x^{\ast })-\nabla f_k]\\ & ={\nabla }^2{f}_k^{-1}[{\nabla }^2{f}_k(x_k-x^{\ast })-(\nabla f_k-\nabla f_{\ast })].\end{array}$$由Taylor定理, 有$$\nabla f_k-\nabla f_{\ast }={\int }_0^1{\nabla }^{2}f(x_k+t(x_k-x^{\ast }))(x_k-x^{\ast })\mathrm{d}t.$$代入有$$\begin{array}{rl}& \Vert {\nabla }^2{f}_k(x_k-x^{\ast })-(\nabla f_k-\nabla f_{\ast })\Vert \\ & =\Vert {\int }_0^1[{\nabla }^2{f}_k-{\nabla }^{2}f(x_k+t(x_k-x^{\ast }))](x_k-x^{\ast })\mathrm{d}t\Vert \\ & \le \Vert x_k-x^{\ast }{\Vert }^2{\int }_0^{1}Lt\mathrm{d}t=\frac{1}{2}L\Vert x_k-x^{\ast }{\Vert }^2,\end{array}$$其中$L$为${\nabla }^{2}f(x)$在$x^{\ast }$附近的Lipschitz常数. 由于${\nabla }^{2}f(x^{\ast })$非奇异, 因此$\exists r>0$, 对$\forall x_k:\Vert x_k-x^{\ast }\Vert \le r$, 有$\Vert {\nabla }^2{f}_k^{-1}\Vert \le 2\Vert {\nabla }^{2}f(x^{\ast }{)}^{-1}\Vert$. 从而$$\Vert x_k+p_k-x^{\ast }\Vert \le \frac{1}{2}L\Vert {\nabla }^{2}f(x^{\ast }{)}^{-1}\Vert \Vert x_k-x^{\ast }{\Vert }^2\triangleq \tilde{L}\Vert x_k-x^{\ast }{\Vert }^2,$$其中$\tilde{L}=L\Vert {\nabla }^{2}f(x^{\ast }{)}^{-1}\Vert$. 因此选取$x_0$使得$\Vert x_0-x^{\ast }\Vert \le \min \left(r,1/\tilde{L}\right)$即得序列 { x k } \{x_k\} {xk​}二次收敛于$x^{\ast }$.
再讨论梯度范数序列 { ∥ ∇ f k ∥ } \{\Vert\nabla f_k\Vert\} {∥∇fk​∥}的收敛性. $$\begin{array}{rl}\Vert \nabla f_{k+1}\Vert & =\Vert \nabla f_{k+1}-\nabla f_k-{\nabla }^{2}f(x_k)p_k\Vert (\because \nabla f_k+{\nabla }^2{f}_k{p}_k=0)\\ & =\Vert {\int }_0^1({\nabla }^{2}f(x_k+t{p}_k)-{\nabla }^{2}f(x_k))p_k\mathrm{d}t\Vert \\ & \le {\int }_0^1\Vert {\nabla }^{2}f(x_k+t{p}_k)-{\nabla }^{2}f(x_k)\Vert \Vert p_k\Vert \mathrm{d}t\\ & \le \frac{1}{2}L\Vert p_k{\Vert }^2\\ & \le \frac{1}{2}L\Vert {\nabla }^2{f}_k^{-1}{\Vert }^2\Vert \nabla f_k{\Vert }^2\\ & \le 2L\Vert {\nabla }^{2}f(x^{\ast }{)}^{-1}{\Vert }^2\Vert \nabla f_k{\Vert }^2.\end{array}$$因此梯度范数序列 { ∥ ∇ f k ∥ } \{\Vert\nabla f_k\Vert\} {∥∇fk​∥}二次收敛性于0.

3.2.2 带Hessian修正的牛顿法 (Cont’d)

本小节主要讨论Hessian非正定情形下的修正牛顿法. 事实上, 对于Hessian非正定情形, 有两种可以获得全局收敛的处理方案:

1. 基于修正牛顿法的线搜索. 此时主要对Hessian矩阵${\nabla }^2{f}_k$作必要的修正使它正定, 从而能产生下降方向;
2. 信赖域-牛顿法. 与线搜索不同的是, 信赖域可直接纳入非正定的Hessian进行计算. 此法将放在下一章讨论信赖域时提到.

修正牛顿算法的一般步骤为:
算法2 (基于修正牛顿法的线搜索)
给定初始点$x_0$;
for $k=0,1,2,\dots$
$分解矩阵B_k={\nabla }^{2}f(x_k)+E_k,其中$
$若{\nabla }^{2}f(x_k)充分正定,则E_k=0;否则选取E_k使得B_k充分正定.$
$求解B_k{p}_k=-\nabla f_k$
$x_{k+1}\leftarrow x_k+{\alpha }_k{p}_k$
$其中{\alpha }_k满足Wolfe条件、Goldstein条件或Armijo回溯条件$
end

算法2在$E_k$选取满足有界修正分解性质时, 可以证明是全局收敛的, 即: 只要Hessian矩阵序列 { ∇ 2 f ( x k ) } \{\nabla^2f(x_k)\} {∇2f(xk​)}有界, 就有序列 { B k } \{B_k\} {Bk​}的条件数有界:$$\kappa (B_k)=\Vert B_k\Vert \Vert B_k^{-1}\Vert \le C,C>0,k=0,1,2,\dots .$$再来考虑算法2的收敛速度. 分两种情形.

1. 序列 { x k } \{x_k\} {xk​}收敛到$x^{\ast }$, 且${\nabla }^{2}f(x^{\ast })$充分正定. 则由${\nabla }^{2}f(x)$的连续性 (或Lipschitz) , 在后期$E_k=0$从而恢复到原本的牛顿法, 达到二次收敛.
2. ${\nabla }^{2}f(x^{\ast })$ (接近) 奇异. 此时$E_k$将不一定消逝, 且算法的收敛速度可能只有线性. 我们首先需要$B_k$的条件数适当 (从而满足有界修正分解性质), 其次还需要修正项$E_k$尽可能的小从而尽量保留二阶信息. 同时计算的耗费也要在可接受的范围.

3.2.2.1 基于特征值的修正

让我们从一个例子说起. 设当前迭代点为$x_k,\nabla f_k=(1,-3,2{)}^T,{\nabla }^2{f}_k=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(10,3,-1{)}^T$. 显然二阶信息是不定的. 由特征值分解, 我们有$Q=I,\Lambda =\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}({\lambda }_1,{\lambda }_2,{\lambda }_3)$, 使得$${\nabla }^2{f}_k=Q\Lambda Q^T=\sum _{i=1}^n{\lambda }_i{q}_i{q}_i^T.$$一个直接的想法是, 将对角元中的负元用较小的正数$\delta$替代, 譬如说$\delta =\sqrt{\mathrm{u}}$, 其中$\mathrm{u}$是机器精度 (设为$10^{-16}$). 则此例中$$B_k=\sum _{i=1}^2{\lambda }_i{q}_i{q}_i^T+\delta q_3{q}_3^T=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}(10,3,10^{-8}).$$不过, 基于此修正Hessian的搜索方向为$$p_k=-B_k^{-1}\nabla f_k=-\sum _{i=1}^2\frac{1}{{\lambda }_i}{q}_i(q_i^T\nabla f_k)-\frac{1}{\delta }{q}_3(q_3^T\nabla f_k)\approx -(2\times 10^8)q_3.$$对于很小的$\delta$, 这一搜索方向将近与$q_3$平行, 并且长度惊人. 这有悖牛顿法的想法.
不过受此启发, 很多人提出了其他的修正方案. 例如干脆将负元变个符号. 在我们的例子中就是$\delta =1$. 或者干脆把负曲率的项去掉, 在我们的例子中就是$B_k={\sum }_{i=1}^2{\lambda }_i{q}_i{q}_i^T$. 再比如, 选取$\delta$使得$p_k$看起来合理一些.

基于最小F-范数的修正是说, 给定$\delta >0$, 求 Δ A = arg ⁡ min ⁡ B { ∥ B ∥ F ∣ λ min ⁡ ( A + B ) ≥ δ } . \Delta A=\arg\min_B\{\Vert B\Vert_F|\lambda_{\min}(A+B)\ge\delta\}. ΔA=argBmin​{∥B∥F​∣λmin​(A+B)≥δ}.最后得到的结果是$$\Delta A=Q(\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}({\tau }_i))Q^T,{\tau }_i=\{\begin{array}{ll}0,& {\lambda }_i\ge \delta ,\\ \delta -{\lambda }_i,& {\lambda }_i<\delta .\end{array}$$这与第一种想法契合.

基于最小2-范数的修正则会给出对角修正. 给定$\delta >0$, 求 Δ A = arg ⁡ min ⁡ B { ∥ B ∥ 2 ∣ λ min ⁡ ( A + B ) ≥ δ } . \Delta A=\arg\min_B\{\Vert B\Vert_2|\lambda_{\min}(A+B)\ge\delta\}. ΔA=argBmin​{∥B∥2​∣λmin​(A+B)≥δ}.最后得到的结果是$$\Delta A=\tau I,\tau =\max (0,\delta -{\lambda }_{\min }(A)).$$实际应用时往往采用试探法, 不断增大$\tau$直至Cholesky分解能够完成. 这与信赖域的修正有些类似.

实际中, 往往不直接做特征值分解, 因为这一步骤的计算量过于庞大. 利用Gauss消去间接地修正往往更受欢迎.

3.2.2.2 修正Cholesky分解

另一种修正Hessian矩阵的方案, 是在Cholesky分解的过程中进行的. 具体说, 如果分解的过程中计算出了负的 (准确说是虚的) 对角元, 则增加一定量保证充分正定. 修正Cholesky分解的应当做到:

1. 保证修正的Cholesky因子存在;
2. 修正后的矩阵条件数一致有界;
3. 若Hessian本身充分正定则不做修正.

设A的LDL分解为$$A=LD{L}^T.$$其中$L$为下三角矩阵, $D$为对角阵. 相应的算法为

算法3 (LDL分解)
for $j=1,2,\dots ,n$
$c_{jj}\leftarrow a_{jj}-{\sum }_{s=1}^{j-1}{d}_s{l}_{js}^2$;
$d_j\leftarrow c_{jj}$;
\quad\quad for $i=j+1,\dots ,n$
$c_{ij}\leftarrow a_{ij}-{\sum }_{s=1}^{j-1}{d}_s{l}_{is}{l}_{js}$;
$l_{ij}\leftarrow c_{ij}/d_j$;
\quad\quad end
end

当$A$不定时, LDL分解是不稳定的. 下面对其加修正. 为了控制修正的质量, 选取两正数$\delta ,\beta$, 要求在计算$L,D$的第$j$列时, 满足:$$d_j\ge \delta ,|m_{ij}|\le \beta ,i=j+1,\dots ,n,$$其中$m_{ij}=l_{ij}\sqrt{d_j}$. 为保证上述不等式, 我们只需在算法3中修改一步: 计算对角元$d_j$的公式改为$$d_j=\max \left(|c_{jj}|,{\left(\frac{{\theta }_j}{\beta }\right)}^2,\delta \right),{\theta }_j=\underset{j<i\le n}{\max }|c_{ij}|.$$下面验证不等式成立. 注意$c_{ij}=l_{ij}{d}_j$, 因此$$|m_{ij}|=|l_{ij}\sqrt{d_j}|=\frac{|c_{ij}|}{\sqrt{d_j}}\le \frac{|c_{ij}|\beta }{{\theta }_j}\le \beta ,\forall i>j.$$注意${\theta }_j$可以先于$d_j$得到, 这就是算法中$c_{ij}$的作用.

3.2.2.3 修正对称不定分解

任一对称矩阵$A$, 不论正定或不正定, 都有$$PA{P}^T=LB{L}^T,$$其中$L$为单位下三角矩阵, $B$为块对角矩阵, 其中块为1阶或2阶, $P$为排列矩阵. 我们之前提到过直接计算不定矩阵的LDL分解不稳定, 这是因为$L,D$中可能含有远大于$A$中元素的项. 但$B$为块对角矩阵, 其对角元为1阶或2阶矩阵, 这样的分解总是存在且计算稳定.

3.3 拟牛顿法的收敛速度

拟牛顿法的搜索方向为$$p_k=-B_k^{-1}\nabla f_k,$$其中对称正定矩阵$B_k$在每步都被更新. 我们假定${\alpha }_k$满足Wolfe条件或强Wolfe条件, 且步长线搜索必首先尝试$\alpha =1$. 这一点在实现快速收敛时将起到关键的作用.
下面的定理说明, 若拟牛顿法的搜索方向与牛顿法的搜索方向足够接近, 则随着迭代点向解逼近, 单位步长将最终满足Wolfe条件. 定理同时还给出了拟牛顿法搜索方向必须遵循的条件, 以达成超线性收敛.

定理5 设$f:{\mathbb{R}}^n\to \mathbb{R}$二阶连续可微. 考虑迭代$$x_{k+1}=x_k+{\alpha }_k{p}_k,$$其中$p_k$为下降方向, ${\alpha }_k$满足Wolfe条件 ($c_1\le 1/2$). 若序列 { x k } \{x_k\} {xk​}收敛到$x^{\ast }$, 且$\nabla f^{(}{x}^{\ast })=0,{\nabla }^{2}f(x^{\ast })$正定, 并且搜索方向满足$$\underset{k\to \infty }{\lim }\frac{\Vert \nabla f_k+{\nabla }^2{f}_k{p}_k\Vert }{\Vert p_k\Vert }=0,$$则

1. 步长${\alpha }_k=1$对充分大的$k$($\exists k_0:k>k_0$)满足Wolfe条件;
2. 若${\alpha }_k=1,\forall k>k_0$, 则 { x k } \{x_k\} {xk​}超线性收敛于$x^{\ast }$.

易得若$c_1>1/2$, 则线搜索可能会排除二次函数的最小点, 且算法将不会采纳单位步长. 若$p_k$是拟牛顿搜索方向, 则上面的极限等价于$$\underset{k\to \infty }{\lim }\frac{\Vert (B_k-{\nabla }^{2}f(x^{\ast }))p_k\Vert }{\Vert p_k\Vert }=0.$$从这个结论中我们知道: 想要达成超线性收敛速度并不需要拟牛顿矩阵序列 { B k } \{B_k\} {Bk​}收敛于Hessian矩阵${\nabla }^{2}f(x^{\ast })$, 而只需要$B_k$沿着搜索方向$p_k$越来越靠近${\nabla }^{2}f(x^{\ast })$.

定理6 设$f:{\mathbb{R}}^n\to \mathbb{R}$二阶连续可微. 考虑迭代式$$x_{k+1}=x_k+p_k$$(即${\alpha }_k$一致选为1), $p_k=-B_k^{-1}\nabla f_k$. 设 { x k } \{x_k\} {xk​}收敛于$x^{\ast }$, $f$在$x^{\ast }$上满足二阶充分条件. 若$$\underset{k\to \infty }{\lim }\frac{\Vert (B_k-{\nabla }^{2}f(x^{\ast }))p_k\Vert }{\Vert p_k\Vert }=0,$$则有 { x k } \{x_k\} {xk​}超线性收敛于$x^{\ast }$.

证明: 我们用$p_k^N$表示牛顿步. 下面证明极限式等价于$$p_k-p_k^N=o(\Vert p_k\Vert ).$$设极限式成立, 于是$$\begin{array}{rl}{p}_k-p_k^N& ={\nabla }^2{f}_k^{-1}({\nabla }^2{f}_k{p}_k+\nabla f_k)\\ & ={\nabla }^2{f}_k^{-1}({\nabla }^2{f}_k-B_k)p_k\\ & =O(\Vert ({\nabla }^2{f}_k-B_k)p_k\Vert )(\because {\nabla }^{2}f(x)连续且在x^{\ast }处正定)\\ & =o(\Vert p_k\Vert ).(\because 极限式)\end{array}$$反向的证明是显然的. 结合牛顿步的二次收敛与上式, 有$$\Vert x_k+p_k-x^{\ast }\Vert \le \Vert x_k+p_k^N-x^{\ast }\Vert +\Vert p_k-p_k^N\Vert =O(\Vert x_k-x^{\ast }{\Vert }^2)+o(\Vert p_k\Vert ).$$由于$$\Vert x_k-x^{\ast }\Vert -\Vert p_k\Vert \le \Vert x_k+p_k-x^{\ast }\Vert ,$$$$\Vert p_k\Vert -\Vert x_k-x^{\ast }\Vert \le \Vert x_k+p_k-x^{\ast }\Vert ,$$经移项就有$$\Vert p_k\Vert +o(\Vert p_k\Vert )\ge \Vert x_k-x^{\ast }\Vert -O(\Vert x_k-x^{\ast }{\Vert }^2),$$$$\Vert p_k\Vert -o(\Vert p_k\Vert )\le O(\Vert x_k-x^{\ast }{\Vert }^2)+\Vert x_k-x^{\ast }\Vert .$$因此$\Vert p_k\Vert =O(\Vert x_k-x^{\ast }\Vert )$. 于是$$\Vert x_k+p_k-x^{\ast }\Vert \le O(\Vert x_k-x^{\ast }{\Vert }^2)+o(\Vert p_k\Vert )=o(\Vert x_k-x^{\ast }\Vert ),$$即有超线性收敛.

注意上面的结论是局部收敛性质. 我们将在后面看到, 拟牛顿法应用时极限式一般是满足的.

4. 步长选取算法

在本章的第一节, 我们谈到了选取步长的两个阶段: 确定区间、二分或插值. 我们将在本节讨论具体的步长选取算法. 考虑最小化函数$$\varphi (\alpha )=f(x_k+\alpha p_k),$$其中$\alpha >0$, $p_k$为下降方向. 若原函数$f(x)=\frac{1}{2}{x}^{T}Qx-b^{T}x$为凸二次函数, 则最优步长${\alpha }_k$具有解析表达式$${\alpha }_k=-\frac{\nabla f_k^T{p}_k}{p_k^{T}Q{p}_k}.$$而对于一般的非线性函数, 则有必要使用迭代算法确定${\alpha }_k$. 线搜索步长的算法需要精心设计, 因为步长的选取会影响非线性优化的强健性与效率.
线搜索步长的方法可根据提供的函数信息加以分类. 譬如有只需函数值的算法. 这种算法往往并不高效. 再如提供梯度信息的算法. 这种算法就能够判定是否选取到了合适的步长. 这里我们仅讨论带导数信息的算法.

4.1 插值

我们先就基于已知$\varphi$函数值和导数值进行插值的方法讨论. 这种方法可以看做是Armijo回溯算法的加强版. 目的就是寻找$\alpha$既能够带来充分的函数下降, 也不至于太小. 注意充分下降条件又可写作$$\varphi ({\alpha }_k)\le \varphi (0)+c_1{\alpha }_k{\varphi }^{\prime }(0).$$我们判定算法的高效性在于, 尽可能少地获取导数信息, 即假定获取导数信息要比获取函数信息要昂贵的多.

假设初始${\alpha }_0$给定. 若$$\varphi ({\alpha }_0)\le \varphi (0)+c_1{\alpha }_0{\varphi }^{\prime }(0),$$则${\alpha }_0$恰好满足要求, 终止搜索. 否则, $[0,{\alpha }_0]$中就应当有合适的步长. 下面根据$\varphi (0),{\varphi }^{\prime }(0),\varphi ({\alpha }_0)$建立$\varphi$的二次近似模型${\varphi }_q(\alpha )$:$${\varphi }_q(\alpha )=\left(\frac{\varphi ({\alpha }_0)-\varphi (0)-{\alpha }_0{\varphi }^{\prime }(0)}{{\alpha }_0^2}\right){\alpha }^2+{\varphi }^{\prime }(0)\alpha +\varphi (0).$$注意这个近似模型满足${\varphi }_q(0)=\varphi (0),{\varphi }_q^{\prime }(0)={\varphi }^{\prime }(0),{\varphi }_q({\alpha }_0)=\varphi ({\alpha }_0)$. 新的${\alpha }_1$则定义为此二次模型的最小值点, 即$${\alpha }_1=-\frac{{\varphi }^{\prime }(0){\alpha }_0^2}{2[\varphi ({\alpha }_0)-\varphi (0)-{\varphi }^{\prime }(0){\alpha }_0]}.$$若在${\alpha }_1$满足充分下降条件, 则终止搜索. 否则继续以$\varphi (0),{\varphi }^{\prime }(0),\varphi ({\alpha }_0),\varphi ({\alpha }_1)$建立$\varphi$的三次近似模型 (三次近似模型对捕捉函数曲率变化有较好的效果, 常带来二次收敛效果)$${\varphi }_c(\alpha )=a{\alpha }^3+b{\alpha }^2+\alpha {\varphi }^{\prime }(0)+\varphi (0),$$其中$$\left(\begin{array}{c}a\\ b\end{array}\right)=\frac{1}{{\alpha }_0^2{\alpha }_1^2({\alpha }_1-{\alpha }_0)}\left(\begin{array}{cc}{\alpha }_0^2& -{\alpha }_1^2\\ -{\alpha }_0^3& {\alpha }_1^3\end{array}\right)\left(\begin{array}{c}\varphi ({\alpha }_1)-\varphi (0)-{\varphi }^{\prime }(0){\alpha }_1\\ \varphi ({\alpha }_0)-\varphi (0)-{\varphi }^{\prime }(0){\alpha }_0\end{array}\right).$$对${\varphi }_c(x)$求导, 可得${\varphi }_c$的最小值点${\alpha }_2$位于区间$[0,{\alpha }_1]$中, 且为$${\alpha }_2=\frac{-b+\sqrt{b^2-3a{\varphi }^{\prime }(0)}}{3a}.$$若需要, 此重复这一步骤, 以$\varphi (0),{\varphi }^{\prime }(0)$和最近的两个$\varphi$函数值建立三次近似模型直至找到$\alpha$满足充分下降条件. 若算得的${\alpha }_i$要比前一个${\alpha }_{i-1}$太近或者小很多, 则置${\alpha }_i={\alpha }_{i-1}/2$以保证迭代的收敛速度以及最终的$\alpha$不会太小.

事实上, 后面会提到使用差分的方法获得导数, 从而降低不少耗费. 相应地, 我们有以$\varphi ,{\varphi }^{\prime }$在最近的两个$\alpha$值上构建某个区间$[\bar{a},\bar{b}]$上的的三次近似模型. 此时三次模型的最小值点要么在区间端点, 要么在区间内部. 内部最小值点的解析表达为$${\alpha }_{i+1}={\alpha }_i-({\alpha }_i-{\alpha }_{i-1})\left[\frac{{\varphi }^{\prime }({\alpha }_i)+d_2-d_1}{{\varphi }^{\prime }({\alpha }_i)-{\varphi }^{\prime }({\alpha }_{i-1})+2d_2}\right],$$其中$$d_1={\varphi }^{\prime }({\alpha }_{i-1})+{\varphi }^{\prime }({\alpha }_i)-3\frac{\varphi ({\alpha }_{i-1})-\varphi ({\alpha }_i)}{{\alpha }_{i-1}-{\alpha }_i},$$$$d_2=\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}({\alpha }_i-{\alpha }_{i-1})[d_1^2-{\varphi }^{\prime }({\alpha }_{i-1}){\varphi }^{\prime }({\alpha }_i){]}^{1/2}.$$对于下一次计算时改舍弃${\alpha }_{i-1}$还是${\alpha }_i$, 需要根据具体的终止条件确定. 我们将在后文对Wolfe条件进行说明.

4.2 初始步长

对于牛顿法与拟牛顿法, 显然应当以${\alpha }_0=1$为初始步长. 随着迭代点靠近满足二阶充分条件的极小点, $\alpha =1$将满足Wolfe条件. 根据前述定理, 算法将达到超线性收敛甚至二次收敛速度.
而对于那些会产生尺度不一的搜索方向的算法, 比如最速下降法与共轭梯度法, 我们就需要利用问题本身的信息来决定初始步长. 一个常用的策略是, 假设在迭代点$x_k$的一阶改变量与上一步的相同, 即$${\alpha }_0\nabla f_k^T{p}_k={\alpha }_{k-1}\nabla f_{k-1}^T{p}_{k-1},$$从而$${\alpha }_0={\alpha }_{k-1}\frac{\nabla f_{k-1}^T{p}_{k-1}}{\nabla f_k^T{p}_k}.$$另一种实用的策略是, 以$f(x_{k-1}),f(x_k),\nabla f_{k-1}^T{p}_{k-1}$构建二次模型, 并定义${\alpha }_0$为其最小值点. 这种方案给出$${\alpha }_0=\frac{2(f_k-f_{k-1})}{{\varphi }^{\prime }(0)}.$$若进一步设${\alpha }_0\leftarrow \min (1,1.01{\alpha }_0),$ 则我们会发现最终${\alpha }_0=1$总会被采纳, 从而保证了牛顿法和拟牛顿法的超线性收敛.

4.3 Wolfe条件下的线搜索步长算法

算法4 (线搜索步长算法)
设${\alpha }_0\leftarrow 0$, 选取${\alpha }_{\max }>0,{\alpha }_1\in (0,{\alpha }_{\max })$;
$i\leftarrow 1$;
repeat
$获取\varphi ({\alpha }_i)$;
\quad\quad if $\varphi ({\alpha }_i)>\varphi (0)+c_1{\alpha }_i{\varphi }^{\prime }(0)$ or $[\varphi ({\alpha }_i)\ge \varphi ({\alpha }_{i-1})$ and $i>1]$
${\alpha }_{\ast }\leftarrow$zoom(${\alpha }_{i-1},{\alpha }_i)$ and stop;
$获取{\varphi }^{\prime }({\alpha }_i)$;
\quad\quad if $|{\varphi }^{\prime }({\alpha }_i)|\le -c_2{\varphi }^{\prime }(0)$
${\alpha }_{\ast }\leftarrow {\alpha }_i$ and stop;
\quad\quad if ${\varphi }^{\prime }({\alpha }_i)\ge 0$
${\alpha }_{\ast }\leftarrow$zoom$({\alpha }_i,{\alpha }_{i-1})$ and stop;
$选取{\alpha }_{i+1}\in ({\alpha }_i,{\alpha }_{\max })$;
$i\leftarrow i+1$;
end (repeat)

下面明确zoom, 调用形式为zoom$({\alpha }_{lo},{\alpha }_{hi})$, 其中

1. 由${\alpha }_{lo},{\alpha }_{hi}$围住的区间中有满足强Wolfe条件的点;
2. ${\alpha }_{lo}$满足充分下降条件, 且给出了目前最低的函数值;
3. ${\alpha }_{hi}$满足${\varphi }^{\prime }({\alpha }_{lo})({\alpha }_{hi}-{\alpha }_{lo})<0$.

每一步迭代zoom都在${\alpha }_{lo},{\alpha }_{hi}$中间产生${\alpha }_j$, 之后再根据1,2,3提出某个端点.

算法5 (zoom)
repeat
$插值求得试探步{\alpha }_j$;
$获取\varphi ({\alpha }_j)$;
\quad\quad if $\varphi ({\alpha }_j)>\varphi (0)+c_1{\alpha }_j{\varphi }^{\prime }(0)$ or $\varphi ({\alpha }_j)\ge \varphi ({\alpha }_{lo})$
${\alpha }_{hi}\leftarrow {\alpha }_j$;
\quad\quad else
$获取{\varphi }^{\prime }({\alpha }_j)$;
\quad\quad\quad\quad if $|{\varphi }^{\prime }({\alpha }_j)|\le -c_2{\varphi }^{\prime }(0)$
${\alpha }_{\ast }\leftarrow {\alpha }_j$ and stop;
\quad\quad\quad\quad if $|{\varphi }^{\prime }({\alpha }_j)({\alpha }_{hi}-{\alpha }_{lo})\ge 0$
${\alpha }_{hi}\leftarrow {\alpha }_{lo}$;
${\alpha }_{lo}\leftarrow {\alpha }_j$;
end (repeat)

在实施过程中, 需要注意:

1. ${\alpha }_j$不应离区间端点太近;
2. 当接近最优时, $f(x_k),f(x_{k-1})$可能在机器精度下无法区分. 因此线搜索必须在多次 (譬如, 10) 无法得到函数下降后停止搜索. 类似地还可根据$x$的变化停止.

---

1. 前者对应于牛顿或拟牛顿法, 后者对应非线性共轭梯度法. ↩︎

2. 凸性是能够反映曲率的, 一方面可以从微分几何中得到解释, 另一方面凸性往往与函数的二阶信息有关. 设函数$f$的Hessian矩阵为$H$, 则在后面我们会提到, 正曲率方向就是满足$p^{T}Hp>0$的$p$, 负曲率方向就是满足$p^{T}Hp<0$的$p$. ↩︎

3. 图示可参考用“人话”解释不精确线搜索中的Armijo-Goldstein准则及Wolfe-Powell准则. ↩︎