无约束优化方法

本文可能存在纰漏，欢迎大家指出错误！非常感谢！

本文讲解的是无约束优化方法，主要包括梯度下降法、牛顿法、BFGS 三种。

后续会添加 LBFGS 方法的说明。

基本问题

对于无约束优化问题，变量$x\in R^n$，目标函数为：

$$ \min_x f(x) $$

泰勒级数

基于梯度的问题都会涉及泰勒级数问题。简单来说，泰勒级数就是在级数$f(x)$在点$x_0$的邻域内具有$n+1$阶导数的前提下，可以使用$n$阶泰勒级数去近似表示邻域内的函数$f(x)$，如下所示：

$$ f(x)=f(x_0)+\nabla f(x_0)(x-x_0)+\frac{\nabla^2 f(x_0)}{2!}(x-x_0)^2+\cdots+\frac{\nabla^n f(x_0)}{n!}(x-x_0)^n $$

梯度下降法

梯度下降法是一种通过迭代来求取极值的方法，它使用目标函数的一阶导数信息，每次迭代取目标函数下降最快的方向作为搜索方向，即梯度的反方向。下面首先说明为什么目标函数下降最快的方向为梯度的反方向，然后再继续说明梯度下降法。

搜索方向

假设目前已经经过了$k-1$次迭代，到达第$k$次迭代。选择自变量$x_k$处，对函数$f(x)$进行一阶泰勒展开，可以得到：

$$ f(x)=f(x_k)+\nabla f(x_k)(x-x_k)+o(x) $$

其中$o(x)$相对一阶泰勒为无穷小量，可以近似忽略，而$\nabla f(x_k)$是一个具体数值，表示函数$f(x)$在点$x_k$处的梯度值。结合上式，可以看作对于$f(x)$，是从$f(x_k)$出发，$x-x_k$是下一步的方向，记作$\alpha g$，$\alpha g$具体含义可以视作$x$沿着单位向量$g\in R^n$的方向迭代，$\alpha$表示步长参数。此时可以结合前文提到的含义，对一阶泰勒展开进行变形：

$$ f(x_k)-f(x)=-\nabla f(x_k)(x-x_k)=-\alpha \nabla f(x_k)g $$

参考无约束优化的基本问题，是希望函数$f(x)$的取值尽可能小，即希望每一次迭代$f(x_k)-f(x)$尽可能大。上式中，步长参数$\alpha>0$，可以认为极大化$f(x_k)-f(x)$等价于：

$$ \min \nabla f(x_k)\cdot g $$

其中$\nabla f(x_k)\cdot g$表示两个向量的乘积，显然，两个向量方向相反时便能取得最小值，如果希望迭代后的目标函数能够以最快速度下降，迭代的方向$g$应该是$\nabla f(x_k)$的反方向，也即梯度的反方向，此时$x$沿着$-\nabla f(x)$的方向前进便能得到极小值。

综上所述，梯度下降法的算法如下：

$$ x^{k+1}\coloneqq x^k-\alpha\nabla f(x_k) $$

牛顿法

牛顿法不同于梯度下降法，在于牛顿法额外引入了二阶导数信息。假设目前已经经过了$k-1$次迭代，到达第$k$次迭代，应用二阶泰勒展开对目标函数在自变量$x_k$处进行展开：

$$ f(x)=f(x_k)+\nabla f(x_k)(x-x_k)+(x-x_k)^T\frac{\nabla^2 f(x_k)}{2}(x-x_k)+o(x^2) $$

同理，无穷小量$o(x^2)$相对于二阶泰勒展开可以近似忽略。当上式导数为 0 时，得到的点即可作为下一次自变量的值$x_{k+1}$。

对等式两侧同时对$x$进行求导，令导数为 0 求得的极值点即为下一次迭代的取值$x_{k+1}$：

$$ \nabla f(x)=\nabla f(x_k)+\nabla^2 f(x_k)(x-x_k)=0 $$

公式中的$\nabla f$表示梯度，$\nabla^2 f$表示 Hessian 矩阵。为了表示方便，记$\nabla f$为$g$，$\nabla^2 f$为$H$。$g$与$H$的具体形式如下所示：

