拟牛顿法推导

·4795·10 分钟·
AI摘要: 本文介绍了拟牛顿法的两种主要形式:BFGS法和L-BFGS法。BFGS法利用曲率信息来预处理梯度,从而避免了传统方法中对Hessian矩阵进行完整计算的需要。L-BFGS法进一步优化了这种预处理,通过仅保存最近m次迭代的曲率信息来计算Hessian矩阵的近似值,显著减少了内存使用和计算量。

拟牛顿方程

f(xk)=f(xk+1)+Gk+1(xkxk+1)感觉可以当成割线法 \nabla f(x_k) = \nabla f(x_{k + 1}) + G_{k+1} (x_k - x_{k+1})\\ 感觉可以当成割线法

BFGS法

BFGS是拟牛顿法的一种,主要是利用曲率信息来对梯度进行预处理来确定下降方向

曲率信息则是通过维护一个使用广义的割线法逐步近似的关于损失函数的 Hessian矩阵 来获得。

近似与校正

Gk+1Bk+1G_{k+1} \approx B_{k+1},也就是说,用“割线法”算出来的数据近似Hessian矩阵,既然是近似必定有误差,下面尝试纠正这个误差

Bk+1=Bk+Bk(Bk是割线法算出来的) B_{k+1} = B_k + \nabla B_k\\ (B_{k}是割线法算出来的)

当然,这只是假设出来的,是否真的能通过这个方式校正假设还需要进一步计算。

一个很有意思的特征是:构造出来的近似Hessian矩阵BkB_k一定是正定的;

牛顿法中,如果Hessian矩阵不是正定的,那目标函数就不一定是下降的,只有正定才一定下降。拟牛顿法正好避免了这个问题,由于一定正定,所以目标函数一定下降

为了方便描述,记:

sk=xk+1xkyk=gk+1gk s_k = x_{k + 1} - x_k\\ y_k = g_{k + 1} - g_k

那么有:

yk=Bk+1sk y_k = B_{k+1}s_k

Bk\nabla B_k 分解:

Bk=aukukT+βvkvkT \nabla B_k = a u_ku_k^T + \beta v_kv_k^T

只要解出了 uvu和v 两个向量,那么就能算出 Bk\nabla B_k, 那么就能算出Bk+1B_{k + 1} .

容易知道:

yk=Bk+1sk=(Bk+Bk)sk=(Bk+aukukT+βvkvkT)sk=Bksk+(aukTsk)uk+(βvkTsk)vkykBksk=(aukTsk)uk+(βvkTsk)vk \begin{align*} y_k &= B_{k + 1}s_k\\ &= (B_{k} + \nabla B_k)s_k\\ &= (B_k + au_ku_k^T + \beta v_kv_k^T)s_k\\ &= B_ks_k + (au_k^Ts_k)u_k + (\beta v_k^T s_k) v_k\\ \Rightarrow & y_k - B_ks_k = (au_k^T s_k)u_k +(\beta v_k^T s_k)v_k \end{align*}

可以发现,该方程十分宽泛,能解出很多解,不妨取个特殊情况(真的是随机猜的条件)试试:

aukTsk=1,βvkTsk=1 au_k^T s_k = 1, \beta v_k^T s_k = -1

那么就能推出:

uk=yk,vk=Bksk u_k = y_k , v_k = B_k s_k

只要能解出个解就行,这里引入了一些额外条件,解出了个还算漂亮的递推结果。继续推下去:

a=1ukTsk=1ykTskβ=1vkTsk=1skTBksk a = \frac{1}{u_k^T s_k} = \frac{1}{y_k^T s_k}\\ \beta = - \frac{1}{v_k^T s_k} = -\frac{1}{s_k^T B_k s_k}

推完了,合并下得到最终的迭代公式:

Bk+1=Bk+ykykTykTskBkskskTBkskBksk B_{k+1} = B_k + \frac{y_k y_k^T}{y_k^T s_k} - \frac{B_k s_k s_k^T B_k}{s_k B_k s_k}

由于牛顿法的方向时Hk+11gk-H_{k+1}^{-1}g_k , 所以,最好能直接求出Hessian矩阵的逆,根据Sherman-Morrison公式得到:

Bk+11=Bk1+(1skTyk+ykTBk1yk(skTyk)2)skskT1skTyk(Bk1ykskT+skykTBk1) B_{k+1}^{-1} = B_{k}^{-1} + \left(\frac{1}{s_k^T y_k} + \frac{y_k^T B_k^{-1}y_k}{(s_k^T y_k)^2}\right)s_k s_k^T - \frac{1}{s_k^T y_k}(B_k^{-1} y_k s_k^T + s_k y_k^T B_{k}^{-1})

伪代码

