牛顿法与拟牛顿法

牛顿法和拟牛顿法是求解无约束最优化的常用方法，有收敛速度快的优点. 牛顿法属于迭代算法，每一步需要求解目标函数的海赛矩阵的逆矩阵，计算复杂. 拟牛顿法通过正定矩阵近似海赛矩阵的逆矩阵，简化了这个过程.

牛顿法

对于无约束优化
$$\underset{x\in R^n}{\min }f(x)$$
$x^{\ast }$是目标的极小值点.

假设$f(x)$有二阶连续偏导数，第k次迭代值为$x^{(}k)$，将$f(x)$在$x^{(}k)$附近进行二阶展开：
$$\begin{array}{rl}f(x)& =f(x^{(k)})+f^{\prime }(x^{(k)})\Delta x+\frac{1}{2}{f}^{\prime \prime }(x^{(k)}\Delta x^2\\ & =f(x^{(k)})+f^{\prime }(x^{(k)})(x-x^{(k)})+\frac{1}{2}{f}^{\prime \prime }(x^{(k)})(x-x^{(k)}{)}^2\end{array}$$
当这里的$x$是高维时
$$f(x)=f(x^{(k)})+g_k^T(x-x^{(k)})+\frac{1}{2}(x-x^{(k)}{)}^{T}H(x^{(k)})(x-x^{(k)})$$
其中，$g_k=g(x^{(k)})=▽f(x^{(k)})$是$f(x)$的梯度向量在点$x^{(k)}$的值，$H(x^{(k)})$是$f(x)$的海赛矩阵在点$x^{(k)}$的值.

海赛矩阵式二阶导数矩阵
$$H(x)=[\frac{{\partial }^{2}f}{\partial x_i\partial x_j}{]}_{n\times n}$$

$f(x)$有极值的必要条件是在极值点处一阶导数为0，即梯度向量为0. 当$H(x^{(k)})$为正定矩阵时，$f(x)$的极值为极小值.

假设$▽f(x^{(k+1)})=0$，我们有
$$▽f(x)=f^{\prime }(x)+f^{\prime \prime }(x)(x-x^{(k)})$$
$$▽f(x)=g_k+H_k(x-x^{(k)})\text{\hspace{1em}(1.1)}$$
所以
$$\begin{gathered} g_k+H_k(x-x^{(k)})=0 \\ {x}^{(k+1)}=x^{(k)}-H_k^{-1}{g}_k \end{gathered}$$
我们假设
$$H_k{p}_k=-g_k$$
则有
$$x^{(k+1)}=x^{(k)}+p_k$$

算法过程

输入：目标函数$f(x)$，梯度$g(x)=▽f(x)$，海赛矩阵$H(x)$，精度要求$ϵ$.

输出：$f(x)$的极小值点$x^{\ast }$.

1. 取初始点$x^{(0)}$，令$k=0$.
2. 计算$g_k$.
3. 若$||g_k||<ϵ$，则停止计算，得近似解$x^{\ast }=x^{(k)}$.
4. 计算$H_k$，并求$p_k$
5. $x^{(k+1)}=x^{(k)}+p_k$
6. $k=k+1$，转2

拟牛顿法

基本思想：

海赛矩阵的逆矩阵计算困难，考虑用一个n阶矩阵$G_k=G(x^{(k)})$来近似替代$H_k^{-1}$.

在$(1.1)$中代入$x=x^{(k+1)}$，即有如下：
$$\begin{gathered} ▽f(x^{(k+1)})=g_k+H_k(x^{(k+1)}-x^{(k)}) \\ {g}_{k+1}-g_k=H_k(x^{(k+1)}-x^{(k)}) \end{gathered}$$

记$y_k=g_{(k+1)}-g_k,{\delta }_k=x^{(x+1)}-x^{(k)}$，则有
$$\begin{gathered} y_k=H)k{\delta }_k \\ {H}_k^{-1}{y}_k={\delta }_k\text{\hspace{1em}(2.1)} \end{gathered}$$
我们称$(2.1)$为拟牛顿条件.

如果$H_k$是正定的，那么可以保证牛顿法搜索方向$p_k$是下降方向：

因为搜索方向是$p_k=-\lambda g_k$
$$x=x^{(k)}+\lambda p_k=x^{(k)}-\lambda H_k^{(-1)}{g}_k$$

所以$f(x)$在$x^{(k)}$的泰勒展开式可以近似写成：
$$f(x)=f(x^{(k)})-\lambda g_k^T{H}_k^{-1}{g}_k$$
由于$H_k^{-1}$正定，所以有$g_k^T{H}_k^{-1}{g}_k>0$. 当$\lambda$为一个充分小的正数时，总有$f(x)<f(x^{(k)})$，也就是说$p_k$是下降方向.

拟牛顿法将$G_k$作为$H_k^{-1}$的近似，要求矩阵$G_k$满足同样的条件，每次迭代矩阵$G_k$都是正定的，且$G_k$要满足拟牛顿条件：
$$G_{k_1}{y}_k={\delta }_k$$

按照拟牛顿条件选择$G_k$作为$H_k^{-1}$的近似或选择$B_k$作为$H_k$的近似的算法称为拟牛顿法.

按照拟牛顿条件，每次迭代都可以更新矩阵：
$$G_{k+1}=G_k+\Delta G_k$$

有多种具体实现方法

DFP(Davidon-Fletcher-Powell)算法

