数值优化方法:梯度下降，最速梯度下降与共轭梯度下降

数值优化方法:梯度下降，最速梯度下降与共轭梯度下降

本文对常见的数值优化方法进行总结，重点关注优化方法的原理、迭代过程和计算复杂度。常见的数值优化方法包括梯度下降法（最速梯度下降法）、牛顿法、共轭梯度下降法等。本笔记的算法思路只是便于对各种优化方法进行理解，并不是完整的逻辑推导，如果想了解其中的数学推导，建议还是查看相关教材。（文章里公式渲染有点慢，需要等一会）

本文相关代码在https://github.com/flztiii/numerical_optimization_test.git

梯度下降算法

算法思路

考虑到以下问题,找到一个$x^{\ast }$使$f(x^{\ast })<f(x),\forall x\in R^n$。梯度下降算法的思路很简单：每一次迭代的点要比上一次次迭代点的值更小。用公式表达就是

$f(x_{k+1})$ < $f(x_k)(1)$

接下来就是如何根据这个思路得到从$x_k$到$x_{k+1}$的增量$\Delta x$。可以对$f(x)$进行一阶泰勒展开，表达式如下所示

$f(x+\Delta x)=f(x)+f^{{}^{\prime }}(x)\cdot \Delta x(2)$

再根据公式(2)，可以得到

$f(x+\Delta x)-f(x)=f^{{}^{\prime }}(x)\cdot \Delta x$ < $0(3)$

那么，只要满足公式(4)的条件，就可以得到$f(x+\Delta x)$ < $f(x)$。为了使公式(4)成立，可以令(总所周知，平方大于0)

$\Delta x=-\alpha \cdot f^{{}^{\prime }}(x),\alpha$ > $0(4)$

所以，可以得到从$x_k$到$x_{k+1}$的迭代过程满足

$x_{k+1}=x_k+\Delta x=x_k-\alpha \cdot f^{{}^{\prime }}(x_k)(5)$

时，始终存在$f(x_{k+1})<f(x_k)$。也就是一开始谈到的思路，只要不断寻找值更小的点，就可以最终找到最小值。

算法过程

1. 给出初始点$x_0$，学习率$\alpha$
2. 计算$x_0$处的梯度$f^{{}^{\prime }}(x_0)$
3. 判断$|f^{{}^{\prime }}(x_0)|$是否小于阈值$\delta$，如果小于，停止迭代，算法结束
4. 得到下一个迭代点$x_1=x_0-\alpha \cdot f^{{}^{\prime }}(x_0)$
5. 重复2-4过程

算法总结

梯度下降算法是最为基础的优化方法之一，它只需要知道梯度方法，不需要计算海森矩阵。计算复杂度较低。但是此方法容易陷入局部极小值，收敛结果对初始点和学习率的要求较高。

最速梯度下降

算法思路

最速梯度下降算法与梯度下降算法思路相似，都是希望每一次迭代的点要比上一次次迭代点的值更小。但是最速梯度下降算法与梯度下降算法的不同之处在于，它对学习率进行了调整，使学习率并不是一个固定值，而是不断变化的。

最速梯度下降算法希望在$x_k$处沿梯度方法$f(x_k)$下降时，找到一个学习率${\alpha }_k$，使得沿$f(x_k)$下降量比选择其他$\alpha$都要大。用公式表达就是

${\alpha }_k=argmi{n}_{\alpha }f(x_k-\alpha \cdot f^{{}^{\prime }}(x_k))(6)$

以上就是最速梯度下降算法的核心。接下来，使用二次型为例子，继续进行说明。（因为二次型具有比较好的性质）我们假设

$f(x)=\frac{1}{2}{x}^{T}Ax-b^{T}x(7)$

则可以计算得到f(x)的梯度和海森矩阵分别为

$f^{{}^{\prime }}(x)=Ax-b(8)$
$f^{{}^{\prime \prime }}(x)=A(9)$

