拟牛顿法

拟牛顿法的基本思想是用不包含二阶导数的矩阵近似牛顿法中的 Hessian 矩阵的逆矩阵，从而避免计算二阶导。拟牛顿法具有二次终止性，且对于一般情形具有 n 步二级收敛速率。缺点是所需存储量较大。是求解无约束最优化问题最有效的一类方法。

一、拟牛顿条件

  牛顿法的迭代公式为：
$$
\boldsymbol{x}^{(k+1)}=\boldsymbol{x}^{(k)}+\lambda_k\boldsymbol{d}^{(k)}
$$
其中， $\lambda_k$ 是沿着牛顿方向搜索的步长， $\boldsymbol{d}^{(k)}$ 是在点 $\boldsymbol{x}^{(k)}$ 处的牛顿方向：
$$
\boldsymbol{d}^{(k)}=-\nabla^2f(\boldsymbol{x}^{(k)})^{-1}\nabla f(\boldsymbol{x}^{(k)})
$$
  为构造 $\nabla^2f(\boldsymbol{x}^{(k)})^{-1}$ 的近似矩阵 $\bold{H}k$ ，先分析 $\nabla^2f(\boldsymbol{x}^{(k)})^{-1}$ 与一阶导数的关系。假设在第 $k$ 次迭代后得到点 $\boldsymbol{x}^{(k+1)}$ ，在该点处将目标函数进行二阶 Taylor 近似：
$$
\begin{aligned}
f(\boldsymbol{x}) \approx & f\left({\boldsymbol{x}}^{(k+1)}\right)+\nabla f\left({\boldsymbol{x}}^{(k+1)}\right)^{\mathrm{T}}\left({\boldsymbol{x}}-{\boldsymbol{x}}^{(k+1)}\right) \
&+\frac{1}{2}\left({\boldsymbol{x}}-{\boldsymbol{x}}^{(k+1)}\right)^{\mathrm{T}} \nabla^{2} f\left({\boldsymbol{x}}^{(k+1)}\right)\left({\boldsymbol{x}}-{\boldsymbol{x}}^{(k+1)}\right)
\end{aligned}
$$
因此在 $\boldsymbol{x}^{(k+1)}$ 附近有（移项然后按求导定义变换）：
$$
\nabla f(\boldsymbol{x}) \approx \nabla f\left(\boldsymbol{x}^{(k+1)}\right)+\nabla^{2} f\left(\boldsymbol{x}^{(k+1)}\right)\left(\boldsymbol{x}-\boldsymbol{x}^{(k+1)}\right)
$$
令 $\boldsymbol{x}=\boldsymbol{x}^{(k)}$ ，则
$$
\nabla f\left({\boldsymbol{x}}^{(k)}\right) \approx \nabla f\left({\boldsymbol{x}}^{(k+1)}\right)+\nabla^{2} f\left({\boldsymbol{x}}^{(k+1)}\right)\left({\boldsymbol{x}}^{(k)}-{\boldsymbol{x}}^{(k+1)}\right)
$$
记：
$$
\begin{aligned}
&\boldsymbol{p}^{(k)}={\boldsymbol{x}}^{(k+1)}-{\boldsymbol{x}}^{(k)} \
&\boldsymbol{q}^{(k)}=\nabla f\left({\boldsymbol{x}}^{(k+1)}\right)-\nabla f\left({\boldsymbol{x}}^{(k)}\right)
\end{aligned}
$$
则有：
$$
\boldsymbol{q}^{(k)} \approx \nabla^{2} f\left(\boldsymbol{x}^{(k+1)}\right) \boldsymbol{p}^{(k)}
$$
假设 Hessian 矩阵 $\nabla^{2} f\left(\boldsymbol{x}^{(k+1)}\right) $ 可逆，则有（Hessian矩阵的逆与一阶导的关系式）：
$$
\boldsymbol{p}^{(k)} \approx \nabla^{2} f\left({\boldsymbol{x}}^{(k+1)}\right)^{-1} \boldsymbol{q}^{(k)}
$$
因此当计算出 $\boldsymbol{p}^{(k)}$ 和 $\boldsymbol{q}^{(k)}$ 后，就可以根据上式估计出点 $\boldsymbol{x}^{(k+1)}$ 处的 Hessian 矩阵的逆，因此可用不包含二阶导数的矩阵 $\bold{H}{k+1}$ 替代 $\nabla^{2} f\left({\boldsymbol{x}}^{(k+1)}\right)^{-1}$ ，需要满足：
$$
\boldsymbol{p}^{(k)}=\bold{H}_{k+1}\boldsymbol{q}^{(k)}
$$
称之为拟牛顿条件。