$$ g=\nabla f=\begin{bmatrix}\frac{\partial f}{\partial x_1} \\ \frac{\partial f}{\partial x_2} \\ \vdots \\ \frac{\partial f}{\partial x_n}\end{bmatrix},H=\nabla^2f=\begin{bmatrix}\frac{\partial^2 f}{\partial x_1^2} & \frac{\partial^2 f}{\partial x_1\partial x_2} & \cdots & \frac{\partial^2 f}{\partial x_1\partial x_n} \\ \frac{\partial^2 f}{\partial x_2\partial x_1} & \frac{\partial^2 f}{\partial x_2^2} & \cdots & \frac{\partial^2 f}{\partial x_2\partial x_n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial^2 f}{\partial x_n\partial x_1} & \frac{\partial^2 f}{\partial x_n\partial x_2} & \cdots & \frac{\partial^2 f}{\partial x_n^2}\end{bmatrix} $$

将符号$g$与$H$带入上面的等式，得到：

$$ g_k+H_k(x-x_k)=0 $$

对等式进行变换，得到：

$$ x-x_k=-H_k^{-1}g_k $$

因此 x 在下一次迭代的取值为：

$$ x_{k+1}=x_k-H_k^{-1}g_k $$

其中，$-H_k^{-1}g_k$即为搜索方向向量，函数取值需要在此方向上下降，因此需要与梯度$g_k$反向，即$-g_kH_k^{-1}g_k<0$，显而易见 Hessian 矩阵要求是正定的。如果 Hessian 矩阵非正定，可能导致目标函数不一定下降，从而牛顿法不收敛。

对于二次型函数，牛顿法迭代一次即可到达极值点，但是对于非二次型，牛顿法缺少步长因子，因此迭代过程中可能发散或者存在不稳定的现象。故迭代时，需要在每一步迭代时给出步长参数$\lambda_k$：

$$ \lambda_k=\arg \min_\lambda f(x_k-\lambda H_k^{-1}g_k) $$

综上所述，总结下牛顿法的特点：

• 牛顿法每次计算都需要计算目标函数的 Hessian 矩阵，计算量比较大，且需要 Hessian 矩阵正定。
• 梯度下降法是一阶收敛，牛顿法是二阶收敛，由于二阶导可能违负甚至为 0（非正定），则迭代后的$x_{k+1}$ 可能不在下降方向上。
• 牛顿法是一个二次曲面去你和当前所处位置的局部曲面，而梯度下降法是用一个平面去拟合当前的局部平面，所以牛顿法选择的下降路径会更加符合真实的最优下降路径。

拟牛顿法

牛顿法部分提到，牛顿法计算量大，且需要保证 Hessian 矩阵正定。为了克服这些缺点，拟牛顿法被提出。拟牛顿法不需要计算目标函数的 Hessian 矩阵，而是构造一个近似 Hessian 矩阵的对称正定矩阵，利用近似矩阵替代 Hessian 矩阵优化目标函数。而不同的构造 Hessian 矩阵的方法决定了不同的拟牛顿方法。

拟牛顿条件

构造 Hessian 矩阵需要满足拟牛顿条件，下面进行说明。首先将目标函数$f(x)$在$x_{k+1}$处进行二阶泰勒展开：

$$ f(x)=f(x_{k+1})+\nabla f(x_{k+1})(x-x_{k+1})+(x-x_{k+1})^T\frac{\nabla^2f(x_{k+1})}{2}(x-x_{k+1}) $$

两边同时对$x_{k+1}$求导，得到：

$$ g=g_{k+1}+H_{k+1}(x-x_{k+1}) $$

令$x=x_k$，整理可得：

$$ g_{k+1}-g_k=H_{k+1}(x_{k+1}-x_k) $$

上式即为拟牛顿条件，迭代过程中对$H_{k+1}$做出约束，根据约束构造一个近似矩阵$B_{k+1}$替代 Hessian 矩阵即可。为了简便，引入记号$s_k$与$y_k$：

$$ \begin{matrix} s_k & = & x_{k+1}-x_k \\ y_k & = & g_{k+1}-g_k \end{matrix} $$

因此得到简化符号后的拟牛顿条件：

$$ y_k=B_{k+1}\cdot s_k $$

由于牛顿法中的迭代方向为$-H^{-1}$，所以令$D_{k+1}=H_{k+1}^{-1}$，拟牛顿条件还可以写为：

$$ s_k=D_{k+1}\cdot y_k $$

BFGS

