数值优化（7）——限制空间的优化算法：LBFGS，LSR1

学弱猹

发布于 2021-08-09 16:43:44

8790

发布于 2021-08-09 16:43:44

————————————————————————————————————

大家好！首先跟大家说一件事情，就是受到疫情的影响，我今年将不前往CMU攻读硕士，而开始入职从事偏DS的算法工程师的工作～最近因为这些事情也比较忙，在更新的进度上有了较大的滞后，这里向大家表示抱歉！在这一年专栏的更新会有什么样的一个倾斜，目前我还没有一个头绪。但我相信，这丝毫不会影响大家对于专栏更新的期待哈哈哈。

这一节我们主要会关注用于求解大规模问题的优化算法。这一节也是有关无约束优化的最后一节。之后我们会开始介绍一些更加偏应用和实际的带约束规划的内容。

需要强调的是，我们的这一个系列关注的更多的是最优化，不是凸优化，会更加偏重于数值算法，而不是凸分析的理论。意味着大家可能最关心的有关拉格朗日乘数法和对偶形式的内容，其实不会在这个系列出现（但是KKT条件是会有的，因为我们很多的算法依赖于这个条件）。但是没有关系，凸优化的相关内容必然也会在专栏中出现，只是可能需要很长一段时间，所以大家倒也不必太担心。如果对于对偶理论等内容比较好奇，可以参考一下这个

https://zhuanlan.zhihu.com/p/52544386

好的不说废话了，我们开始吧。

Source

R. H. Byrd, H. F. Khalfan, and R. B. Schnabel. Analysis of a symmetric rank-one trust region method.
J. E. Dennis and R. B. Schnabel. Numerical methods for unconstrained optimization and nonlinear equations
J. Nocedal and S. J. Wright. Numerical Optimization
课堂笔记，教授主页：https://www.math.fsu.edu/~whuang2/index.html

引入：什么是大规模问题

提到大规模问题必然要提到机器学习这个领域。比方说对于回归模型，如果我们利用梯度下降法去优化函数

f(\hat y) = \frac 12 (y - \hat y)^2

那么必然会依赖到所有的数据点（因为

\hat y - y

就是它的梯度，

\hat y

是一个与所有样本都有关的量）。那么如果你每一次的迭代，都需要计算到所有的样本点，当样本数到达足够大的量级的时候，这个时间的开销就会难以承受。

要解决这个问题，工业界一般使用的是SGD（随机梯度下降）和mini-batch（小批量下降），相当于说每一次更新只会采取一个样本/一小部分样本。这样的思想除了在梯度下降，在聚类中心的计算中（也就是K-means）也会有用到，具体的可以参考吴恩达Coursera的课程。

但是大规模问题并不是只有一种的，事实上，如果我们说优化模型长下面这个样子

\min_{x \in \mathbb{R}^n}\sum_{i=1}^N f_i(x)

那么之前的例子就是

会特别大。但是如果是

特别大呢？这种情况下，对于每一个数据点，它都涵盖很多信息，这个时候很多方法就不再奏效了。比方说在第5节所提到的牛顿CG方法，这个方法的实现可以依赖海塞矩阵的作用（就是矩阵与向量的相乘，这样子可以避免直接计算海塞矩阵，比方说

x^T A x

这个量，如果你知道了

和

就可以计算，并不需要知道

是什么）。所以如果说海塞矩阵比较稀疏，计算作用的复杂度就会很低，那么这个方法就可行。反过来，如果海塞矩阵很稠密，那就完全的不可行了。

事实上，按照这个思维来看，对于上一节我们刚提到的拟牛顿法，很多时候就不奏效了，比方说BFGS方法就需要依赖我们每一步计算的对称正定矩阵，如果说

比较大，那么储存这个矩阵需要依赖

O(n^2)

的一个储存空间，这在现实中是做不到的。

到此我们算是介绍了大规模问题的主要的两种情况，我们后面的方法主要关注的问题，都是第二种问题。

限制空间的拟牛顿方法

拟牛顿方法算是在一次收敛速度和二次收敛速度之间找了一个平衡，换句话说，我们希望对方法做一些修改，使得新的方法依然能够保证收敛速度，且可以应用于大规模问题。这就产生了限制空间的拟牛顿方法（Limited Memory Quasi-Newton Methods）。

