写在前面，本系列主要是对下面这本书做的学习笔记

常用数学符号的 LaTeX 表示方法

Markdown 常用数学符号和公式

基本概念

前言

最优化方法在机器学习领域处于中心地位，绝大多数算法最后都归结于求解最优化问题，从而确定模型参数，或直接获得预测结果

有监督学习，通过最小化损失函数或优化其他类型的目标函数确定模型的参数
数据降维算法，通过优化某种目标函数，确定降维后的结果，如主成分分析

按照优化变量的类型，可以将优化问题分为连续型优化问题与离散型优化问题

连续型优化问题的优化变量是连续变量，一般可借助导数求解

离散型优化问题的优化变量则为离散值

连续型优化问题又可分为凸优化问题与非凸优化问题，凸优化问题可以保证求得全局最优解

按照目标函数的数量

单目标优化：只有一个目标函数,
多目标优化：有多个目标函数

按照是否带有约束条件

不带约束的优化
带约束的优化

按照求解方式可分为

数值优化：求问题的近似解
解析解：求精确的公式解

基于概率的优化算法是优化算法家族中一种特殊的存在，典型的是遗传算法与贝叶斯优化

通常情况下最优化问题是求函数的极值，其优化变量是整数或实数

有一类特殊的优化问题，其目标是寻找某一函数，使得泛函的值最大化或最小化，变分法是求解此类问题的经典方法

基本概念

最优化问题的目标是求函数或泛函的极值(Extrema)，在基础数学、计算数学、应用数学以及工程、管理、经济学领域均有应用

最优化算法是求解最优化问题的方法

确定优化目标函数之后，需要根据问题的特点以及现实条件的限制选拝合适的算法

在机器学习与深度学习库中，最优化算法通常以优化器(Optimizer)或求解器(Solver)的形式出现

问题定义

接下来考虑的最优化问题是求解函数极值的问题，包括极大值与极小值

要计算极值的函数称为目标函数，其自变量称为优化变量，对于函数

$f(x)=(x-2)^{2}+5$

其极小值在$x=2$点处取得，此时函数值为$5$，$x=2$ 为该问题的解

一般将最优化问题统一表述为极小值问题，对于极大值问题，只需要将目标函数反号，即可转化为极小值问题

要求$f(x)$的极大值，等价于求$-f(x)$的极小值，最优化问题可以形式化地定义为

$\min _{\boldsymbol{x}} f(\boldsymbol{x}) \quad \boldsymbol{x} \in X$

其中$x$为优化变量，$f$为目标函数，$X \subseteq \mathbb{R}^{n}$为优化变量允许的取值集合，称为可行域(Feasible Set)，它由目标函数的定义域、等式及不等式约束(Constraint Function)共同确定

可行域之内的解称为可行解(Feasible Solution)，下面是一个典型的最优化问题

$\min _{x}\left(x^{3}-4 x^{2}+5\right) \quad x \in[-10,10]$

该问题的可行域为区间$[-10,10]$，如不进行特殊说明，这里目标函数均指多元函数，一元函数为其特例，无须单独讨论

线性规划：目标函数为一次函数(线性函数)
非线性规划：目标函数是非线性函数

非线性规划的一种特例是目标函数为二次函数，称为二次规划

很多实际应用问题可能带有等式和不等式约束条件，可以写成

$\min _{\boldsymbol{x}} f(\boldsymbol{x}) \quad g_{i}(\boldsymbol{x})=0, i=1, \cdots, p \quad h_{i}(\boldsymbol{x}) \leqslant 0, i=1, \cdots, q$

这里将不等式约束统一写成小于或等于号的形式，满足等式和不等式约束条件的解称为可行解，否则称为不可行解

下面是一个带有等式和不等式约束的优化问题

$\min _{x, y}\left(x^{2}+2 y^{2}-3 x y+4 y\right) \quad x+y=1 \quad x^{2}+y^{2} \leqslant 4$

等式和不等式约束定义的区域与目标函数定义域的交集为可行域

此问题的可行域为约束条件$x+y=1$与$x^{2}+y^{2} \leqslant 4$所定义的交集，是直线位于圆内的部分，如下图所示

xxx图可行域示意图

在很多实际问题中出现的二次规划可以写成下面的形式

$\begin{array}{l} \min _{\boldsymbol{x}}\left(\frac{1}{2} \boldsymbol{x}^{\mathrm{T}} \boldsymbol{Q} \boldsymbol{x}+\boldsymbol{c}^{\mathrm{T}} \boldsymbol{x}\right) \\ \boldsymbol{A x} \leqslant \boldsymbol{b} \end{array}$

其中$x \in \mathbb{R}^{n}$，$Q$是$n \times n$的二次项系数矩阵，$c \in \mathbb{R}^{n}$是一次项系数向量

$\boldsymbol{A}$是$m \times n$的不等式约束系数矩阵，$\boldsymbol{b} \in \mathbb{R}^{m}$是不等式约束的常数向量

局部最优解与全局最优解

假设 $x^{*}$ 是一个可行解，如果对可行域内所有点 $x$ 都有 $f\left(x^{*}\right) \leqslant f(x)$ ，则称 $x^{*}$ 为全局极小值

类似地可以定义全局极大值

全局极小值是最优化问题的解

对于可行解 $x^{*}$ ，如果存在其 $\delta$ 邻域，使得该邻域内的所有点即所有满足 $\left\|x-x^{*}\right\| \leqslant \delta$ 的点 $\boldsymbol{x}$ ，都有 $f\left(\boldsymbol{x}^{*}\right) \leqslant f(\boldsymbol{x})$ ，则称 $\boldsymbol{x}^{*}$ 为局部极小值

类似地可以定义局部极大值

局部极小值可能是最优化问题的解，也可能不是

最优化算法的目标是寻找目标函数的全局极值点而非局部极值点。下图为全局最优解与局部最优解的示意图，目标函数为

xxx图全局最优解与局部最优解

上图中的目标函数有两个局部极大值点、1个局部极小值点，区间$[0,2]$内的局部极大值点也是全局极大值点

$f(x) = (-x^4 + 10x^3 + 100x^2+10)e^{-x^2}$

迭代法的基本思想

如果目标函数可导，那么可以利用导数信息确定极值点

微积分为求解可导函数极值提供了统一的方法，即寻找函数的驻点