BFGS 是一种拟牛顿方法，通过迭代构建近似的 Hessian 矩阵，省去了求解 Hessian 的复杂步骤，同时使用 BFGS 构造出来的近似 Hessian 矩阵一定是正定的。虽然搜索方向不一定最优，但可以保证一定朝着最优的方向前进。

首先初始化近似 Hessian 矩阵 $B_0=I$，接下来每次迭代对矩阵$B_k$进行更新即可：

$$ B_{k+1}=B_k+\Delta B,k=1,2,\cdots $$

构造近似 Hessian 矩阵的关键即为矩阵$\Delta B$的构造，令：

$$ \Delta B_k=\alpha uu^T+\beta vv^T $$

此处的向量$u$与$v$是待定向量，$\alpha$与$\beta$为待定系数，知道了这些待定参数后即可构造$\Delta B_k$，且这样构造出的矩阵可以满足对称条件。根据拟牛顿条件：

$$ \begin{array}{ccl} y_k & = & B_{k+1}s_k \\ & = & (B_k+\Delta B_k)s_k \\ & = & (B_k+\alpha uu^T+\beta vv^T)s_k \\ & = & B_ks_k+(\alpha uu^T)s_k+(\beta vv^T)s_k \\ & = & B_ks_k+u(\alpha u^Ts_k)+v(\beta v^Ts_k) \end{array} $$

其中，$\alpha u^Ts_k$与$\beta v^Ts_k$均为实数，代表了$u$与$v$方向的拉伸程度，为了计算简单，做如下赋值运算：

$$ \begin{array}{ccc} \alpha u^Ts_k & = & 1 \\ \beta v^Ts_k & = & -1 \end{array} $$

带入上式可得：

$$ u-v=y_k-B_ks_k $$

可以得到对$u$与$v$的一个近似：

$$ \begin{array}{ccc} u & = & y_k \\ v & = & B_ks_k \end{array} $$

将上式结合近似 Hessian 矩阵$B_k$对称的特性($B_k=B_k^T$)，进而可以得到系数$\alpha$与$\beta$的值：

$$ \begin{array}{cccl} \alpha & = & \frac{1}{u^Ts_k} & =\frac{1}{y_k^Ts_k} \\ \beta & = & -\frac{1}{v^Ts_k}=-\frac{1}{(B_ks_k)^Ts_k}=-\frac{1}{s_k^TB_k^Ts_k} & =-\frac{1}{s_k^TB_ks_k} \end{array} $$

此时，$\alpha$、$\beta$、$u$与$v$均已求得，可得最终的$\Delta B_k$的更新公式：

$$ \begin{array}{ccc} \Delta B_k & = & \alpha uu^T+\beta vv^T \\ & = & \frac{y_ky_k^T}{y_k^Ts_k}-\frac{B_ks_ks_k^TB_k}{s_k^TB_ks_k} \end{array} $$

因此，迭代计算的公式为：

$$ B_{k+1}=B_k+\frac{y_ky_k^T}{y_k^Ts_k}-\frac{B_ks_ks_k^TB_k}{s_k^TB_ks_k} $$

Sherman-Morrison 公式

由于牛顿法的方向是$-H_k^{-1}g_k$的，所以最好可以直接计算出$B_k^{-1}$，这样就不需要进行求逆计算。针对 BFGS 方法，采用 Sherman-Morrison 公式进行矩阵求逆计算。下面首先对 Sherman-Morrision 公式进行说明。

Sherman-Morrision 公式：假设$A\in R^{n\times n}$为可逆矩阵，$u,v\in R^n$为列向量，则$A+uv^T$可逆当且仅当$1+v^TA^{-1}u\neq 0$，且当$A+uv^T$可逆时，逆矩阵由以下公式给出：

$$ (A+uv^T)^{-1}=A^{-1}-\frac{A^{-1}uv^TA^{-1}}{1+v^TA^{-1}u} $$

证明

($\Leftarrow$)当$1+v^TA^{-1}u\neq 0$时，令$X=A+uv^T$，$Y=A^-1-\frac{A^{-1}uv^TA^{-1}}{1+v^TA^{-1}u}$，则此时只需要证明$XY=YX=I$即可说明$A+uv^T$可逆，其中$I$为$n$阶单位矩阵。

$$ \begin{array}{ccc} XY & = & (A+uv^T)(A^{-1}-\frac{A^{-1}uv^TA^{-1}}{1+v^TA^{-1}u}) \\ & = & \end{array} $$

未完待续。。。