我们先回顾一下拟牛顿法。对于一般的拟牛顿法，其核心思想可以写成

B_{k + 1} = f(B_k, s_k, y_k)

其中

s_k = x_{k + 1} - x_k

，

y_k = \nabla f(x_{k + 1}) - \nabla f(x_k)

。也就是说，通过给定了这几个量，我们会希望推出下一个矩阵

B_{k + 1}

是什么。有了这个思路，我们自然会想，既然是限制空间了，那么可以产生的信息一定会受限，也就是说，如果我们以迭代的思维把这个函数写开，正常情况下就是

B_{k + 1} = f(f(f(\cdots (f(B_0, s_0, y_0)))))

现在我们要限制空间，其含义就是保存部分数据。那么这个时候，迭代就会变成

B_{k + 1} = f(f(f(\cdots (f(B_{k - m+ 1}, s_{k - m + 1}, y_{k - m + 1})))))

你可以看出来，我们相当于保留了最近的

组数据。对于更加靠后的数据，我们就丢弃不要了。这样的想法其实也不难理解，相信你知道牛顿法的局部二次收敛速度只有在迭代点足够靠近极值点的时候才会存在。而优化算法的迭代过程其实刚开始一般都会偏离的很远，所以丢弃一些足够远的数据无伤大雅。

丢弃数据之后必然会产生一个问题就是这里的

B_{k - m + 1}

应该是什么？很显然，我们需要找一个已知的矩阵代替它。这就是限制空间拟牛顿法的思路。

限制空间的BFGS方法

还是一样，我们先提一下BFGS方法。如果变成了限制空间的情况，这个方法就会变成LBFGS。

在上一节我们提过一般情况下的BFGS方法，它利用的是线搜索方法的框架。那么在这里我们依然考虑的是这个框架，所以我们要选取的搜索方向其实就是

p_k = -B_k^{-1}\nabla f(x_k) = -H_k \nabla f(x_k)

有一件非常离谱的事情在于，通过两次循环即可计算出这个搜索方向，并且不需要计算出

H_k

。我们简单提一下这个思路和算法，但不会给出详细的计算过程。

回顾一下上一节所提到的带逆形式的BFGS更新公式

H_{k + 1} = (I_n - \frac{s_ky_k^T}{y_k^Ts_k})H_k(I_n- \frac{y_ks_k^T}{y_k^Ts_k}) + \frac{s_ks_k^T}{y_k^Ts_k}

我们令

\rho_k = \frac 1{y_k^T s_k}

，

V_k = I_n - \frac{y_ks_k^T}{y_k^T s_k}

，我们就可以得到

H_{k + 1} = V_k^T H_k V_k + \rho_k s_k s_k^T

我们如果代入

H_k

，就可以得到

H_{k + 1} = V_k^T (V_{k - 1}^T H_{k - 1} V_{k - 1} + \rho_{k - 1}s_{k - 1}s_{k - 1}^T) V_k + \rho_k s_k s_k^T

所以可以一直按照这个思路写下去，一直写到最后的

H_{k - m + 1}

，再靠后的数据就丢弃了。而这个矩阵我们会用一个简单的

H_{k + 1}^{(0)}

代替。

好的，接下来我们祭出这个计算的方法，这个方法就可以很好地计算出

H_{k + 1}\nabla f(x_{k + 1})

这个量

它的计算复杂度为

O(8mn)

。我们建议大家代入表达式一步一步地计算一遍。如果你没有亲手一步一步的算一遍，你可能也不会相信这个方法和我们上面写的东西居然是完全一样的……

那么现在，问题就落在了我们的对于

H_{k+1}^{(0)}

，也就是

H_{k - m + 1}

的选取上。一般来说，取单位阵的倍数（也就是数量矩阵）是方便计算的，所以一般来说，我们设置为

H_{k + 1}^{(0)} = \gamma_{k + 1}I_n

其中

\gamma_{k + 1}

为BB步长（见第3节）。使用BB步长的原因在于，如果我们将其理解为我们正常的最速下降法，取

H_{k + 1}^{(0)}

其实本质上就是在取初始步长。因此取成BB步长是有依据的，当然了它在之后的证明中也起到了作用。

一般来说，BB2是使用的比较多的，也就是一般来说，我们有