根据费马引理，对于一元函数局部极值点必定是导数为0的点；对于多元函数则是梯度为0的点(在数值计算中，也称为静止点(Stationary Point)

机器学习中绝大多数目标函数可导，因此这种方法是适用的

通过求解驻点来寻找极值虽然在理论上可行，但实现时却存在困难

实际问题中目标函数梯度为0的方程组通常难以求解，对于下面的二元目标函数

$f(x，y)=x^{3}-2 x^{2}+\mathrm{e}^{x y}-y^{3}+10 y^{2}+100 \sin (x y)$

对$x$和$y$分别求偏导数并令它们为0，得到如下方程组

$\left\{\begin{array}{l} 3 x^{2}-4 x+y \mathrm{e}^{x y}+100 y \cos (x y)=0 \\ x \mathrm{e}^{x y}-3 y^{2}+20 y+100 x \cos (x y)=0 \end{array}\right.$

显然，这个方程组很难求解，含有指数函数、对数函数、三角函数以及反三角函数的方程一般情况下没有公式解，称为超越方程

即使是代数方程(多项式方程)，4次以上的方程没有求根公式

方程系数的有限次加减乘除以及开方运算均不可能是方程的根

因此，直接解导数为0的方程组不是一种可行的方法

对于大多数最优化问题通常只能近似求解，称为数值优化

一般采用迭代法，从一个初始可行点$x_{0}$开始，反复使用某种规则迭代直至收敛到最优解

具体地，在第 $k$ 次迭代时，从当前点 $x_{k-1}$ 移动到下一个点 $x_{k}$

如果能构造一个数列${x_k}$，保证它收敛到梯度为0的点，即下面的极限成立

$\lim \limits_{k \rightarrow+\infty} \nabla f\left(x_{k}\right)=\mathbf{0}$

则能找到函数的极值点

算法的核心

1️⃣如何定义从上一个点移动到下一个点的规则，这些规则一般利用一阶导数(梯度)或二阶导数(黑塞矩阵)

因此，迭代法的核心是得到如下式形式的迭代公式

$\boldsymbol{x}_{k+1}=h\left(\boldsymbol{x}_{k}\right)$

梯度下降法，牛顿法及拟牛顿法均采用了此思路，区别在于构造迭代公式的方法不同，迭代法的原理如图所示

xxx图迭代法的原理

2️⃣初始点$\boldsymbol{x}_{0}$的选择，通常用常数或随机数进行初始化

算法要保证对任意可行的$x_{0}$均收敛到极值点处

一阶优化算法

一阶优化算法利用目标函数的一阶导数构造迭代公式，典型代表是梯度下降法及其变种

本节介绍基本的梯度下降法、最速下降法、梯度下降法的其他改进版本(包括动量项、 AdaGrad、RMSProp、AdaDelta、Adam 算法等)以及随机梯度下降法

梯度下降法

梯度下降法(Gradient Descent Method)由数学家柯西提出，它沿着当前点 $\boldsymbol{x}_{k}$ 处梯度相反的方向进行迭代，得到 $x_{k+1}$ ，直至收致到梯度为0的点

其理论依据：在梯度不为0的任意点处，梯度正方向是函数值上升的方向，梯度反方向是函数值下降的方向

下面先通过例子说明，然后给出严格的证明

1️⃣首先考虑一元函数的情况，如图所示

xxx图一元函数的导数值与函数单调性的关系

对于一元函数，梯度是一维的，只有两个方向：沿着$x$轴向右和向左

如果导数为正，则梯度向右；否则向左

当导数为正时，是增函数，$x$变量向右移动时(即沿着梯度方向)函数值增大，否则减小，对于上图所示的函数

当$x<x_{0}$时，导数为正，此时向左前进函数值减小，向右则增大
当$x>x_{0}$时，导数为负，此时向左前进函数值增大，向右则减小

2️⃣接下来考虑二元函数，二元函数的梯度有无穷多个方向

对于函数$x^{2}+y^{2}$，其在$(x，y)$点处的梯度值为$(2 x，2 y)$，函数在$(0,0)$点处有极小值，在任意点$(x，y)$处，从$(0,0)$点指向$(x，y)$方向(即梯度方向)的函数值都是单调递增的

该函数的形状如下图所示

xxx图4.6$x^{2}+y^{2}$的形状

下图为该函数的等高线，在同一条等高线上的所有点处函数值相等

xxx图4.6$x^{2}+y^{2}$的等高线

在任意点处，梯度均为从原点指向该点，是远离原点的方向，函数值单调增

下面考虑函数$-x^{2}-y^{2}$，其在$(x，y)$点处的梯度值为$(-2 x,-2 y)$，函数在$(0,0)$点处有极大值，在任意点$(x，y)$处，从$(x，y)$点指向$(0,0)$方向的函数值都是单调递增的，$(-x,-y)$即梯度方向，下图为该目标函数的形状

xxx图$-x^{2}-y^{2}$的形状

下面给出严格的证明，将函数在$\boldsymbol{x}$点处作一阶泰勒展开

$f(\boldsymbol{x}+\Delta \boldsymbol{x})=f(\boldsymbol{x})+(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x}+o(\|\Delta \boldsymbol{x}\|)$

对上式变形，函数的增量与自变量增量、函数梯度的关系为

$f(\boldsymbol{x}+\Delta \boldsymbol{x})-f(\boldsymbol{x})=(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x}+o(\|\Delta \boldsymbol{x}\|)$

如果令$\Delta \boldsymbol{x}=\nabla f(\boldsymbol{x})$，则有

$f(\boldsymbol{x}+\Delta \boldsymbol{x})-f(\boldsymbol{x})=(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \nabla f(\boldsymbol{x})+o(\|\Delta \boldsymbol{x}\|)$

如果$\Delta x$足够小，则可以忽略高阶无穷小项，有

$f(\boldsymbol{x}+\Delta \boldsymbol{x})-f(\boldsymbol{x}) \approx(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \nabla f(\boldsymbol{x}) \geqslant 0$

如果在$\boldsymbol{x}$点处梯度不为$\boldsymbol{0}$，则能保证移动到$\boldsymbol{x}+\Delta \boldsymbol{x}$时函数值增大，相反地，如果令$\Delta \boldsymbol{x}=-\nabla f(\boldsymbol{x})$，则有

$f(\boldsymbol{x}+\Delta \boldsymbol{x})-f(\boldsymbol{x}) \approx-(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \nabla f(\boldsymbol{x}) \leqslant 0$

即函数值减小，事实上，只要确保

$(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x} \leqslant 0$

则有

$f(\boldsymbol{x}+\Delta \boldsymbol{x}) \leqslant f(\boldsymbol{x})$

因此，选择合适的增量$\Delta x$就能保证函数值下降，负梯度方向是其中的一个特例

接下来证明：增量的模一定时，在负梯度方向，函数值是下降最快的

由于

$(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x}=\|\nabla f(\boldsymbol{x})\| \cdot\|\Delta \boldsymbol{x}\| \cdot \cos \theta$

其中$\theta$为$\nabla f(\boldsymbol{x})$与$\Delta \boldsymbol{x}$之间的夹角，因此，如果$\theta<\frac{\pi}{2}$, 则$\cos \theta>0$，从而有

$(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x} \geqslant 0$

此时函数值增大

$f(\boldsymbol{x}+\Delta \boldsymbol{x}) \geqslant f(\boldsymbol{x})$

$\Delta \boldsymbol{x}$沿着正梯度方向是其特例，如果$\theta>\frac{\pi}{2}$，则$\cos \theta<0$，从而有

$(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x} \leqslant 0$

此时函数值下降

$f(\boldsymbol{x}+\Delta \boldsymbol{x}) \leqslant f(\boldsymbol{x})$

$\Delta x$沿着负梯度方向即$\theta=\pi$是其特例，由于$-1 \leqslant \cos \theta \leqslant 1$，因此，如果向量$\Delta x$的模大小一定，则$\Delta x=-\nabla f(\boldsymbol{x})$，即在梯度相反的方向函数值下降最快，此时$\cos \theta=-1$。

梯度下降法每次的迭代增量为

$\Delta \boldsymbol{x}=-\alpha \nabla f(\boldsymbol{x})$

其中$\alpha$为人工设定的接近于0的正数，称为步长或学习率，其作用是保证$\boldsymbol{x}+\Delta \boldsymbol{x}$在$\boldsymbol{x}$的邻域内、从而可以忽略泰勒公式中的$o(|\Delta x|)$项，否则不能保证每次达代时函数值下降，使用该增量则有

$(\nabla f(\boldsymbol{x}))^{\mathrm{T}} \Delta \boldsymbol{x}=-\alpha(\nabla f(\boldsymbol{x}))^{\mathrm{T}}(\nabla f(\boldsymbol{x})) \leqslant 0$

函数值下降，由此得到梯度下降法的迭代公式，从初始点$x_{0}$开始，反复使用如下迭代公式

$\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}-\alpha \nabla f\left(\boldsymbol{x}_{k}\right)$

只要没有到达梯度为0的点，函数值会沿序列$\boldsymbol{x}_{k}$递减，最终收敛到梯度为0的点

从 $\boldsymbol{x}_{0}$ 出发，用上式进行迭代，会形成一个函数值递减的序列 $\left\{x_{i}\right\}$

$f\left(x_{0}\right) \geqslant f\left(x_{1}\right) \geqslant f\left(x_{2}\right) \geqslant \cdots \geqslant f\left(x_{k}\right)$

迭代终止的条件是函数的梯度值为$\mathbf{0}$(实际实现时是接近于 0 即可)，此时认为已经达到极值点。梯度下降法的流程如算法$4.1$所示

🧮算法梯度下降法

初始化$x_{0}$，$k=0$

while $\left|\nabla f\left(\boldsymbol{x}_{k}\right)\right|>$eps and $k<N$ do

迭代 $x_{k+1}=x_{k}-\alpha \nabla f(x_k)$

$k=k+1$

end while

$x_{0}$可初始化为固定值，如$\mathbf{0}$，或随机数(通常为均匀分布或正态分布)，后者在训练神经网络时经常被采用

eps为人工指定的接近于0的正数，用于判定梯度是否已经接近于$0$；$N$为最大迭代次数，防止死循环的出现

梯度下降法在每次迭代时只需要计算函数在当前点处的梯度值，具有计算量小、实现简单的优点

只要未到达驻点处且学习率设置恰当，每次达代时均能保证函数值下降

下图为用梯度下降法求解$x^{2}+y^{2}$极值的过程，迭代初始值(图中的大圆点)设置为$(0,4)$，学习率设置为$0.1$

每次迭代时的值$\boldsymbol{x}_{i}$以小圆点显示

学习率$\alpha$的设定也是需要考虑的问题，一般情况下设置为固定的常数，如$10^{-5}$

在深度学习中，采用了更复杂的策略，可以在迭代时动态调整其值

最速下降法

梯度下降法中步长$\alpha$是固定的，或者根据某种人工指定的策略动态调整

最速下降法(Steepest Descent Method)是对梯度下降法的改进，它用算法自动确定步长值

最速下降法同样沿着梯度相反的方向进行迭代，但每次需要计算最佳步长$\alpha$

最速下降法的搜索方向与梯度下降法相同，也是负梯度方向

$\boldsymbol{d}_{k}=-\nabla f\left(\boldsymbol{x}_{k}\right)$

在该方向上寻找使得函数值最小的步长，通过求解如下一元函数优化问题实现

$\alpha_{k}=\arg \min _{\alpha} f\left(\boldsymbol{x}_{k}+\alpha d_{k}\right)$

优化变量是$\alpha$，实现时有两种方案

将 $\alpha$ 的取值离散化，取典型值 $\alpha_{1},\cdots,\alpha_{n}$ ，分别计算取这些值的目标函数值，然后确定最优值
直接求解上式目标函数的驻点，对于有些情况可得到解析解

这类方法也称为直线搜索(Line Search)，它沿着某一确定的方向在直线上寻找最优步长

梯度下降法的改进

梯度下降法在某些情况下存在收敛速度慢、收敛效果差的同题，因此出现了大量改进方案

标准的梯度下降法可能存在振荡问题，具体表现为在优化变昰的某些分量方向上来回振荡，导致收敛速度慢

下图显示了用梯度下降法求解 $0.1 x_{1}^{2}+2 x_{2}^{2}$ 的极值时的逄代过程，可以看到，迭代序列在 $x_{2}$ 方向来回振荡

xxx图梯度下降法振荡问题

动量项梯度下降法

动量项梯度下降法通过引入动量项解决此问题，类似于物理中的动量，依靠惯性保持迭代时的前进方向

动量项的计算公式为

$\boldsymbol{v}_{k}=-\alpha \nabla f\left(\boldsymbol{x}_{k}\right)+\mu \boldsymbol{v}_{k-1}$

它是上次迭代时的动量项与本次负梯度值的加权和，其中$\alpha$是学习率，其作用与标准的梯度下降法相同，$\mu$是动量项系数

如果按照时间线展开，则第$k$次迭代时使用了从1到$k$次迭代时的所有负梯度值，且负梯度值按系数$\mu$指数级衰减，即使用了移动指数加权平均

反复利用上式，展开之后的动量项为

$\begin{aligned} v_{k}= & -\alpha \nabla f\left(\boldsymbol{x}_{k}\right)+\mu v_{k-1}=-\alpha \nabla f\left(\boldsymbol{x}_{k}\right)+\mu\left(-\alpha \nabla f\left(\boldsymbol{x}_{k-1}\right)+\mu \boldsymbol{v}_{k-2}\right) \\ = & -\alpha \nabla f(\boldsymbol{x})-\alpha \mu \nabla f\left(\boldsymbol{x}_{k-1}\right)+\mu^{2} \boldsymbol{v}_{k-2} \\ = & -\alpha \nabla f(\boldsymbol{x})-\alpha \mu \nabla f\left(\boldsymbol{x}_{k-1}\right)+\mu^{2}\left(-\alpha \nabla f\left(\boldsymbol{x}_{k-2}\right)+\mu \boldsymbol{v}_{k-3}\right) \\ = & -\alpha \nabla f(\boldsymbol{x})-\alpha \mu \nabla f\left(\boldsymbol{x}_{k-1}\right)-\alpha \mu^{2} \nabla f\left(\boldsymbol{x}_{k-2}\right)+\mu^{3} v_{k-3} \\ & \vdots \\ = & -\alpha \nabla f(\boldsymbol{x})-\alpha \mu \nabla f\left(\boldsymbol{x}_{k-1}\right)-\alpha \mu^{2} \nabla f\left(\boldsymbol{x}_{k-2}\right)-\alpha \mu^{3} \nabla f\left(\boldsymbol{x}_{k-3}\right) \cdots \end{aligned}$

更新优化变量值时使用动量项代替负梯度项，梯度下降更新公式为

$\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}+\boldsymbol{v}_{k}$

动量项加快了梯度下降法的收敛速度，它使用历史信息对当前梯度值进行修正，消除病态条件问题上的来回振荡

下图显示了用动量项梯度下降法求解 $0.1 x_{1}^{2}+2 x_{2}^{2}$ 极值时的迭代过程，与上图相比，迭代代序列更为平滑，且收敛更快

xxx图使用动量项后的迭代轨迹

AdaGrad

标准梯度下降法的步长值难以确定，且优化变量的各个分量采用了相同的步长，AdaGrad(Adaptive Gradient)算法根据前几轮迭代时的历史梯度值动态计算步长值，且优化向量的每一个分量都有自已的步长，梯度下降迭代公式为

$\left(\boldsymbol{x}_{k+1}\right)_{i}=\left(\boldsymbol{x}_{k}\right)_{i}-\alpha \frac{\left(\boldsymbol{g}_{k}\right)_{i}}{\sqrt{\sum_{j=1}^{k}\left(\left(\boldsymbol{g}_{j}\right)_{i}\right)^{2}+\varepsilon}}$

其中 $\alpha$ 是人工设定的全局学习率，$g_{k}$是第$k$次迭代时的梯度向量，$\varepsilon$是为避免除0操作而增加的接近于0的正数，$i$为向量的分量下标，这里的计算针对向量的每个分量分别进行

与标准梯度下降法相比，上式多了分母项，它累积了到本次迭代为止的梯度的历史值信息，用于计算步长值

历史导数值的绝对值越大，在该分量上的学习率越小，反之越大

虽然实现了自适应学习率，但这种算法还存在问题：需要人工设置全局学习率$\alpha$；随着时间的累积，上式中的分母会越来越大，导致学习率趋向于0，优化变量无法有效更新

RMSProp

RMSProp算法是对AdaGrad的改进，避免了长期累积梯度值所导致的学习率趋向于0的问题

算法维持一个梯度平方累加值的向量$E\left[\boldsymbol{g}^{2}\right]$，其初始值为$\mathbf{0}$，更新公式为

$E\left[\boldsymbol{g}^{2}\right]_{k}=\delta E\left[\boldsymbol{g}^{2}\right]_{k-1}+(1-\delta) \boldsymbol{g}_{k}^{2}$

这里的$g^{2}$是对梯度向量的每个分量分别进行平方，$\delta$是人工设定的衰减系数

不同于AdaGrad直接累加所有历史梯度的平方和，RMSProp将历史梯度平方值按照系数$\delta$指数级衰减之后再累加，即使用了移动指数加权平均

梯度下降法更新公式为

$\left(\boldsymbol{x}_{k+1}\right)_{i}=\left(\boldsymbol{x}_{k}\right)_{i}-\alpha \frac{\left(\boldsymbol{g}_{k}\right)_{i}}{\sqrt{\left(E\left[\boldsymbol{g}^{2}\right]_{k}\right)_{i}+\varepsilon}}$

$\alpha$是人工设定的全局学习率，标准梯度下降方法相比，这里也只多了一个分母项

AdaDelta

AdaDelta算法也是对AdaGrad的改进，避免了长期累积梯度值所导致的学习率趋向于0的问题，还去掉了对人工设置全局学习率的依赖

算法定义了两个向量，初始值均为0

$E\left[\boldsymbol{g}^{2}\right]_{0}=\mathbf{0} \quad E\left[\Delta \boldsymbol{x}^{2}\right]_{0}=\mathbf{0}$

$E\left[g^{2}\right]$是梯度平方(对每个分量分别平方)的累计值，与RMSProp算法相同，更新公式为

$E\left[\boldsymbol{g}^{2}\right]_{k}=\rho E\left[\boldsymbol{g}^{2}\right]_{k-1}+(1-\rho) \boldsymbol{g}_{k}^{2}$

$g^{2}$是向量每个元素分别计算平方，后面所有的计算公式都是对向量的每个分量分别进行计算，接下来计算RMS向量

$\operatorname{RMS}[\boldsymbol{g}]_{k}=\sqrt{E\left[\boldsymbol{g}^{2}\right]_{k}+\varepsilon}$

然后计算优化变量的更新值

$\Delta \boldsymbol{x}_{k}=-\frac{\operatorname{RMS}[\Delta \boldsymbol{x}]_{k-1}}{\operatorname{RMS}[\boldsymbol{g}]_{k}} \boldsymbol{g}_{k}$

这里 $\mathrm{RMS}[\Delta x]_{k-1}$ 根据 $E\left[\Delta \boldsymbol{x}^{2}\right]$ 计算，计算公式与 $RMS[\boldsymbol{g}]_{k}$ 相同

这个更新值同样通过梯度来构造，但学习率是通过梯度的历史值确定的

$E\left[\Delta x^{2}\right]$是优化变量更新值的平方累加值，它们的更新公式为

$E\left[\Delta x^{2}\right]_{k}=\rho E\left[\Delta x^{2}\right]_{k-1}+(1-\rho) \Delta x_{k}^{2}$

在这里， $\Delta x_{k}^{2}$ 是对 $\Delta x_{k}$ 的每个分量进行平方。梯度下降的达代公式为

$\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}+\Delta \boldsymbol{x}_{k}$

Adam

Adam(Adaptive Moment Estimation)算法整合了自适应学习率与动量项

算法用梯度构造了两个向量$m$和$v$，初始值为$0$，更新公式为

$\begin{aligned} \left(m_{k}\right)_{i} & =\beta_{1}\left(m_{k-1}\right)_{i}+\left(1-\beta_{1}\right)\left(g_{k}\right)_{i} \\ \left(v_{k}\right)_{i} & =\beta_{2}\left(\boldsymbol{v}_{k-1}\right)_{i}+\left(1-\beta_{2}\right)\left(\boldsymbol{g}_{k}\right)_{i}^{2} \end{aligned}$

这里 $\beta_{1}$ 、 $\beta_{2}$ 是人工置顶的参数，梯度下降的迭代公式为

$\left(\boldsymbol{x}_{k+1}\right)_{i}=\left(\boldsymbol{x}_{k}\right)_{i}-\alpha \frac{\sqrt{1-\beta_{2}^{k}}}{1-\beta_{1}^{k}} \frac{\left(\boldsymbol{m}_{k}\right)_{i}}{\sqrt{\left(\boldsymbol{v}_{k}\right)_{i}}+\varepsilon}$

$m$的作用相当于于动量项，$v$用于构造学习率

随机梯度下降法

在机器学习中，目标函数通常定义在一个训练样本集上

假设训练样本集有$N$个样本，机器学习模型在训练时优化的目标是这个数据集上的平均损失函数

$L(\boldsymbol{w})=\frac{1}{N} \sum_{i=1}^{N} L\left(\boldsymbol{w},\boldsymbol{x}_{i},y_{i}\right)$

其中 $L\left(\boldsymbol{w}，\boldsymbol{x}_{i}，y_{i}\right)$ 是对单个训练样本 $\left(\boldsymbol{x}_{i}，y_{i}\right)$ 的损失函数， $\boldsymbol{w}$ 是机哭学习模型需要学习的参数，是优化变量

显然 $\nabla L(w)=\frac{1}{N} \sum_{i=1}^{N} \nabla L\left(w，x_{i}，y_{i}\right)$ ，因此计算目标函数梯度时需要计算对每个训练样本损失函数的梯度，然后求均值

如果训练时每次达代都用所有样本，那么计算成本太高

作为改进，可以在每次迭代时选取一批样本，将损失函数定义在这些样本上，作为整个样本集的损失函数的近似值

小批量随机梯度下降法(Mini Batch Gradient Descent Method在每次迭代时使用上面目标函数的随机逼近值，只使用$M \ll N$个样本来近似计算损失函数，在每次迭代时，要优化的目标函数变为

$L(w) \approx \frac{1}{M} \sum_{i=1}^{M} L\left(w,x_{i},y_{i}\right)$

随机梯度下降法在数学期望的意义下收敛，随机采样产生的梯度的期望值是真实的梯度

在具体实现时，每次先对所有训练样本进行随机洗牌，打乱顺序；然后将其均匀分成多份，每份$M$个样本；接下来依次用每一份执行梯度下降法迭代

一种特殊情况是$M=1$，每次迭代只使用一个训练样本

随机梯度下降法并不能保证每次迭代后目标函数值下降，事实上，每次迭代时使用的是不同的目标函数

但通常情况下目标函数的整体趋势是下降的，能够收敛到局部极值点处

下图是用随机梯度下降法训练神经网络时损失函数的曲线，横轴为迭代次数，纵轴为损失函数的值

可以看到，迭代时函数的值会出现振荡，但整体趋势是下降的，最后收敛

xxx图用随机梯度下降法训练神经网络时的损失函数曲线

除具有实现效率高的优点之外，随机梯度下降法还会影响收敛的效果

对于深度神经网络，随机梯度下降法比批量梯度下降法更容收敛到一个好的极值点处

应用一工神经网络

假设有 $N$ 个训练样本 $\left(x_{i}，y_{i}\right)$ ，其中 $x_{i}$ 为输入向量， $y_{i}$ 为标签向量

训练的目标是最小化样本标签值与神经网络预测值之间的误差，如果使用均方误差，则优化的目标为

$L(\boldsymbol{w})=\frac{1}{2 N} \sum_{i=1}^{N}\left\|h\left(\boldsymbol{x}_{i}\right)-\boldsymbol{y}_{i}\right\|^{2}$

其中$w$为神经网络所有参数的集合，包括各层的权重和偏置，$h(x)$是神经网络实现的映射

这个最优化问题是一个不带约束条件的问题，可以用梯度下降法求解

如果计算出了损失函数对参数的梯度值，梯度下降法第$k+1$次迭代时参数的更新公式为

$\boldsymbol{w}_{k+1}=\boldsymbol{w}_{k}-\alpha \nabla L\left(\boldsymbol{w}_{k}\right)$

梯度值的计算通过反向传播算法实现

参数的初始化是一个需要考虑的问题，一般用随机数进行初始化

如果训练样本数很大，那么通常采用随机梯度下降法，通常情况下，随机梯度下降法有很好的收敛效果

二阶优化算法

梯度下降法只利用了一阶导数信息，收敛速度慢，通常情况下，利用二阶导数信息可以加快收敘速度，典型代表是牛顿法和拟牛顿法

牛顿法在每个迭代点处将目标函数近似为二次函数，然后通过求解梯度为0的方程得到迭代方向

牛顿法在每次迭代时需要计算梯度向量与黑塞矩阵，并求解一个线性方程组，计算量大且面临黑塞矩阵不可逆的问题

拟牛顿法是对它的改进，算法构造出一个矩阵作为黑塞矩阵或其逆矩阵的近似

牛顿法

牛顿法(Newton Method)寻找目标函数作二阶近似后梯度为$\mathbf{0}$的点，逐步逼近极值点

根据费马引理，函数在点$x$处取得极值的必要条件是梯度为0

$\nabla f(\boldsymbol{x})=\mathbf{0}$

直接计算函数的梯度然后解上面的方程组通常很困难，和梯度下降法类似，可以采用迭代法近似求解

对目标函数在$x_{0}$处作二阶泰勒展开

$f(\boldsymbol{x})=f\left(\boldsymbol{x}_{0}\right)+\nabla f\left(\boldsymbol{x}_{0}\right)^{\mathrm{T}}\left(\boldsymbol{x}-\boldsymbol{x}_{0}\right)+\frac{1}{2}\left(\boldsymbol{x}-\boldsymbol{x}_{0}\right)^{\mathrm{T}} \nabla^{2} f\left(\boldsymbol{x}_{0}\right)\left(\boldsymbol{x}-\boldsymbol{x}_{0}\right)+o\left(\left\|\boldsymbol{x}-\boldsymbol{x}_{0}\right\|^{2}\right)$

忽略二次以上的项，将目标函数近似成二次函数，等式两边同时对$\boldsymbol{x}$求梯度，可得

$\nabla f(\boldsymbol{x}) \approx \nabla f\left(\boldsymbol{x}_{0}\right)+\nabla^{2} f\left(\boldsymbol{x}_{0}\right)\left(\boldsymbol{x}-\boldsymbol{x}_{0}\right)$

其中 $\nabla^{2} f\left(x_{0}\right)$ 为在 $x_{0}$ 处的黑塞矩阵

从上面可以看出，这里至少要展开到二阶，如果只有一阶，那么无法建立梯度为$\mathbf{0}$的方程组，因为此时一次近似函数的梯度值为常数

令函数的梯度为$\mathbf{0}$，有

$\nabla f\left(x_{0}\right)+\nabla^{2} f\left(x_{0}\right)\left(\boldsymbol{x}-x_{0}\right)=\mathbf{0}$

解这个线性方程组可以得到

$\boldsymbol{x}=\boldsymbol{x}_{0}-\left(\nabla^{2} f\left(\boldsymbol{x}_{0}\right)\right)^{-1} \nabla f\left(\boldsymbol{x}_{0}\right)$

如果将梯度向量简写为$\boldsymbol{g}$，黑塞矩阵简记为$\boldsymbol{H}$，上式可以简写为

$x=x_{0}-H^{-1} g$

由于在泰勒公式中忽略了高阶项将函数进行了近似，因此这个解不一定是目标函数的驻点，需要反复用上式进行迭代

从初始点$x_{0}$处开始，计算函数在当前点处的黑塞矩阵和梯度向量，然后用下面的公式进行迭代

$\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}-\alpha \boldsymbol{H}_{k}^{-1} \boldsymbol{g}_{k}$

直至收敛到驻点处

即在毎次迭代之后，在当前点处将目标函数近似成二次函数，然后寻找梯度为0的点

$-H^{-1} g$称为牛顿方向，迭代终止的条件是梯度的模接近于0，或达到指定的达代次数

牛顿法的流程如下算法所示

🧮牛顿法

初始化$x_{0}$，$k=0$

while $k<N$ do

计算当前点处的梯度值$g_k$以及黑塞矩阵$H_k$

if $\left|g_{k}\right|<$ eps then

停止迭代

end if

迭代 $d_{k}=-H_{k}^{-1} g_{k}$

迭代 $\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}+\alpha \boldsymbol{d}_{k}$

$k=k+1$

end while

其中$\alpha$是人工设置的学习率，需要学习率的原因与梯度下降法相同，是为了保证能够忽略泰勒公式中的高阶无穷小项

如果目标函数是二次函数，则黑塞矩阵是一个常数矩阵，且泰靿公式中的高阶项为0

对于任意给定的初始点$\boldsymbol{x}_{0}$，牛顿法只需要一次迭代即可收敛到驻点

与梯度下降法不同，牛顿法无法保证每次迭代时目标函数值下降。为了确定学习率的值，通常使用直线搜索技术

具体做法是让$\alpha$取一些典型的离散值，如下面的值

$0.0001,0.001,0.01$

选择使得 $f\left(\boldsymbol{x}_{k}+\alpha \boldsymbol{d}_{k}\right)$ 最小化的步长值作为最优步长，保证迭代之后的函数值充分下降

与梯度下降法相比，牛顿法有更快的收敛速度，但每次迭代的成本也更高

1️⃣每次迭代时需要计算梯度向量与黑塞矩阵，并计算黑塞矩阵的逆矩阵，最后计算矩阵与向量乘积

实现时通常不直接求黑塞矩阵的逆矩阵，而是求解如下方程组

$\boldsymbol{H}_{k} \boldsymbol{d}=-\boldsymbol{g}_{k}$

求解线性方程组可使用迭代法，如共轭梯度法

2️⃣另一个问题是黑塞矩阵可能不可逆，从而导致其失效

拟牛顿法

拟牛顿法(Quasi-Newton Methods)对牛顿法存在的问题进行了改进

其核心思路是不精确计算目标函数的黑塞矩阵然后求逆矩阵，而是通过其他手段得到黑塞矩阵的逆

具体做法是构造一个近似黑塞矩阵或其逆矩阵的正定对称矩阵，用该矩阵进行牛顿法迭代

由于要推导下一个迭代点 $x_{k+1}$ 的黑塞矩阵需要满足的条件，并建立与上一个迭代点 $x_{k}$ 处的函数值、导数值之间的关系，以指导近似矩阵的构造，因此需要在 $\boldsymbol{x}_{k+1}$ 点处作泰勒展开，并将 $x_{k}$ 的值代入泰勒公式

将函数在$x_{k+1}$点处作二阶泰勒展开，忽略高次项，有

$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)$