二、秩 1 校正

原理

  当 $\nabla^{2} f\left({\boldsymbol{x}}^{(k)}\right)^{-1}$ 是 $n$ 阶对称正定矩阵时，满足拟牛顿条件的矩阵 $\bold{H}{k}$ 也应该是 $n$ 阶对称正定矩阵。这种近似矩阵的一般构造策略是：初始的 $\bold{H}{k}$ （$k=1$）取任意一个 $n$ 阶对称正定矩阵（一般取单位矩阵 $I$ ），然后通过修正 $\bold{H}{k}$ 给出 $\bold{H}{k+1}$ ，令
$$
\bold{H}{k+1}=\bold{H}{k}+\Delta\bold{H}{k}
$$
其中，$\Delta\bold{H}{k}$ 称为校正矩阵。

  确定 $\Delta\bold{H}{k}$ 的方法之一是令
$$
\Delta\bold{H}{k}=\alpha_{k} \boldsymbol{z}^{(k)}\left(\boldsymbol{z}^{(k)}\right)^{\mathrm{T}}
$$
其中，$\alpha_k$ 是一个常数，$z^{(k)}$ 是 $n$ 维列向量。这样得到的 $\Delta\bold{H}{k}$ 是一个秩为1的对称矩阵，$\boldsymbol{z}^{(k)}$ 的选择应使得拟牛顿条件得到满足，因此
$$
\boldsymbol{p}^{(k)}=\bold{H}{k} \boldsymbol{q}^{(k)}+\alpha_{k} \boldsymbol{z}^{(k)} \boldsymbol{z}^{(k) \mathrm{T}} \boldsymbol{q}^{(k)}
$$
等式两端同乘 $\boldsymbol{q}^{(k)\mathrm{T}}$ ，整理可得到
$$
\boldsymbol{q}^{(k) \mathrm{T}}\left(\boldsymbol{p}^{(k)}-\bold{H}{k} \boldsymbol{q}^{(k)}\right)=\alpha{k}\left(\boldsymbol{z}^{(k) \mathrm{T}} \boldsymbol{q}^{(k)}\right)^{2}
$$
另外
$$
\boldsymbol{z}^{(k)}=\frac{\boldsymbol{p}^{(k)}-\bold{H}{k} \boldsymbol{q}^{(k)}}{\alpha{k} \boldsymbol{z}^{(k) \mathrm{T}} \boldsymbol{q}^{(k)}}
$$
从而可得到 秩1校正公式：
$$
\bold{H}{k+1}=\bold{H}{k}+\frac{\left(\boldsymbol{p}^{(k)}-\bold{H}{k} \boldsymbol{q}^{(k)}\right)\left(\boldsymbol{p}^{(k)}-\bold{H}{k} \boldsymbol{q}^{(k)}\right)^{\mathrm{T}}}{\boldsymbol{q}^{(k) T}\left(\boldsymbol{p}^{(k)}-\bold{H}_{k} \boldsymbol{q}^{(k)}\right)}
$$
用法

在利用该公式进行极小化目标函数 $f(x)$ 时，第 $k$ 次迭代中令搜索方向为
$$
\boldsymbol{d}^{(k)}=-\bold{H}{k} \nabla f\left(\boldsymbol{x}^{(k)}\right)
$$
然后沿着该搜索方向搜索，计算步长 $\lambda_k$ ，从而确定后继点 $x^{(k+1)}$ 。进而计算后继点 $x^{(k+1)}$ 处的梯度 $\nabla f(\boldsymbol{x}^{(k+1)})$ 以及 $\boldsymbol{p}^{(k)}$ 和 $\boldsymbol{q}^{(k)}$ ，再利用 秩1校正公式 计算 $\bold{H}{k+1}$ ，接着继续计算搜索方向 $\boldsymbol{d}^{(k+1)}$。以此类推直至满足算法的终止条件（最大迭代次数或精度要求）。