\gamma_{ k + 1} = \frac{s_k ^T y_k}{y_k ^T y_k}

这里还有个细节是我们使用的是最近一步的BB步长。这里虽然我们对应的矩阵

H_{k + 1}^{(0)}

不是最近一步，但是我们希望这个我们代替的矩阵，可以更好地反映当前我们迭代的函数的海塞矩阵特征值的情况，这就是我们如此定义的原因。

那么到此为止，我们就可以放出它的算法部分了。

这里的two-loop recursion，就是我们上面所写的那个循环。

最后，我们来看一下LBFGS方法的收敛性分析。我们主要还是想说的是它的全局收敛性。体现在下面这个定理中。

Theorem 1: 设

x_0

为初值，

f \in C^2

，

\mathcal{N}_{x_0} = \{x: f(x) \le f(x_0)\}

凸，并且存在

m, M > 0

，使得

L\|z\|^2 \le z^T \nabla^2 f(x) z \le U\|z\|^2

对任意的

z \in \mathbb{R}^n, x \in \mathcal{N}_{x_0}

成立。那么LBFGS更新会使得

\{x_k\}

收敛到

的极小值点

x^*

。

我们证明一下这个结论。事实上，主要的证明过程都和上一节证明BFGS的过程相似。

我们考虑同样的BFGS的原版更新（也就是说，我们考虑的是

B_k

的公式，不是

H_k

的），这个时候自然就就会有

\operatorname{tr}(B_{k + 1}) = \operatorname{tr}(\bar B_k) - \frac{\|\bar B_ks_k\|^2}{s_k^T\bar B_ks_k} + \frac{\|y_k\|^2}{y_k^Ts_k}

|B_{k + 1}| = |\bar B_k| \frac{y_k^T s_k}{s_k^T \bar B_k s_k}

这个地方要注意的是我们把

B_k

改成了

\bar B_k

。这是因为我们的更新公式不再是更新到

B_0

，而是更新到

H_{k+1}^{(0)}

就停止了。所以为了体现这个变化，我们用

\bar B_k

表示它。

设

\cos\theta_k = \frac{s_k^T \bar B_k s_k}{\|s_k\|\|\bar B_ks_k\|}

，

q_k = \frac{s_k^T \bar B_ks_k}{s_k^Ts_k}

（除了改了下标记以外，都和之前的证明没有任何本质区别），并设

\Psi(B) = \operatorname{tr} (B) - \ln (\det(B))

，那么同样的思路我们可以推导出

\Psi(B_{k + 1}) = \operatorname{tr}(\bar B_k) + M_k - \frac{q_k}{\cos ^2\theta_k} - \ln(\det(\bar B_k)) - \ln m_k + \ln q_k \\ \le \Psi(B_k) + (M_k - \ln m_k - 1) + \ln \cos^2 \theta_k

因此我们迭代

步（只保留

步的信息），就可以得到

0 < \Psi(B_{k + 1}) \le \Psi(B_{k + 1}^{(0)}) + m (U - \ln L - 1) + \sum_{j = k - m + 1}^k \ln \cos^2 \theta_j

这里

B_{k + 1}^{(0)} = H_{k + 1}^{(0) ^{-1}}

。

到此情况开始有一些变化了。我们思考一下之前的证明是如何做的？它是利用迭代不可能无限做下去，使得最终推出

\cos \theta_j \not \to 0

。那么在这里，我们自然也是希望说明这一个条件是成立的，那会怎么利用上这个定理呢？

首先我们注意到其实这个条件的本质就是说明搜索方向不能取得太差。换句话说，我们的搜索方向对应的对称正定的矩阵不能够太过于病态，也就是条件数不能没有上界。这个思路会引导我们去使用第2节的Theorem 5

数值优化（2）——线搜索：步长选取条件的收敛性

也就是要证明它的条件数存在上界。而我们看到，其实

\Psi(\cdot)

这个函数本质就是一个条件数的估计，因为

\Psi(X)

这个函数如果有界，那么

的各个特征值就一定有界（想想为什么？），所以基于这个思考，只要证明

\Psi(B_{k + 1})

有界就可以了，而这是非常容易的，因为我们使用的

H_{k + 1}^{(0)}