上式两边同时对$x$取梯度，可以得到

$\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)$

如果令$x=x_{k}$，则有

$\nabla f\left(\boldsymbol{x}_{k+1}\right)-\nabla f\left(\boldsymbol{x}_{k}\right) \approx \nabla^{2} f\left(\boldsymbol{x}_{k+1}\right)\left(\boldsymbol{x}_{k+1}-\boldsymbol{x}_{k}\right)$

将梯度向量与黑塞矩阵简写，则有

$\boldsymbol{g}_{k+1}-\boldsymbol{g}_{k} \approx \boldsymbol{H}_{k+1}\left(\boldsymbol{x}_{k+1}-\boldsymbol{x}_{k}\right)$

如果令

$\boldsymbol{s}_{k}=\boldsymbol{x}_{k+1}-\boldsymbol{x}_{k} \qquad \boldsymbol{y}_{k}=\boldsymbol{y}_{k+1}-\boldsymbol{g}_{k}$

则上式可简写为

$y_{k} \approx \boldsymbol{H}_{k+1} s_{k}$

这里的 $s_{k}$ 和 $y_{k}$ 都可以根据之前的迭代结果直接算出，如果 $H_{k+1}$ 可逆，那么上式等价于

$\boldsymbol{s}_{k} \approx \boldsymbol{H}_{k+1}^{-1} \boldsymbol{y}_{k}$