优缺点

• 优点

  1. 不需要计算 Hessian 矩阵并求逆；
  2. 具有二次终止性，收敛性好。

• 缺点

  1. 需要满足一定的条件；

     仅当 $\boldsymbol{q}^{(k) \mathrm{T}}\left(\boldsymbol{p}^{(k)}-\bold{H}{k} \boldsymbol{q}^{(k)}\right)>0$ ，所求 $\bold{H}{k+1}$ 才具有正定性。

  2. 并且上一点是无法保证的，即使满足条件，也可能由于舍入误差的影响，导致 $\Delta\bold{H}_k$ 无界，数值计算困难。

三、DFP算法（变尺度法）

校正矩阵
$$
\Delta \bold{H}{k}=\frac{\boldsymbol{p}^{(k)} \boldsymbol{p}^{(k) T}}{\boldsymbol{p}^{(k) T} \boldsymbol{q}^{(k)}}-\frac{\bold{H}{k} \boldsymbol{q}^{(k)} \boldsymbol{q}^{(k) \mathrm{T}} \bold{H}{k}}{\boldsymbol{q}^{(k) T} \bold{H}{k} \boldsymbol{q}^{(k)}}
$$
DFP公式

这样定义的校正矩阵能够使得 $\bold{H}{k+1}$ 满足拟牛顿条件，记为 DFP 公式，表示为
$$
\bold{H}{k+1}=\bold{H}{k}+\frac{\boldsymbol{p}^{(k)} \boldsymbol{p}^{(k) T}}{\boldsymbol{p}^{(k) T} \boldsymbol{q}^{(k)}}-\frac{\bold{H}{k} \boldsymbol{q}^{(k)} \boldsymbol{q}^{(k) \mathrm{T}} \bold{H}{k}}{\boldsymbol{q}^{(k) T} \bold{H}{k} \boldsymbol{q}^{(k)}}
$$
其计算步骤与秩1校正公式近似：

1. 选取初始点 $x^{(1)}$ ，给定终止精度要求 $\varepsilon>0$ 。

2. 置 $\bold{H}{1}=\bold{I}{n}$ （单位矩阵），计算在点 $x^{(1)}$ 处的梯度
   $$
   \boldsymbol{g}_1=\nabla f(\boldsymbol{x}^{(1)})
   $$
   置 $k=1$ 。

3. 搜索方向 $\boldsymbol d^{(k)}=-\bold H_{k} \boldsymbol g_{k}$ 。

4. 从 $\boldsymbol x^{(k)}$ 出发，沿着 $\boldsymbol d^{(k)}$ 方向搜索，求出步长 $\lambda_k$ ，使得 $\lambda_k=\arg\min _\limits{\lambda \geq 0} f\left(\boldsymbol x^{(k)}+\lambda \boldsymbol d^{(k)}\right)$ 。然后迭代搜索点 $\boldsymbol x^{(k+1)}=\boldsymbol x^{(k)}+\lambda_k\boldsymbol d^{(k)}$ 。

5. 检验是否达到终止条件（$\left|\nabla f\left(\boldsymbol{x}^{(k+1)}\right)\right| \leqslant \varepsilon$），若达到则终止迭代；否则进入第6步。

6. 若 $k=n$ ，则令 $\boldsymbol x^{(1)}=\boldsymbol x^{(k+1)}$ ，返回步骤2；否则进入步骤7。

7. 令 $\boldsymbol g_{(k+1)}=\nabla f(\boldsymbol{x}^{(k+1)}),\boldsymbol p^{(k)}=\boldsymbol x^{(k+1)}-\boldsymbol x^{(k)}, \boldsymbol q^{(k)}=\boldsymbol g_{k+1}-\boldsymbol g_{k}$ 。利用 DFP 公式计算 $\bold{H}_{k+1}$ ，置 $k=k+1$ ，返回步骤3。

优缺点

• 用 DFP 方法构造出来的 $\bold{H}_{k}$ 均为对称正定矩阵。
• 具有二次终止性，且对于一般情形具有 n 步二级收敛速率。

四、BFGS公式及Broyden族

另一种形式的拟牛顿条件（DFP公式的对偶公式）。

参考《最优化理论与算法》P314页。

参考

• 《最优化理论与算法》陈宝林著。