是一个数量矩阵，这个系数是BB步长，它的本质是一个特征值的加权平均（要理解这个结论，可以考虑二次凸问题的情况，这个情况下

y_k = A s_k

，代入即可发现结论）。因此我们有理由相信存在

，使得

\Psi(B_{k + 1}) \le C

通过这个，我们就可以推出来

\sum_{i=1}^n (\lambda_i - \ln \lambda_i) \le C

，其中

\lambda_i

为这个矩阵的各个特征值。也就可以推出来，存在

\delta_1, \delta_2

，使得

\lambda_{min} > \delta_1, \lambda_{max} < \delta_2

这也就是证明了存在

，使得

\|B_{k+1}\|\|B_{k+1}^{-1}\| \le T

，那么自然全局收敛性就得到了保证。

这个部分的证明其实很多部分都是上一节内容。而后半部分的推论思考则需要花时间消化一下。关于它的局部收敛性，其实结论和上一节的结论无异，也就是说局部超线性收敛速度依然存在。

最后提一下它的实操。其实我们可以看出，要求海塞矩阵凸的条件依然是存在的，也就是说非凸问题它是不保证全局收敛性的。因此上一节提到的各种修正，到这一节依然是需要继续使用的，具体的可以参考上一节的相关内容。

到此，我们就介绍好了LBFGS的全部内容。

限制空间的SR1方法

我们在上一节有提到过SR1方法，事实上限制空间的SR1方法（也即LSR1方法）也在近期被提出，但是远不如LBFGS历史悠久了。它也是存在一个很有趣的展开表示使得它可以通过修正，在大规模的情况下被使用。

我们回顾一下SR1方法的更新公式，其实就是

B_{k + 1} = B_k + \frac{(y_k - B_k s_k)(y_k - B_k s_k)^T}{(y_k - B_ks_k)^T s_k}

如果我们设

\phi_k = y_k - B_k s_k

，那么会有

B_{k + 1} = B_{k- m + 1} + \frac{\phi_{k - m + 1}\phi_{k - m + 1}^T}{\phi_{k - m + 1}^T s_{k - m + 1}} + \cdots + \frac{\phi_k\phi_k^T}{\phi_k^T s_k}

非常牛逼的地方在于这个公式也被人发现了一个通用的限制空间版本的更新公式，就是下面的这个定理。

Theorem 2: 设

Y_{k + 1} = \begin{bmatrix}y_{k - m + 1} & \ldots & y_k\end{bmatrix}

，

S_{k + 1} = \begin{bmatrix}s_{k - m + 1} & \ldots & s_k\end{bmatrix}

，且

D_{k + 1} = \operatorname{diag}(s_{k -m + 1}^T y_{k - m + 1}, \ldots, s_{k }^T y_{k })

，且

L_{k + 1} \in \mathbb{R}^{m \times m}

定义为

(L_{k + 1})_{i,j} = \begin{cases}s_{k - m + i}^T y_{k - m + j} & i >j \\ 0 & \operatorname{otherwise}\end{cases}

设

B_{k - m + 1}

为对称阵，且

\phi_i^T s_i \ne 0, i = k - m + 1, \ldots, k

，那么有下面的公式成立

B_{k + 1} = B_{k - m + 1} + (Y_ {k + 1} - B_{k - m + 1}S_{k + 1})(D_{k + 1} + L_{k + 1} + L_{k + 1}^T - S_{k + 1}^T B_{k - m + 1}S_{k + 1})^{-1}(Y_{k + 1} - B_{k -m +1}S_{k + 1})^T

并且

D_{k + 1} + L_{k + 1} + L_{k + 1}^T - S_{k + 1}^T B_{k - m + 1}S_{k + 1}

可逆。

这个定理我们不证明，具体的可以参考Source上的BKS96（对应年份和名字首字母）（好吧我似乎应该很早就把每一个地方要用的source文献标记一下位置……）。

我们可以看出，这个方法与我们上面写的LBFGS方法不同的地方在于，我们这里给出了一个直接可以一步到位，从

B_{k - m + 1}

到

B_{k + 1}

的公式。所以我们的目标自然也就是修改

B_{k - m + 1}

为

B_{k + 1}^{(0)}

。有了LBFGS的经验，我们的修改方法也应运而出，考虑BB步长。现实中一般采用