上两式称为拟牛顿条件，用于近似代替黑塞矩阵和它的逆矩阵的矩阵需要满足该条件

利用该条件，根据上一个迭代点 $x_{k}$ 和当前迭代点 $x_{k+1}$ 的值以及这两点处的梯度值，就可以近似计算出当前点的黑塞矩阵或其逆矩阵

由于黑塞矩阵与它的逆矩阵均对称，因此它们的近似矩阵也要求是对称的

此外，通常还要求近似矩阵正定，拟牛顿法通过各种方法构造出满足上述条件的近似矩阵

下面介绍典型的实现：DFP算法以及BFGS算法

DFP算法

问题的核心是构造黑塞矩阵或其逆矩阵的近似矩阵$H_{k}$，保证满足拟牛顿条件

首先为该矩阵设定初始值，然后在每次迭代时更新此近似矩阵

$\boldsymbol{H}_{k+1}=\boldsymbol{H}_{\boldsymbol{k}}+\boldsymbol{E}_{k}$

其中 $E_{k}$ 称为校正矩阵，现在的任务变为寻找该矩阵，根据上式，如果以 $H_{k}$ 充当黑塞矩阵逆矩阵的近似，有

$\left(\boldsymbol{H}_{k}+\boldsymbol{E}_{k}\right) \boldsymbol{y}_{k}=\boldsymbol{s}_{k}$