同时我们可以令

$g(\alpha )=f(x_{k+1})=f(x_k-\alpha f^{{}^{\prime }}(x_k))(10)$

则公式(6)的必要条件为$\frac{dg(\alpha )}{d\alpha }=0$。因此，我们通过公式(10)中可以得到

$g^{{}^{\prime }}(\alpha )=-f^{{}^{\prime }}(x_{k+1}{)}^T{f}^{{}^{\prime }}(x_k)=0(11)$

将公式(8)代入公式(11)可以得到

$f^{{}^{\prime }}(x_{k+1}{)}^T{f}^{{}^{\prime }}(x_k)=(A{x}_{k+1}-b{)}^T{f}^{{}^{\prime }}(x_k)$

$f^{{}^{\prime }}(x_{k+1}{)}^T{f}^{{}^{\prime }}(x_k)=(A(x_k-\alpha f^{{}^{\prime }}(x_k))-b{)}^T{f}^{{}^{\prime }}(x_k)$

$f^{{}^{\prime }}(x_{k+1}{)}^T{f}^{{}^{\prime }}(x_k)=(A{x}_k-b{)}^T{f}^{{}^{\prime }}(x_k)-\alpha f^{{}^{\prime }}(x_k{)}^{T}A{f}^{{}^{\prime }}(x_k)$

$f^{{}^{\prime }}(x_{k+1}{)}^T{f}^{{}^{\prime }}(x_k)=f^{{}^{\prime }}(x_k{)}^T{f}^{{}^{\prime }}(x_k)-\alpha f^{{}^{\prime }}(x_k{)}^{T}A{f}^{{}^{\prime }}(x_k)=0$

$\alpha =\frac{f^{{}^{\prime }}(x_k{)}^T{f}^{{}^{\prime }}(x_k)}{f^{{}^{\prime }}(x_k{)}^{T}A{f}^{{}^{\prime }}(x_k)}(12)$

所以当待优化函数为二次型时，每一次的学习率的计算公式为公式(12)。梯度下降方法就变成最速梯度下降算法。

算法过程

1. 给出初始点$x_0$
2. 计算$x_0$处的梯度$f^{{}^{\prime }}(x_0)$
3. 判断$|f^{{}^{\prime }}(x_0)|$是否小于阈值$\delta$，如果小于，停止迭代，算法结束
4. 计算当前梯度f^{’}(x_0)下的学习率，计算公式为${\alpha }_0=argmi{n}_{\alpha }f(x_0-\alpha \cdot f^{{}^{\prime }}(x_0))$(如果$f(x)$为二次型，则计算公式为${\alpha }_0=\frac{f^{{}^{\prime }}(x_0{)}^T{f}^{{}^{\prime }}(x_0)}{f^{{}^{\prime }}(x_0{)}^{T}A{f}^{{}^{\prime }}(x_0)}$)(也可以使用Armijo-Goldstein准则或Wolfe-Powell准则)
5. 得到下一个迭代点$x_1=x_0-{\alpha }_0\cdot f^{{}^{\prime }}(x_0)$
6. 重复2-5过程

算法总结

最速梯度下降算法是在梯度下降算法基础上的改进，通过修改学习率，可以提高算法的收敛速度。但是仍然无法避免梯度下降容易被困于局部极小值的问题。

牛顿法

算法思路

牛顿法的思想与梯度下降不同，它并不是寻找比当前点具有更小值的点，而是从极小值的必要条件出发。我们知道，极小值的必要条件为梯度为0，也就是

$f^{{}^{\prime }}(x)=0(13)$

从这一点出发，我们可以对$f(x)$进行二阶泰勒展开，如下所示

$f(x+\Delta x)=f(x)+f^{{}^{\prime }}(x)\Delta x+\Delta x^T{f}^{{}^{\prime \prime }}(x)\Delta x(14)$

如果满足公式(13)，则可以得到$f(x+\Delta x)=f(x)$。将其代入公式(14)可以得到