\tilde \gamma_{k + 1} = \frac 1 {\alpha_{k + 1}^{B1}} = \frac{y_k^T s_k}{s_k^T s_k}

这里要注意的是，因为我们关注的问题一般都不是纯凸的问题，所以很多情况下，好的，它的更新公式介绍完了，但是它的算法过程不是很简单。上一节我们说过，SR1方法一般是使用信赖域框架解决的。换句话说我们如何解对应的信赖域子问题？

我们看一下我们的信赖域子问题，其实本质上就是下面这个问题

\min_{ \| p \| \le \Delta} p^T g + \frac 12 p^T B p

但是这里的

B = \tilde \gamma I_n + \Phi M^{-1} \Phi^T

，其中

M \in \mathbb{R}^{m \times m}

对称非奇异，

\Phi \in \mathbb{R}^{n \times m}

（我们修改了标记，所以要对应上面的那个定理的表达式来看）。所以要解这个子问题，联想到第5节的Proposition 1，我们其实就需要找到这样的

\lambda \ge 0

，使得

(B + \lambda I)p^* = -\nabla f(x)

\lambda (\Delta - \|p^*\|) = 0

B + \lambda I

为对称半正定的矩阵。

事实上，这就相当于一方面找到一个合适的

p^*

，也就是求

p(\lambda) = -(B + \lambda I)^{-1} g

。再进一步找到一个合适的

\lambda

，也就是求解

h(\lambda) = \frac 1{\|p(\lambda)\|} - \frac 1 {\Delta}

。你也能看出来，总的来说就是希望

h(\lambda) = 0

。

同时要保证

\|p^*\|

是可行解。

要化简上面的式子，这个地方的处理确实需要一点技巧。我们注意到

的一部分是一个半正定的矩阵，所以可以做的一件事情是对它做一个合同变换，使得它的特征值被暴露出来。为什么要做这件事我们到之后再解释。

为了找到合同变换，我们需要想尽办法构造一些正交矩阵，所以我们的第一件事就是做QR分解。这里想到QR分解，是因为我们的矩阵

\Phi

并不是一个方阵。再然后，我们保留

，对剩下的部分再做一次谱分解。具体来说就是

\tilde \gamma I_n + \Phi M^{-1} \Phi^{T} = \tilde \gamma I_n + QR M^{-1} R^T Q^T = \tilde \gamma I_n + QUSU^TQ^T

然后考虑令

P = QU

即可。

但是对于

来说，构造并没有结束。我们有说过，如果说我们有

X^TX = ODO^T

为一般的谱分解，那么

X^TX + kI = O(D + kI)O^T

，但是这是要

的大小和

X^TX

的相同才可以做到。在这里因为

的大小并不是

而是

（想想为什么），所以必须要对特征值补0，才可以进行我们上面的操作。在这个时候，

因为只是列正交矩阵，要想做到谱分解那一步，还需要补上正交补空间的那一部分。总体来说我们可以写成

B = [P \quad P_\perp] (\tilde \gamma I + \begin{bmatrix}S & \\ & 0_{(n - m) \times (n - m)}\end{bmatrix})[P \quad P_\perp]^T

这个时候，新的矩阵

[P \quad P_\perp]

就是一个成熟的正交矩阵啦！为了讨论的方便，后面我们把

(\tilde \gamma I + \begin{bmatrix}S & \\ & 0_{(n - m) \times (n - m)}\end{bmatrix})

的这一部分拉出来，假设

的

个特征值是按从小到大重新排列的，重新写成一个特征值矩阵

\operatorname{diag}(\sigma_1, \cdots, \sigma_m, \sigma_{m+1}, \sigma_{m + 1}, \cdots, \sigma_{m + 1})

。

在做了特征值分解之后，我们所关注的函数的书写都不再抽象复杂了。因为这个时候，我们有

p(\lambda) = [P \quad P_\perp] \operatorname{diag}(\frac 1{\sigma_1 + \lambda}, \cdots, \frac 1{\sigma_m + \lambda}, \frac 1{\sigma_{m + 1} + \lambda}, \cdots, \frac 1{\sigma_{m + 1} + \lambda})\begin{bmatrix}P^T g \\ P^T_\perp g \end{bmatrix} \\ = P \operatorname{diag}(\frac 1{\sigma_1 + \lambda}, \cdots, \frac 1{\sigma_m + \lambda}, \frac 1{\sigma_{m + 1} + \lambda}, \cdots, \frac 1{\sigma_{m + 1} + \lambda})P^T g + \frac 1{\sigma_{m + 1} + \lambda} (I - PP^T)g