上式变形后得到

$\boldsymbol{E}_{k} \boldsymbol{y}_{k}=\boldsymbol{s}_{k}-\boldsymbol{H}_{k} \boldsymbol{y}_{k}$

DFP算法采用了这种思路，DFP(Davidon-Fletcher-Powell)算法以其3位发明人的名字命名

算法构造黑塞矩阵逆矩阵的近似(Inverse Hessian Approximation)，其初始值为单位矩阵$I$，每次迭代时按照下式更新该矩阵

$\boldsymbol{H}_{k+1}=\boldsymbol{H}_{k}+\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}$

即校正矩阵为

$\boldsymbol{E}_{k}=\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}$

其中 $\boldsymbol{u}_{k}$ 和 $\boldsymbol{v}_{k}$ 为待定的 $n$ 维向量， $\alpha_{k}$ 和 $\beta_{k}$ 为待定的系数

显然，按照上式构造的$\boldsymbol{H}_{k}$是一个对称矩阵根据上几个式子，校正矩阵必须满足

$\left(\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}\right) \boldsymbol{y}_{k}=\boldsymbol{s}_{k}-\boldsymbol{H}_{k} \boldsymbol{y}_{k}$

即

$\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}=\boldsymbol{s}_{k}-\boldsymbol{H}_{k} \boldsymbol{y}_{k}$

此方程的解不唯一，可以取某些特殊值从而简化问题的求解，这里令

$\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}=\boldsymbol{s}_{k} \quad \beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}=-\boldsymbol{H}_{k} \boldsymbol{y}_{k}$

同时令

$v_k=s_k \qquad v_k = H_ky_k$

将这两个解代人上面的两个方程，可以得到

$\alpha_{k} \boldsymbol{s}_{k} \boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}=\alpha_{k} \boldsymbol{s}_{k}\left(\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}\right)=\alpha_{k}\left(\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}\right) \boldsymbol{s}_{k}=\boldsymbol{s}_{k}$

以及

$\begin{aligned} \beta_{k} \boldsymbol{H}_{k} \boldsymbol{y}_{k}\left(\boldsymbol{H}_{k} \boldsymbol{y}_{k}\right)^{\mathrm{T}} \boldsymbol{y}_{k} & =\beta_{k} \boldsymbol{H}_{k} \boldsymbol{y}_{k} \boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}=\beta_{k} \boldsymbol{H}_{k} \boldsymbol{y}_{k} \boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{\boldsymbol{k}} \boldsymbol{y}_{k} \\ & =\beta_{k} \boldsymbol{H}_{\boldsymbol{k}} \boldsymbol{y}_{k}\left(\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k} \boldsymbol{y}_{k}\right)=\beta_{k}\left(\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k} \boldsymbol{y}_{k}\right) \boldsymbol{H}_{k} \boldsymbol{y}_{k}=-\boldsymbol{H}_{k} \boldsymbol{y}_{k} \end{aligned}$