设置起始点x0{\displaystyle \mathbf {x} _{0}}和初始的Hessian矩阵B0{\displaystyle B_{0}},重复以下步骤,xk{\displaystyle \mathbf {x} _{k}} 会收敛到问题的解:

  1. 通过求解方程pk=Bkf(xk){\mathbf {p} _{k}=-\displaystyle B_{k}\nabla f(\mathbf {x} _{k})},获得下降方向pk{\displaystyle \mathbf {p} _{k}}

  2. pk{\displaystyle \mathbf {p} _{k}}方向上进行一维的优化线搜索,找到合适的步长αk{\displaystyle \alpha _{k}}。如果这个搜索是完全的,则αk=argminf(xk+αpk){\displaystyle \alpha _{k}=\arg \min f(\mathbf {x} _{k}+\alpha \mathbf {p} _{k})} 。在实际中,不完全的搜索一般就足够,此时只要求αk{\displaystyle \alpha _{k}}满足Wolfe条件。

  3. sk=αkpk{\displaystyle \mathbf {s} _{k}=\alpha _{k}\mathbf {p} _{k}},并且令xk+1=xk+sk{\displaystyle \mathbf {x} _{k+1}=\mathbf {x} _{k}+\mathbf {s} _{k}}。 (十分正常普通的梯度迭代)

  4. yk=f(xk+1)f(xk){\displaystyle \mathbf {y} _{k}={\nabla f(\mathbf {x} _{k+1})-\nabla f(\mathbf {x} _{k})}} (保存迭代后的梯度变化量)

  5. Bk+1=Bk+(1skTyk+ykTBkyk(skTyk)2)skskT1skTyk(BkykskT+skykTBk) B_{k+1} = B_{k} + \left(\frac{1}{s_k^T y_k} + \frac{y_k^T B_ky_k}{(s_k^T y_k)^2}\right)s_k s_k^T - \frac{1}{s_k^T y_k}(B_k y_k s_k^T + s_k y_k^T B_{k}) (这里的BkB_k是近似Hessian矩阵的)

f(x){\displaystyle f(\mathbf {x} )}表示要最小化的目标函数。可以通过检查梯度的范数 f(xk){\displaystyle ||\nabla f(\mathbf {x} _{k})||}来判断收敛性。如果B0{\displaystyle B_{0}}初始化为B0=I{\displaystyle B_{0}=I},第一步将等效于梯度下降,但接下来的步骤会受到近似于Hessian矩阵的Bk{\displaystyle B_{k}}的调节。

L-BFGS

BFGS的好处

容易知道,牛顿法实际上是要求初始迭代点和最优点相距比较近才行,否则可能收敛到其他极小值点,所以才有牛顿分形这种东西。

那么拟牛顿法由于不是完全的牛顿法,所以也就是说,它的收敛情况不是很需要限定成“离最优点比较近”,实际上,“离最有点很近”这种概念本来就很虚,所以拟牛顿法既有牛顿法的收敛速度,又没有那么限制使用区域,这就是BFGS的好处。

BFGS的坏处与L-BFGS

BFGS作为牛顿法,每次迭代都要保存一个Hessian矩阵,空间复杂度O(n2)O(n^2),这让牛顿法(BFGS之类的)无法应用于具有百万参数的现代神经网络模型

有困难就有解决办法,于是就有L-BFGS

L-BFGS 全称 储存受限的BFGS,L是指 limited memory

L-BFGS中非常重要的一点是:迭代初值的设置为单位矩阵

迭代公式和BFGS一样:

Bk+11=Bk1+(1skTyk+ykTBk1yk(skTyk)2)skskT1skTyk(Bk1ykskT+skykTBk1) B_{k+1}^{-1} = B_{k}^{-1} + \left(\frac{1}{s_k^T y_k} + \frac{y_k^T B_k^{-1}y_k}{(s_k^T y_k)^2}\right)s_k s_k^T - \frac{1}{s_k^T y_k}(B_k^{-1} y_k s_k^T + s_k y_k^T B_{k}^{-1})

但是,在BFGS中,每次都只计算skyks_k和y_k,每次都只保留Bk1B_{k}^{-1}, 可是这个Bk1B_{k}^{-1}是非常大的,难以保存。

(上面的说法可能合理、可能不合理,现在脑子有点乱)

为了计算Bk1B_{k}^{-1},不再是以往所有的sk,yks_k, y_k的运算和,而是最近m次的 (这个m是人为指定的)

在第 k 次迭代,算法求得了 xkx_k ,并且保存的曲率信息为si,yi(k1)...(km) s_i,yi_{(k−1)...(k−m)}。为了得到Hk H_k,算法每次迭代均需选择一个初始的矩阵 Hk0H_k^0,这是不同于 BFGS 算法的一个地方,接下来只用最近的 m 个向量对该初始矩阵进行修正,实践中 Hk0H^0_k 的设定通常如下:

Hk0=rkErk=sk1Tyk1yk1Tyk1 H_k^0 = r_k E\\ r_k = \frac{s_k - 1^T y_{k-1}}{y_{k -1}^T y_{k-1}}

其中 rkr_k 表示比例系数,它利用最近一次的曲率信息来估计真实Hessian矩阵的大小,这就使得当前步的搜索方向较为理想,而不至于跑得“太偏”,这样就省去了步长搜索的步骤,节省了时间。

在L-BFGS算法中,通过保存最近 m 次的曲率信息来更新近似矩阵的这种方法在实践中是很有效的,虽然 L-BFGS 算法是线性收敛,但是每次迭代的开销非常小,因此 L-BFGS 算法执行速度还是很快的,而且由于每一步迭代都能保证近似矩阵的正定,因此算法的鲁棒性还是很强的。

总结

总结下 BFGS 与 L-BFGS 的:

  • BFGS算法在运行的时候,每一步迭代都需要保存一个 n×nn×n 的矩阵,现在很多机器学习问题都是高维的,当 n 很大的时候,这个矩阵占用的内存是非常惊人的,并且所需的计算量也是很大的,这使得传统的 BFGS 算法变得非常不适用。而 L-BFGS 则是很对这个问题的改进版,从上面所说可知,BFGS 算法是通过曲率信息 (sk,yk)(sk,yk)来修正 HkH_k 从而得到 Hk+1H_{k+1}

  • L-BFGS 算法的主要思路是:算法仅仅保存最近 m 次迭代的曲率信息来计算 Hk+1H_{k+1} 。这样,我们所需的存储空间就从 n×nn×n 变成了 2m×n2m×n 而通常情况下 m<<n 。

Kaggle学习赛初探