（注意我们有

PP^T + P_\perp P_\perp^T = I

）。有了这个之后，注意到如果我们记

g_{\parallel} = P^T g

且

g_\perp = P^T_\perp g

，就可以利用补空间的正交性和勾股定理，得到

h(\lambda) = \frac 1{\sqrt{(\sum_{i = 1}^m \frac{(g_\parallel)_i^2}{(\sigma_i+\lambda)^2}) + \frac{\|g_\perp\|^2}{(\sigma_{m + 1} + \lambda)^2}}} - \frac 1 {\Delta}

这样的话，我们可以看到，因为我们所有的

\sigma_i

都是非负的，这也意味着，只要我们的

\lambda > \sigma_{min}

，这个函数就会变成一个关于

h(\lambda)

的增函数。在有增函数的条件下，想得到根是极其容易的，二分法，牛顿法都是很好的主意（这里有一个思想是求根和求函数极值在某些情况下是等价的）。

现在，我们来观察一下，我们的一个好的思路是针对信赖域方法的第二个式子，也就是

\lambda(\Delta - \|p^*\|) = 0

来找突破口。如果说我们有

\|p^*\| < \Delta

，也就是说，这个最优的解位于球的内部。那么这个时候，很容易找到

\lambda = 0

，这个时候直接令

p^* = p(0)

即可。

反过来，如果说

\|p^*\| = \Delta

，那么这个时候，如果有

h(-\sigma_{min}) < 0

，那么这个时候因为满足信赖域条件的

\lambda

是存在的，并且非负，因此根据这个论断就不难得出，一定存在一个

\lambda \ge 0

满足

h(\lambda) = 0

，这个

\lambda

就是我们要的值。如果说

h(-\sigma_{min}) = 0

，这个时候，其实相当于说

\lambda = - \sigma_{min}

就是我们要的解了，这个时候直接带入即可，但要注意

B - \sigma_{min}I

不再是一个可逆阵，换句话说要考虑广义逆的存在。也就是

p^* = -(B - \sigma_{min})^+ g

。最后如果说

h(-\sigma_{min}) > 0

，那么这个时候没辙了，我们很显然不能够再通过降低

\lambda

来缩小

h(\lambda)

了。所以没啥办法，我们也只能换一个思路。

注意看一下

h(\lambda)

的写法，如果考虑取

h(-\sigma_{min})

，那么其实，当我们的

\lambda \to -\sigma_{min}

的时候，分母中的某一项的分母会趋于0，使得分母趋于无穷大，整个式子会趋于

-\frac 1{\Delta}

。但注意到这个是在对应的分子不为0的情况下才存在的，换句话说，出现

h(-\sigma_{min})

的反常现象，其实就是

g_{\parallel}

或

g_\perp

的模长为0导致的，换句话说就是我们的

\sigma_{min}

对应的特征向量，与

这个方向正交了。这一点说明我们需要人工添加这个向量，一方面不会影响

p^*

对最小值的收敛方向（因为与梯度正交的方向移动多少，函数值都不会下降），另一方面也可以降低

h(\lambda)

的值，所以这样的话，我们的目标其实就是

p^* = -(B - \sigma_{min})^+ g + t u_{min}

，只要我们设

使得

\|p^*\| = \Delta

。

好的，我们用了非常多的文字描述了我们的结果，相信大家在理解了这些之后，就不会再对下面的算法过程感到陌生。

接下来，我们给出这个算法的合理性，也就是说我们要证明，这个算法是能够达到一定程度上的收敛性的。

Theorem 3: 设

\mathcal{N}_0 = \{x \in \mathbb{R}^n: f(x) \le f(x_0)\}

。并设

的梯度是Lipschitz连续的。也就是说对于任意的

x, y \in \mathcal{N}_0

，我们都有

\|\nabla f(x) - \nabla f(y)\| \le L \|x - y \|

。并且对于任意的

，我们都有