上面两个结果利用了矩阵乘法的结合律以及 $\boldsymbol{H}_{k}$ 是对称矩阵这一条件，在这里 $\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}$ 与 $\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k} y_{k}$ 均为标量，从而解得

$\alpha_{k}=\frac{1}{\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{y}_{k}} \quad \beta_{k}=-\frac{1}{\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k} \boldsymbol{y}_{k}}$

将上面的解代人式 (4.16)，由此得到矩阵$H_{k}$的更新公式

$\boldsymbol{H}_{k+1}=\boldsymbol{H}_{k}-\frac{\boldsymbol{H}_{k} \boldsymbol{y}_{k} \boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k}}{\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{H}_{k} \boldsymbol{y}_{k}}+\frac{\boldsymbol{s}_{k} \boldsymbol{s}_{k}^{\mathrm{T}}}{\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}}$

此更新公式可以保证$H_{k}$的对称正定性

每次迭代时，得到矩阵$H_{k}$之后用牛顿法进行更新，由于构造的是黑塞矩阵逆矩阵的近似，因此可以直接将其与梯度向量相乘从而得到牛顿方向

DFP算法的流程如下算法所示

🧮DFP算法

初始化$x_{0}$，$H_0$，$k=0$

while $k<N$ do

迭代 $d_{k}=-H_{k} g_{k}$

用直线搜索得到步长$\lambda_{k}$

迭代 $s_{k}=\lambda_{k} d_{k}，x_{k+1}=x_{k}+s_{k}$

if $\left|g_{k+1}\right|<eps$ then

结束循环

end if

迭代 $y_{k}=g_{k+1}-g_{k}$

迭代 $\boldsymbol{H}_{k+1}=\boldsymbol{H}_{k}-\frac{\boldsymbol{H}_{k} y_{k} \nu_{k}^{\mathrm{T}} H_{k}}{\boldsymbol{y}_{k}^{\mathrm{T}} H_{k} \boldsymbol{y}_{k}}+\frac{\boldsymbol{s}_{k} \boldsymbol{s}_{k}^{\mathrm{T}}}{\boldsymbol{y}_{k}^{T} \theta_{k}}$

$k=k+1$

end while

如果用单位矩阵初始化$\boldsymbol{H}_{k}$，则第一次迭代时

$\boldsymbol{d}_{0}=-\boldsymbol{H}_{0} \boldsymbol{g}_{0}=-\boldsymbol{I} \boldsymbol{g}_{0}=-\boldsymbol{g}_{0}$

这相当于使用梯度下降法，后面逐步细化$H_{k}$，使其更精确地通近当前点处黑塞矩阵的逆矩阵。

BFGS算法

BFGS(Broyden-Fletcher-Goldfarb-Shanno)算法以其4位发明人的名字命名

算法构造黑塞矩阵的一个近似矩阵$B_{k}$并用下式迭代更新这个矩阵

$\boldsymbol{B}_{k+1}=\boldsymbol{B}_{k}+\Delta \boldsymbol{B}_{k}$

该矩阵的初始值 $B_{0}$ 为单位阵 $I$ ，要解决的问题就是每次的校正矩阵 $\Delta B_{k}$ 的构造

根据前面的式子，黑塞矩阵的近似矩阵$B_{k}$需要满足

$\left(\boldsymbol{B}_{k}+\Delta \boldsymbol{B}_{k}\right) \boldsymbol{s}_{k}=\boldsymbol{y}_{k}$

与DFP算法相同，校正矩阵构造为如下形式

$Delta \boldsymbol{B}_{k}=\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}$

将其代人式 (4.20)，可以得到

$\left(\boldsymbol{B}_{k}+\alpha_{k} \boldsymbol{u}_{k} \boldsymbol{u}_{k}^{\mathrm{T}}+\beta_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}\right) \boldsymbol{s}_{k}=\boldsymbol{y}_{k}$

整理后可得

$\alpha_{k}\left(\boldsymbol{u}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}\right) \boldsymbol{u}_{k}+\beta_{k}\left(\boldsymbol{v}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}\right) \boldsymbol{v}_{k}=\boldsymbol{y}_{k}-\boldsymbol{B}_{k} \boldsymbol{s}_{k}$

同样，可以取这个方程的一组特殊解，这里直接令

$\alpha_{k}\left(\boldsymbol{u}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}\right) \boldsymbol{u}_{k}=\boldsymbol{y}_{k} \qquad \beta_{k}\left(v_{k}^{\mathrm{T}} \boldsymbol{s}_{k}\right) \boldsymbol{v}_{k}=-\boldsymbol{B}_{k} \boldsymbol{s}_{k}$

同时令两个向量为

$\boldsymbol{u}_{k}=\boldsymbol{y}_{k} \quad \boldsymbol{v}_{k}=\boldsymbol{B}_{k} \boldsymbol{s}_{\boldsymbol{k}}$

将它们的代入上面两个方程，可得

$\alpha_{k}\left(\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}\right) \boldsymbol{y}_{k}=\boldsymbol{y}_{k}$

以及

$\beta_{k}\left(\boldsymbol{B}_{k} s_{k}\right)^{\mathrm{T}} \boldsymbol{s}_{k} \boldsymbol{B}_{k} s_{k}=\beta_{k} \boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k}^{\mathrm{T}} \boldsymbol{s}_{k} \boldsymbol{B}_{k} \boldsymbol{s}_{k}=\beta_{k}\left(\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k} \boldsymbol{s}_{k}\right) \boldsymbol{B}_{k} \boldsymbol{s}_{k}=-\boldsymbol{B}_{k} \boldsymbol{s}_{k}$

从而解得两个系数为

$\alpha_{k}=\frac{1}{\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}} \qquad \beta_{k}=-\frac{1}{\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k} \boldsymbol{s}_{k}}$

由此得到校正矩阵为

$\Delta \boldsymbol{B}_{k}=\frac{\boldsymbol{y}_{k} \boldsymbol{y}_{k}^{\mathrm{T}}}{\boldsymbol{y}_{k}^{\mathrm{T}} \boldsymbol{s}_{k}}-\frac{\boldsymbol{B}_{k} \boldsymbol{s}_{k} \boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k}}{\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k} \boldsymbol{s}_{k}}$

如果初始值 $B_{0}$ 是正定矩阵，且在每次迭代时 $y_{k}^{\mathrm{T}} s_{k}>0$ ，则每次更新后得到的 $B_{k}$ 都是正定的

由于BFGS算法构造的是黑塞矩阵的近似，因此还需要求解方程组以得到牛顿方向

而$B_{k}$是正定对称矩阵，可以采用高效的方法求解此线性方程组

比较DFP算法与BFGS算法可以发现，二者的校正矩阵计算公式互为对偶，将 $\boldsymbol{s}_{k}$ 与 $\boldsymbol{y}_{k}$ 的角色进行了对换

BFGS算法的流程如下算法所示

🧮BFGS算法

初始化$x_{0}$，$B_0=I$，$k=0$

while $k<N$ do

迭代 $d_{k}=-B_{k}^{-1} g_{k}$

用直线搜索得到步长$\lambda_{k}$

迭代 $\boldsymbol{s}_{k}=\lambda_{k} \boldsymbol{d}_{k}，\boldsymbol{x}_{k+1}=\boldsymbol{x}_{k}+\boldsymbol{s}_{k}$

if $\left|g_{k+1}\right|<eps$ then

结束循环

end if

迭代 $y_{k}=g_{k+1}-g_{k}$

迭代 $B_{k+1}=B_{k}+\frac{y_{k} y_{k}^{\mathrm{T}}}{\boldsymbol{y}_{k}^{\mathrm{T}} s_{k}}-\frac{\boldsymbol{B}_{k} s_{k} s_{k}^{\mathrm{T}} \boldsymbol{B}_{k}}{\boldsymbol{s}_{k}^{\mathrm{T}} \boldsymbol{B}_{k} \boldsymbol{s}_{k}}$

$k=k+1$

end while

BFGS 算法在每次达代时需要计算$n \times n$的矩阵$\boldsymbol{B}_{k}$，当$n$很大时，存储该矩阵将耗费大量内存

为此，提出了改进方案L-BFGS算法(有限存储的BFGS算法)，其思想是不存储完整的矩阵 $B_{k}$ ，只存储向量 $s_{k}$ 和 $y_{k}$