DFP选择$G_{k+1}$的方法是，假设每一次迭代$G_{k+1}$都是由$G_k$加上两个附加项构成，我们假设这两个附加项分别为$Q_k、P_k$，则有
$$\begin{gathered} G_{k+1}=G_k+Q_k+P_k \\ {G}_{k+1}{y}_k=G_k{y}_k+Q_k{y}_k+P_k{y}_k \end{gathered}$$

使$G_{k+1}$满足拟牛顿条件，可以使
$$\begin{gathered} P_k{y}_k={\delta }_k \\ {Q}_k{y}_k=-G_k{y}_k \end{gathered}$$

取
$$\begin{gathered} P_k=\frac{{\delta }_k{\delta }_k^T}{{\delta }^T{y}_k} \\ Q=-\frac{G_k{y}_k{y}_k^T{G}_k}{y_k^T{G}_k{y}_k} \end{gathered}$$

所以可以得到$G_{k+1}$的迭代公式
$$G_{k+1}=G_k-\frac{G_k{y}_k{y}_k^T{G}_k}{y_k^T{G}_k{y}_k}+\frac{{\delta }_k{\delta }_k^T}{{\delta }^T{y}_k}\text{\hspace{1em}(2.2)}$$

算法过程：

输入：$f(x)$，梯度$g(x)=▽f(x)$，精度$ϵ$.

输出：$f(x)$的极小点$x^{\ast }$

1. 选取初始点$x^{(0)}$，取$G_0$为正定对称矩阵，令$k=0$

2. 计算$g_k$，若$||g_k||<ϵ$，则停止计算，得到近似解$x^{\ast }=x^{(k)}$

3. 令$p_k=-G_k{g}_k$

4. 一维搜索：

   求${\lambda }_k$使得$f(x^{(k)}+{\lambda }_k{p}_k)={\min }_{\lambda \ge 0}f(x^{(k)}+\lambda p_k)$

5. 令$x^{(k+1)}=x^{(k)}+{\lambda }_k{p}_k$

6. 计算$g_{k+1}$，若$||g_{k+1}||<ϵ$，则停止计算，得近似解$x^{\ast }=x^{(k)}$，否则按式$(2.2)$计算$G_{k+1}$.

7. $k=k+1$

BFGS(Boyden-Fletcher-Goldfarb-Shanno)算法

BFGS是最流行得拟牛顿算法.

基本思想：
考虑用$B_k$逼近海赛矩阵$H$，则拟牛顿条件为：

$$B_{k+1}{\delta }_k=y_k$$
令
$$\begin{gathered} B_{k+1}=B_k+P_k+Q_k \\ {B}_{k+1}{\delta }_k=B_k{\delta }_k+P_k{\delta }_k+Q_k{\delta }_k \end{gathered}$$

使$P_k$和$Q_k$满足：
$$\begin{gathered} P_k{\delta }_k=y_k \\ {Q}_k{\delta }_k=-B_k{\delta }_k \end{gathered}$$
取
$$\begin{gathered} P_k=\frac{y_k{y}_k^T}{y_k^T{\delta }_k} \\ {Q}_k=-\frac{B_k{\delta }_k{\delta }_k^T{B}_k}{{\delta }_k^T{B}_k{\delta }_k} \end{gathered}$$
我们得到了$B_{k+1}$的迭代公式：
$$B_{k+1}=B_k-\frac{B_k{\delta }_k{\delta }_k^T{B}_k}{{\delta }_k^T{B}_k{\delta }_k}+\frac{y_k{y}_k^T}{y_k^T{\delta }_k}\text{\hspace{1em}(2.3)}$$
可以证明，如果初始矩阵$B_0$是正定的，那么迭代过程中的每个矩阵$B_k$都是正定的.

算法过程：

输入：$f(x)$，梯度$g(x)=▽f(x)$，精度$ϵ$.

输出：$f(x)$的极小点$x^{\ast }$

1. 选取初始点$x^{(0)}$，取$B_0$为正定对称矩阵，令$k=0$

2. 计算$g_k$，若$||g_k||<ϵ$，则停止计算，得到近似解$x^{\ast }=x^{(k)}$

3. 令$p_k=-B_k{g}_k$

4. 一维搜索：

   求${\lambda }_k$使得$f(x^{(k)}+{\lambda }_k{p}_k)={\min }_{\lambda \ge 0}f(x^{(k)}+\lambda p_k)$

5. 令$x^{(k+1)}=x^{(k)}+{\lambda }_k{p}_k$

6. 计算$g_{k+1}$，若$||g_{k+1}||<ϵ$，则停止计算，得近似解$x^{\ast }=x^{(k)}$，否则按式$(2.3)$计算$B_{k+1}$.

7. $k=k+1$

Broden类算法

可以从BFGS算法的$B_k$矩阵得到关于$G_k$的迭代公式. 若记$G_k=B_k^{-1}、G_{k+1}=B_{k+1}^{-1}$，对$(2.3)$用两次Sherman-Morrison公式可得

$$G_{k+1}=(I-\frac{{\delta }_k{y}_k^T}{{\delta }_k^T{y}_k})G_k(I-\frac{{\delta }_k{y}_k^T}{{\delta }_k^T{y}_k}{)}^T+\frac{{\delta }_k{y}_k^T}{{\delta }_k^T{y}_k}$$

参考文献

李航-统计学习方法