\frac{|(y_k - \tilde B_k s_k)^T s_k |}{\|s_k\| \|y_k - \tilde B_k s_k \|} \ge u

，其中

\mu > 0

为常数，且

B_k

为算法过程中的Hessian的估计。那么我们有

\liminf_{k \to \infty} \| \nabla f(x_k) \| = 0

要想说清楚这个证明，理解之前传统的SR1的证明自然是必不可少的。其实如果大家有观察我们的上一节的Theorem 1，我们可以看到，除去要求我们的

B_k

一致的有一个上界以外，其他的条件都是相同的，因此实际上如果我们能够说明

B_k

的这个性质，就足够完成这里的证明。

首先我们注意到，根据条件我们可以得到

\|y_k \| \le L \|s_k\|

，所以有

|\tilde \gamma_k | = | \frac{y_k^Ts_k}{s_k^Ts_k}| \le \frac {L \|s_k \|^2}{\|s_k\|^2} = L

这个东西就说明我们可以加上这样的条件，因为我们需要使用数学归纳法来解决这样的问题。同时因为我们需要说明一个一致的上界，所以对于之前的几步，也就是“限制空间”还没有填满的时候，我们是不用多过担心的（想想为什么？）。换句话说，我们可以人工的设置

\tilde B_0 = \tilde \gamma I

，然后以递推和归纳的方式证明。

首先我们来计算

\tilde B_1

，注意到

\tilde B_1 = \tilde \gamma I + \frac{(y_i - \gamma_1 s_1)(y_1 - \gamma_1 s_1)^T}{(y_1 - \gamma_1 s_1)^T s_1}

所以我们会有

\|\tilde B_1\| \le |\tilde \gamma| + \frac{\|y_1 - \tilde \gamma_1 s_1 \|^2}{|(y_1 - \tilde \gamma_1 s_1)^T s_1|} \le |\tilde \gamma| + \frac{\|y_1 - \tilde \gamma_1 s_1 \|}{\mu \|s_1\|} \le |\tilde \gamma| + \frac 1u (\frac{\|y_1\|}{\|s_1\|} + | \tilde \gamma|)

\le L + \frac 1u (L + L) = (1+\frac 1u) L + L [(1 + \frac 1u) - 1]

通过这样的变换，我们大胆的做一个推断，就是说我们认为会有

\|\tilde B_j\| \le ( 1 + \frac 1u)^j L + L [(1 + \frac 1u)^j - 1]

如果这个成立的话，要验证它也就是一步归纳法，因为我们有

\|\tilde B_{j + 1}\| \le \|\tilde B_j\| + \frac 1u (\frac{\|y_j\|}{\|s_j\|} + \| \tilde B_j\|)

（这中间其实跳了挺多步的，但是如果你看过前面的推导，你不会对这些感到困难）所以不难推导出我们的推断，换句话说，我们的所有的

\|B\|

最大也不会超过我们所要求的

\|B_{m + 1}\|

，因为我们不管做了多少步，都只会保留很少的几步，这样的操作保证了只要我们的运算中不会让模长爆破，有界的要求就可以达到。所以这样的话，证明就算完成了，因为剩下的部分都和上一节的Theorem 1一模一样了。

好的，到此我们算是介绍好了限制空间的SR1方法，而且更重要的是，我们也算正式的结束了对于无约束优化的方法的介绍。

小结

本节我们主要关注在了两个方法：限制空间的BFGS和SR1方法。而它们性质的证明都充满了浓浓的限制空间的意味。因为限制空间的特性，我们得以在大规模机器学习中看到这些方法（比方说深度学习中的优化算法，就有一个是LBFGS）。事实上这么多方法说完，相信大家对于数值优化也算有了一定的了解，不过这才只是刚刚说完无约束优化的部分，后面我们会开始进入新的方向，介绍一些新的方法与思路。

本文参与腾讯云自媒体同步曝光计划，分享自微信公众号。

原始发表：2020-07-09，如有侵权请联系 cloudcommunity@tencent.com 删除

编程算法

本文分享自学弱猹的精品小屋微信公众号，前往查看

如有侵权，请联系 cloudcommunity@tencent.com 删除。

本文参与腾讯云自媒体同步曝光计划，欢迎热爱写作的你一起参与！

编程算法

登录后参与评论

0 条评论

热度