对于大多数目标函数，BFGS 算法有很好的收敛效果

下图是用L-BFGS算法求解$x^{2}+y^{2}$极值的迭代过程，算法只需要达代4次即可收敛到极小值点处

xxx图L-BFGS算法求解$x^{2}+y^{2}$极值的迭代过程

分治法

分治法是算法设计中常用的思路，它把一个问题拆分成多个子问题，通常情况下，子问题更容易求解

在求得子问题的解之后，将其合并得到整个问题的解

在用于最优化方法时的通行做法是每次只优化部分变量，将高维优化问题分解为低维优化问题

坐标下降法

坐标下降法(Coordinate Descent)是分治法的典型代表

对于多元函数的优化问题，坐标下降法每次只对一个分量进行优化，将其他分量固定不动

算法依次优化每一个变量，直至收敛，假设要求解的优化问题为

$\min _{\boldsymbol{x}} f(\boldsymbol{x})，\boldsymbol{x}=\left(x_{1}，x_{2}，\cdots，x_{n}\right)$

算法在每次迭代时依次选择 $x_{1}, \cdots, x_{n}$ 进行优化，求解单个变量的优化问题

完整的算法流程如下算法所示

🧮坐标下降法

初始化$x_{0}$

while 没有收敛 do

for i=1,2, $\cdots$, n do

求解 $min_{x_i}{f(x)}$

end for

end while

算法每次迭代时在当前点处沿一个坐标轴方向进行一维搜索，固定其他坐标轴方向对应的分量，求解一元函数的极值

在整个过程中依次循环使用不同坐标轴方向对应的分量进行达代，更新这些分量的值，一个周期的一维搜索迭代过程相当于一次梯度下降法造代完成对优化变量每个分量的一次更新

坐标下降法的求解过程如下图所示，这里是二维优化问题

xxx图坐标下降法求解二维优化问题的迭代过程

在每次连代时，首先固定$y$轴分量，优化$x$轴分量；然后固定$x$轴分量，优化$y$轴分量，整个优化过程在各个坐标轴方向之间轮换

坐标下降法具有计算量小、实现效率高的优点，在机器学习领域得到了成功的应用

典型的是求解线性模型的训练问题，在开源库liblinear中有实现，此外，在求解非负矩阵分解问题中也有应用

缺点

坐标下降法的缺点是对非光滑(不可导)的多元目标函数可能无法进行有效处理，以及难以并行化

考虑下面的目标函数

f(x, y)=|x+y|+3|y-x|

如果当前迭代点为$(-2,-2)$，单独在$x$轴和$y$轴方向进行迭代均无法保证目标函数值下降，此时坐标下降法失效

如下图所示，图中画出了目标函数的等高线，$(0,0)$是极小值点，在$(-2,-2)$点处，单独改变$x$或$y$的值均不能使目标函数值下降

xxx图 4.16 坐标下降法失效的例子

SMO算法

SMO(Sequential Minimal Optimization)算法是求解支持向量机对偶问原的高效算法

算法的核心思想是每次从优化变量中挑出两个分量进行优化，让其他分量固定，这样能保证满足等式约束条件

假设训练样本集为 $\left(\boldsymbol{x}_{i}, y_{i}\right), i=1, \cdots, l$ ，其中 $\boldsymbol{x} \in \mathbb{R}^{n}$ 为样本的特征向量， $y=\pm 1$ 为标签值

在使用核函数之后，支持向量机在训练时求解的对偶问题为

$\begin{array}{l} \min _{\alpha} \frac{1}{2} \sum_{i=1}^{l} \sum_{j=1}^{l} \alpha_{i} \alpha_{j} y_{i} y_{j} K\left(x_{i}, x_{j}\right)-\sum_{i=1}^{l} \alpha_{i} \\ 0 \leqslant \alpha_{i} \leqslant C \\ \sum_{j=1}^{l} \alpha_{j} y_{j}=0 \end{array}$

$C$为惩罚因子，是人工设定的正常数，核矩阵的元素为

$K_{i j}=K\left(\boldsymbol{x}_{i}, \boldsymbol{x}_{j}\right)$

这里 $K\left(\boldsymbol{x}_{i},\boldsymbol{x}_{j}\right)$ 为核函数，将两个向量 $\boldsymbol{x}_{\boldsymbol{i}},\boldsymbol{x}_{j}$ 映射为一个实数

该对偶问题是二次规划问题，问题规模由训练样本数$l$决定，其值很大时，常规的求解算法将面临计算效率低和存储空间占用太大的问题

SMO算法每次选择两个变量进行优化，假设选取的两个分量为 $\alpha_{i}$ 和 $\alpha_{j}$ ，其他分量都固定，当成常数

由于 $y_{i} y_{i}=1$ 、 $y_{j} y_{j}=1$ ，对这两个变量的目标函数可以写成

$f\left(\alpha_{i}，\alpha_{j}\right)=\frac{1}{2} K_{i i} \alpha_{i}^{2}+\frac{1}{2} K_{j j} \alpha_{j}^{2}+s K_{i j} \alpha_{i} \alpha_{j}+y_{i} v_{i} \alpha_{i}+y_{j} v_{j} \alpha_{j}-\alpha_{i}-\alpha_{j}+c$

其中$c$是一个常数。这里定义

$s=y_{i} y_{j} \qquad v_{i}=\sum_{k=1，k \neq i，k \neq j}^{l} y_{k} \alpha_{k}^{*} K_{i k} \qquad v_{j}=\sum_{k=1，k \neq i，k \neq j}^{l} y_{k} \alpha_{k}^{*} K_{j k}$

这里的$\alpha^{*}$为$\alpha$在上一轮迭代后的值

子问题的目标函数是二元二次函数，可以直接给出最小值的解析解，这个问题的约束条件为

$0 \leqslant \alpha_{i} \leqslant C \qquad 0 \leqslant \alpha_{j} \leqslant C \qquad y_{i} \alpha_{i}+y_{j} \alpha_{j}=-\sum_{k=1，k \neq i，k \neq j}^{l} y_{k} \alpha_{k}=\xi$

利用上面的等式约束可以消掉 $\alpha_{i}$ ，从而只剩下变量 $\alpha_{j}$ ，目标函数是 $\alpha_{j}$ 的二次函数，直接求得解析解

分阶段优化

AdaBoost算法在训练时同样采取了分治法的策略，每次迭代时先训练弱分类器，然后确定弱分类器的权重系数

AdaBoost算法在训练时的目标是最小化指数损失函数

假设强分类器为$F(x)$，$x \in \mathbb{R}^{n}$为特征向量，$y=\pm 1$为标签值，单个训练样本的指数损失函数为

$L(y，F(\boldsymbol{x}))=\exp (-y F(\boldsymbol{x}))$

强分类器是弱分类器的加权组合，定义为

$F(\boldsymbol{x})=\sum_{i=1}^{M} \beta_{i} f_{i}(\boldsymbol{x})$

其中 $f_{i}(\boldsymbol{x})$ 是第 $i$ 个弱分类器， $\beta_{i}$ 是其权重， $M$ 为弱分类器的数量

训练时依次训练每个弱分类器，将其加人强分类器中

将强分类器的计算公式代人上面的损失函数中，得到训练第$j$个弱分类器时对整个训练样本集的训练损失函数为

$\left(\beta_{j}，f_{j}\right)=\arg \min _{\beta，f} \sum_{i=1}^{l} \exp \left(-y_{i}\left(F_{j-1}\left(\boldsymbol{x}_{i}\right)+\beta f\left(\boldsymbol{x}_{i}\right)\right)\right)$

这里将强分类器拆成两部分

第一部分是之前的迭代已经得到的强分类器$F_{j-1}$
第二部分是当前要训练的弱分类器$f$与其权重$\beta$的乘积对训练样本的损失函数

前者在之前的迭代中已经求出，因此可以看成常数，上式目标函数可以简化为

$\min _{\beta，f} \sum_{i=1}^{l} w_{i}^{j-1} \exp \left(-\beta y_{i} f\left(\boldsymbol{x}_{i}\right)\right)$

其中

$w_{i}^{j-1}=\exp \left(-y_{i} F_{j-1}\left(\boldsymbol{x}_{i}\right)\right)$

它只和前面迭代得到的强分类器有关，与当前的弱分类器、弱分类器权重无关，这就是样本权重

上上式的问题可以分两步求解，首先将$\beta$看成常数

由于 $y_{i}$ 和 $f\left(\boldsymbol{x}_{i}\right)$ 的取值只能为 $+1$ 或 $-1$ ，且样本权重非负，要让上上式的目标函数最小化，必须让二者相等

因此损失函数对$f(x)$的最优解为