$f^{{}^{\prime }}(x)\Delta x+\Delta x^T{f}^{{}^{\prime \prime }}(x)\Delta x=0$

$\Delta x=-\frac{f^{{}^{\prime }}(x)}{f^{{}^{\prime \prime }}(x)}(15)$

这样可以得到每一次迭代的更新量。随着迭代的不断进行，最终可以找到待优化函数的极小值。

算法过程

1. 给出初始点$x_0$
2. 计算$x_0$处的梯度$f^{{}^{\prime }}(x_0)$和海森矩阵$f^{{}^{\prime \prime }}(x_0)$
3. 判断$|f^{{}^{\prime }}(x_0)|$是否小于阈值$\delta$，如果小于，停止迭代，算法结束
4. 得到下一个迭代点$x_1=x_0-\frac{f^{{}^{\prime }}(x_0)}{f^{{}^{\prime \prime }}(x_0)}$
5. 重复2-4过程

算法总结

相比于梯度下降方法，牛顿法收敛速度更快。但是，它需要计算海森矩阵，计算复杂度较高。如果待优化函数为最小二乘，则可以利用高斯-牛顿法，避免计算海森矩阵。

共轭梯度下降算法

算法思路

共轭梯度下降算法也是一种梯度下降。但与普通梯度下降、最速梯度下降不同之处在于，它下降的方向并不是梯度方向。其思路涉及到线性代数相关知识，比普通梯度下降要更加复杂。为了简化后续的说明，我们还是以二次型作为例子进行讲解，即

$f(x)=\frac{1}{2}{x}^{T}Ax-b^{T}x(16)$

在讲解之前，需要先说明一下什么是共轭。设$A$为$n$阶实对称正定矩阵，如果有两个$n$维向量$d_1$和$d_2$满足

$d_1^{T}A{d}_2=0$

则称向量$d_1$和$d_2$对于矩阵$A$共轭。若一组非零向量$d_0,d_1,...,d_{n-1}$满足

$d_i^{T}A{d}_j=0,i\ne j$

则称向量系$d_i,(i=0,1,...,n-1)$为关于矩阵$A$共轭。了解了共轭向量的概念后，可以很容易知道一组共轭向量$d_i,(i=0,1,...,n-1)$是线性无关的，证明如下。我们构造一组参数$a_i,(i=0,1,...,n-1)$，计算以下表达式成立时的$a_i$。只要$a_i$全部为0，则$d_i,(i=0,1,...,n-1)$是线性无关的。

$a_0{d}_0+a_1{d}_1+...+a_{n-1}{d}_{n-1}=0$

我们让等式两边左侧同时乘上$d_i^{T}A,\forall i$，根据共轭的条件可以得到

$a_i{d}_i^{T}A{d}_i=0,\forall i$

由于$A$为正定矩阵，$d_i^{T}A{d}_i>0$必然成立，因此，要让等式成立，必须$a_i=0,\forall i$。则可以证明$d_i,(i=0,1,...,n-1)$是线性无关的。

再回到待优化函数$f(x)$上。其中，$x$是$n$维向量，且公式(16)中的矩阵$A$为正定矩阵。我们可以假设，待优化函数$f(x)$取得最小值时的$x$为$x^{\ast }$，并且存在一组$d_i,(i=0,1,...,n-1)$关于矩阵$A$共轭。由之前的描述可以得到$d_i$是线性无关的。则可以使用$d_i$对$x^{\ast }$进行表达，表达式为

$x^{\ast }-x_0={\alpha }_0{d}_0+{\alpha }_1{d}_1+...+{\alpha }_{n-1}{d}_{n-1}(17)$

那么，可以这样来看，从初始点$x_0$开始，不断在上一个点基础上增加${\alpha }_i{d}_i$，经过$n$次迭代，就可以得到最终解$x^{\ast }$，用公式表达就是

$x_{k+1}=x_k+{\alpha }_k{d}_k(18)$

可以看到，与梯度下降的迭代公式一致。其中${\alpha }_k$就是第$k$次迭代的学习率，而$d_k$是下降方向。那么，下一步就是找到一组关于矩阵$A$共轭的方向向量$d_k$。经过推导可以认为

$d_{k+1}=-g_{k+1}+{\beta }_k{d}_k(19)$

$g_{k+1}=f^{{}^{\prime }}(x_{k+1})$

接下来只要计算出${\beta }_k$即可。在公式(19)左侧同时乘上$d_k^{T}A$可以得到

$d_k^{T}A{d}_{k+1}=-d_k^{T}A{g}_{k+1}+{\beta }_k{d}_k^{T}A{d}_k=0$

${\beta }_k=\frac{-d_k^{T}A{g}_{k+1}}{d_k^{T}A{d}_k}(20)$

下降方向已经得到了，迭代过程中唯一未知的量还有学习率${\alpha }_k$。学习率的计算方法与公式(12)相同，最后计算得到的学习率为

${\alpha }_k=-\frac{g_k^T{d}_k}{d_k^{T}A{d}_k}$

当然，以上公式只有在$f(x)$是二次型的时候才成立，路过不是二次型，需要借助其他方法进行求解。

算法过程

如果$f(x)$是二次型

1. 设定$k:=0$，选择优化初始点$x_0$
2. 计算$x_0$处的梯度$g_0=f^{{}^{\prime }}(x_0)$，判断$g_0$是否为0，如果是，结束迭代；反之，令$d_0=-g_0$
3. 计算学习率${\alpha }_k=-\frac{g_k^T{d}_k}{d_k^{T}A{d}_k}$
4. 得到新的点$x_{k+1}=x_k+{\alpha }_k{d}_k$
5. 计算梯度$g_{k+1}=f^{{}^{\prime }}(x_{k+1})$，判断$g_{k+1}$是否为0，如果是，结束迭代
6. 计算新的下降方向参数${\beta }_k=\frac{g_{k+1}^{T}A{d}_k}{d_k^{T}A{d}_k}$
7. 计算新的梯度下降方向$d_{k+1}=-g_{k+1}+{\beta }_k{d}_k$
8. 重复3-7过程

如果$f(x)$不是二次型

1. 设定$k:=0$，选择优化初始点$x_0$
2. 计算$x_0$处的梯度$g_0=f^{{}^{\prime }}(x_0)$，判断$g_0$是否为0，如果是，结束迭代；反之，令$d_0=-g_0$
3. 计算学习率${\alpha }_k=argmi{n}_{\alpha }f(x_k+\alpha \cdot d_k)$(可以使用Armijo-Goldstein准则或Wolfe-Powell准则)
4. 得到新的点$x_{k+1}=x_k+{\alpha }_k{d}_k$
5. 计算梯度$g_{k+1}=f^{{}^{\prime }}(x_{k+1})$，判断$g_{k+1}$是否为0，如果是，结束迭代
6. (a) 计算新的下降方向参数${\beta }_k=\frac{g_{k+1}^T{g}_{k+1}}{g_k^T{g}_k}$(Fletcher-Reeves Formula) (b) 计算新的下降方向参数${\beta }_k=\frac{g_{k+1}^T(g_{k+1}-g_k)}{g_k^T{g}_k}$(Polak-Ribiere Formula) © 计算新的下降方向参数${\beta }_k=\frac{g_{k+1}^T(g_{k+1}-g_k)}{d_k^T(g_{k+1}-g_k)}$(Hestenes-Stiefel Formula)
7. 计算新的梯度下降方向$d_{k+1}=-g_{k+1}+{\beta }_k{d}_k$
8. 重复3-7过程

算法总结

共轭梯度下降算法对传统梯度下降算法的下降方向进行了优化，其收敛速度介于梯度下降算法和牛顿法之间。并且只需要计算梯度方向，计算复杂度较低。