$f_{j}=\arg \min _{f} \sum_{i=1}^{l} w_{i}^{j-1} \boldsymbol{I}\left(y_{i} \neq f\left(\boldsymbol{x}_{i}\right)\right)$

其中$I$是指标函数，如果括号里的条件成立，其值为1，否则为0

上式的最优解是使得对样本的加权误差最小的弱分类器，在得到弱分类器之后，式子的优化目标可以表示成$\beta$的函数

$L(\beta)=\mathrm{e}^{-\beta} \times \sum_{y_{i}=f_{j}\left(x_{i}\right)} w_{i}^{j-1}+\mathrm{e}^{\beta} \times \sum_{y_{i} \neq J_{j}\left(x_{i}\right)} w_{i}^{j-1}$

上式前半部分是被当前的弱分类器正确分类的样本，此时 $y_{i}$ 与 $f\left(x_{i}\right)$ 同号且 $y_{i} f\left(x_{i}\right)=1$ ， $\exp \left(-\beta y_{i} f\left(x_{i}\right)\right)=\exp (-\beta)$ ，后半部分是被当前的弱分类器错误分类的样本，这种情况有 $y_{i} f\left(\boldsymbol{x}_{i}\right)=-1$ ， $\exp \left(-\beta y_{i} f\left(\boldsymbol{x}_{i}\right)\right)=\exp (\beta)$

目标函数可以进一步写成

$L(\beta)=\left(\mathrm{e}^{\beta}-\mathrm{e}^{-\beta}\right) \times \sum_{i=1}^{l} w_{i}^{j-1} I\left(y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)\right)+\mathrm{e}^{-\beta} \times \sum_{i=1}^{l} w_{i}^{j-1}$

具体推导过程为

$\begin{array}{l} \mathrm{e}^{-\beta} \cdot \sum_{y_{i}=f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1}+\mathrm{e}^{\beta} \cdot \sum_{y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1} \\ =\mathrm{e}^{-\beta} \cdot \sum_{y_{i}=f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1}+\mathrm{e}^{-\beta} \cdot \sum_{y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1}-\mathrm{e}^{-\beta} \cdot \sum_{y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1}+\mathrm{e}^{\beta} \cdot \sum_{y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1} \\ =\mathrm{e}^{-\beta} \cdot \sum_{i=1}^{l} w_{i}^{j-1}+\left(\mathrm{e}^{\beta}-\mathrm{e}^{-\beta}\right) \cdot \sum_{y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)} w_{i}^{j-1} \\ =\mathrm{e}^{-\beta} \cdot \sum_{i=1}^{l} w_{i}^{j-1}+\left(\mathrm{e}^{\beta}-\mathrm{e}^{-\beta}\right) \cdot \sum_{i=1}^{l} w_{i}^{j-1} I\left(y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)\right) \end{array}$

函数在极值点的导数为0，对$\beta$求导并令其为0

$\left(\mathrm{e}^{\beta}+\mathrm{e}^{-\beta}\right) \cdot \sum_{i=1}^{l} w_{i}^{j-1} I\left(y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)\right)-\mathrm{e}^{-\beta} \cdot \sum_{i=1}^{l} w_{i}^{j-1}=0$

上式两边同除以 $\sum_{i=1}^{l} w_{i}^{j-1}$ ，由此得到关于 $\beta$ 的方程

$\left(\mathrm{e}^{\beta}+\mathrm{e}^{-\beta}\right) \cdot \mathrm{err}_{j}-\mathrm{e}^{-\beta}=0$

最后得到最优解为

$\beta=\frac{1}{2} \ln \frac{1-\operatorname{err}_{j}}{\operatorname{err}_{j}}$

其中$err_{j}$为当前弱分类器对训练样本集的加权错误率

$\operatorname{err}_{j}=\left(\sum_{i=1}^{l} w_{i}^{j-1} I\left(y_{i} \neq f_{j}\left(\boldsymbol{x}_{i}\right)\right)\right) /\left(\sum_{i=1}^{l} w_{i}^{j-1}\right)$

在得到当前弱分类器之后，对强分类器进行更新

$F_{j}(\boldsymbol{x})=F_{j-1}(\boldsymbol{x})+\beta_{j} f_{j}(\boldsymbol{x})$

下次迭代时样本的权重为

$w_{i}^{j}=w_{i}^{j-1} \cdot \mathrm{e}^{-\beta_{j} y_{i} f_{j}\left(x_{i}\right)}$

即AdaBoost训练算法中的样本权重更新公式

应用logistic回归

下面介绍坐标下降法在求解logistic回归训练问题中的应用，给定$l$个训练样本 $\left(x_{i}，y_{i}\right)$ ， $i=1，\cdots，l$ ，其中 $x_{i} \in \mathbb{R}^{n}$ 为特征向量， $y_{i}=\pm 1$ 为标签值

L2正则化logistic回归的对偶问题为

$\begin{array}{l} \min _{\boldsymbol{\alpha}} D_{L R}(\boldsymbol{\alpha})=\frac{1}{2} \boldsymbol{\alpha}^{\mathrm{T}} \boldsymbol{Q} \boldsymbol{\alpha}+\sum_{i: \alpha_{i}>0} \alpha_{i} \ln \alpha_{i}+\sum_{i: \alpha_{i}<C}\left(C-\alpha_{i}\right) \ln \left(C-\alpha_{i}\right) \\ 0 \leqslant \alpha_{i} \leqslant C，i=1，\cdots，l \end{array}$

其中$C$为惩罚因子，矩阵$Q$定义为

$Q_{i j}=y_{i} y_{j} x_{i}^{\mathrm{T}} \boldsymbol{x}_{j}$

如果定义

$0 \log 0=0$

它与下面的极限是一致的

$\lim \limits_{x \rightarrow 0^{+}} x \ln x=0$

最小化式可以简化为

$\begin{array}{l} \min _{\alpha} D_{L R}(\boldsymbol{\alpha})=\frac{1}{2} \boldsymbol{\alpha}^{\mathrm{T}} \boldsymbol{Q} \boldsymbol{\alpha}+\sum_{i=1}^{l}\left(\alpha_{i} \ln \alpha_{i}+\left(C-\alpha_{i}\right) \ln \left(C-\alpha_{i}\right)\right) \\ 0 \leqslant \alpha_{i} \leqslant C，i=1，\cdots，l \end{array}$

上式的目标函数中带有对数函数，可以采用坐标下降法求解

与其他最优化方法(如共轭梯度法、拟牛顿法)相比，坐标下降法有更快的迭代速度，更适合大规模问题的求解

下面我们介绍带约束条件的坐标下降法求解思路

在用坐标下降法求解时，采用了一个技巧，不是直接优化一个变量，而是优化该变量的增量

假设本次迭代时要优化 $\alpha_{i}$ ，将其他的 $\alpha_{j},j \neq i$ 固定不动

假设本次迭代之后 $\alpha_{i}$ 的值为 $\alpha_{i}+z$ ，在这里 $\alpha_{i}$ 为一个常数，是 $\alpha_{i}$ 上一次迭代后的值

用 $\alpha_{i}+z$ 替换 $\alpha_{i}$ ，最小化式子的目标函数以及不等式约束可以写成 $z$ 的函数

$\begin{array}{l} \min _{z} g(z)=\left(c_{1}+z\right) \ln \left(c_{1}+z\right)+\left(c_{2}-z\right) \ln \left(c_{2}-z\right)+\frac{a}{2} z^{2}+b z \\ -c_{1} \leqslant z \leqslant c_{2} \end{array}$

其中所有常数定义为

$c_{1}=\alpha_{i},c_{2}=C-\alpha_{i},a=Q_{i i},b=(Q \alpha)_{i}$

因为目标函数含有对数函数，上面的函数是一个超越函数，无法给出极值的公式解

采用牛顿法求解上面的问题，不考虑不等式约束条件 $-c_{1} \leqslant z \leqslant c_{2}$ ，迭代公式为

$z_{k+1}=z_{k}+d \qquad d=-\frac{g^{\prime}\left(z_{k}\right)}{g^{\prime \prime}\left(z_{k}\right)}$

这是牛顿法对一元函数的情况，梯度为一阶导数，黑塞矩阵为二阶导数，子问题目标函数的一阶导数和二阶导数分别为

$g^{\prime}(z)=a z+b+\ln \frac{c_{1}+z}{c_{2}-z} \qquad g^{\prime \prime}(z)=a+\frac{c_{1}+c_{2}}{\left(c_{1}+z\right)\left(c_{2}-z\right)}$

为了保证牛顿法收敛，还需要使用直线搜系，检查迭代之后的函数值是否充分下降