前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >凸优化(9)——近端牛顿方法;矩阵论/数值线性代数基础:浮点数运算

凸优化(9)——近端牛顿方法;矩阵论/数值线性代数基础:浮点数运算

作者头像
学弱猹
发布2021-08-09 17:32:28
7400
发布2021-08-09 17:32:28
举报

上一节笔记:凸优化(8)——内点法中的屏障法与原始-对偶方法,近端牛顿方法

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

大家好!

这一节我们会接着上一节,介绍完近端牛顿方法(Proximal Newton Method),剩下的时间会拿来介绍一些基本的矩阵论和数值计算的知识,用于对之后介绍高阶方法的铺垫~

那么我们开始吧。

目录

  • 近端牛顿方法
  • 矩阵论/数值线性代数基础:浮点数运算
  • 矩阵论/数值线性代数基础:矩阵灵敏性分析

Source

  • CMU 10-725, Convex Optimization
  • Boyd, Vandenberghe, Convex Optimization
  • G. Golub and C. van Loan (1996), Matrix computations, Chapters 1-5, 10
  • N. Vishnoi (2013),
Lx =b

; Laplacian solvers and algorithmic applications

  • J. Friedman, T. Hastie, R. Tibshirani (2009), Regularization paths for generalized linear models via coordinate descent
  • C.J. Hsiesh, M. A. Sustik, et al. (2011), Sparse inverse covariance matrix estimation using quadratic approximation
  • J. Lee, Y. Sun, et al. (2014), Proximal Newton-type methods for minimizing composite functions
  • P. Tseng and S. Yun (2009), A coordinate gradient descent method for nonsmooth separable minimization

近端牛顿方法

在上一节我们给出了近端牛顿方法的更新公式

z^{(k)} = \operatorname{prox}_{H^{(k - 1)}}(x^{(k - 1)} - (H^{(k - 1)}) ^{-1} \nabla g(x^{(k - 1)}))
x^{(k)} = x^{(k - 1)} + t_k (z^{(k)} - x^{(k - 1)})

这里的

H^{(k- 1)} = \nabla^2 g(x^{(k - 1)})

为海塞矩阵。

还是一样,我们要先给出步长的选取条件,再给出收敛性的证明。对于步长选取条件其实没什么特别的,依然是回溯法。设

0 < \alpha \le \frac 12

0 < \beta < 1

,那么只要

f(x + tv) > f(x) + \alpha t \nabla g(x)^T v + \alpha (h(x + tv) - h(x))

那就设

t = \beta t

,否则就停止缩小,采用步长。这里

v = \operatorname{prox}_H (x - H^{-1} \nabla g(x)) - x

。这个不等式其实就是近端梯度方法里所用到的那个不等式,但是这里要注意的是,因为

v

的计算本身比较复杂,所以一般会先计算

v

再做回溯法,这样就可以防止每一次使用不等式判断的时候都要重复计算

v

现在我们来看一下它的收敛性结果。

Theorem 1: 设

mI \preceq \nabla^2 g \preceq LI

,并且

\nabla^2 g

是Lipschitz连续的,且其对应参数为

M

。假设

\operatorname{prox}_H(\cdot)

是可以精确计算的,那么它局部收敛速度为二阶收敛速度,且满足

\|x^{(k)} - x^*\|_2 \le \frac M {2m} \|x^{(k - 1)} - x^* \|_2^2

.

我们证明一下这个结论。首先注意到在

k

足够大的时候,我们会有

t = 1

,这是因为在那种情况下,可以认为函数迭代已经到了局部,即认为函数已经具备了足够的凸性,所以

t = 1

就恰好对应了

g(x)

的Taylor展开,根据

g(x)

海塞矩阵正定就可以得到这个结论,具体的细节交给读者补充。

有了

t = 1

的结论,我们就可以通过投影算子的性质给出一些证明。注意到

\|x^+ - x^* \|_2 \le \frac 1{\sqrt m} \|x^+ - x^* \|_H

这是因为

H \ge mI

,那么接下来呢,我们注意到它是通过一个近端算子得到的,也即

\frac 1{\sqrt m} \|x^+ - x^* \|_H \le \frac 1 {\sqrt m} \|\operatorname{prox}_H (x - H^{-1} \nabla g(x)) - \operatorname{prox}_H (x^* - H^{-1} \nabla g(x^*)) \|_H

虽然在这里近端算子作用之后只会得到

z

,但是因为

t = 1

,所以代入公式可以发现这个算子得到的

z

就是下一步的迭代点

x

。而到了最后一步收敛的时候,和近端梯度法一样,我们也会有

x^* = \operatorname{prox}_H (x^* - H^{-1} \nabla g(x^*))

接下来我们自然会想消去这些近端算子,这里就要用到近端算子的不可延展性(non-expansiveness,这个性质我们放到后面单独证明),也即

\frac 1 {\sqrt m} \|\operatorname{prox}_H (x - H^{-1} \nabla g(x)) - \operatorname{prox}_H (x^* - H^{-1} \nabla g(x^*)) \|_H \le \frac 1 {\sqrt{m}} \|x - x^* + H^{-1} \nabla g(x^*) - H^{-1} \nabla g(x) \|_H

下一步我们希望再做一些处理,消去这个

H^{-1}

。如果你还记得

H

表示的是海塞矩阵的话,你会明白为什么想这么做。

这也需要一点技巧,总结来说就是

\frac 1 {\sqrt m} \|y \|_H \le \frac 1{\sqrt m} \|H^{1/2} y \|_2 \le \frac 1m \|H^{1/2} y \|_H = \frac 1m \|Hy \|_2

把这个不等式代入上面,就可以得到

\frac 1 {\sqrt{m}} \|x - x^* + H^{-1} \nabla g(x^*) - H^{-1} \nabla g(x) \|_H \le \frac 1 m \| \nabla g(x^*) - \nabla g(x) + H(x - x^*) \|_2

关键的地方来了,注意到根据Taylor展开,我们可以得到

\nabla g(x) = \nabla g(x^*) + H(x- x^*) + \frac 12 (x- x^*)^T\nabla^3 g(x^*) (x -x^*)

\nabla^3 g(x^*)

是个啥?不用管它,把每一个梯度的分量对应上一个海塞矩阵,然后组合在一起就好,这个东西严格来说是一个张量(tensor),但不影响我们这里的使用)根据海塞矩阵的Lipschitz连续性我们可以得到

\frac 1 m \| \nabla g(x^*) - \nabla g(x) + H(x - x^*) \|_2 \le \frac M {2m} \|x - x^* \|_2^2

这就证明了结论。

接下来我们来证明一下之前遗留的那个性质——近端算子的不可延展性。这个性质与《数值优化》第8节(数值优化(8)——带约束优化:引入,梯度投影法)介绍的投影的不可延展性是相同的,因为投影算子只是近端算子的一个特例。

Proposition 1: 证明

\operatorname{prox}_{h, t}(x) = u

当且仅当

\forall y, h(y) \ge h(u) + \frac 1 t (x -u)^T (y - u)

,并利用此性质证明

\|\operatorname{prox}_{h, t}(x) - \operatorname{prox}_{h,t} (y)\|_2 \le \| x - y\|_2, \forall x, y

,并证明对于

\operatorname{prox}_{h, H}(\cdot)

H

为正定矩阵)也有相同的结论。

我们先证明一下第一个结论。一方面,注意到

\operatorname{prox}_{h, t}(x) =\arg\min_x \frac 1 {2t} \|x - z \|_2^2 + h(z)

,因此我们设

g(z) = \frac 1 {2t} \|x - z \|_2^2 + h(z)

,那么既然我们希望求解这个算子的精确解,我们自然需要求解次梯度(注意

h

不一定可导,所以不能使用梯度方法),也即

\partial g(z) = - \frac 1t (x - z) + \partial h(z)

那么注意到

u

是极小值点,所以有

0 \in \partial g(u)

,因此有

\frac 1t (x - u) \in \partial h(u)

。这就相当于找到了一个次梯度,因此可以利用次梯度的不等式

\forall y, h(y) \ge h(u) + g^T (y- u), g \in \partial h(u)

g = \frac 1t (x - u)

代入就可以了。

另一方面,其实反过来推就可以得到

0 \in \partial g(u)

,再结合强凸性,极小值唯一,就可以证明完毕。

有了这个结论之后,就可以来证明第二个性质了。我们只需要设

\operatorname{prox}_{h, t}(x) = u, \operatorname{prox}_{h , t}(y) = v

,那么可以得到

\forall z

,有

\forall z, h(z) \ge h(u) + \frac 1 t (x - u)^T (z - u)
\forall z, h(z) \ge h(v) + \frac 1 t (x - v)^T (z - v)

第一个式子的

z

v

,第二个式子的

z

u

,然后两个式子加在一起,就可以得到

0 \ge \frac 1 t (x - u + v - y)^T ( v -u)

因为

t > 0

所以可以消掉,然后移项即可。至于第三个结论其实非常简单,因为只需要把所有地方的

\frac 1t

改成

H

即可,过程完全相同,只需要

H

是正定的就没问题。

所以我们可以看出,近端牛顿方法其实也可以保证二阶收敛速度,但这是在近端算子可以被精确计算的情况下。如果没有的话其实就很难说,一般来说肯定近端算子能够被快速计算出来的话是最好的,比方说在

g

是一个二次函数的情况下。

虽然说它的使用其实在挺多方面是有限制的,但是也有一些它被成功应用的例子,一个是glmnet(可以见Source中的Friedman (2009)),一个是QUIC(可以见Source中的Hsiesh (2011))。具体的细节我们这里就不详谈了。

最后关于这个方法,我们再提一下Source中的Lee (2014)所给出的一些数值实验。下面这两张图展现了不同方法的对graphical lasso问题的比较结果。

可以看出来,在一般的情况下结果还是不错的。但是,在近端算子计算不精确的情况下,结果就不会那么好。下面这张图就说明了这一点。

重点看的是绿色的线,它其实说明了在计算近端算子的时候,如果设置了最大迭代步数为10步,是不足以使得迭代具备二次收敛速度的。这也说明了其实一般情况下,对于近端算子的计算的要求很高,也一定程度上说明了近端牛顿法的限制性。当然了这篇paper里也给出了一些其他的结果,比方说红色的线就是他们自己的方法,他们也给出了局部超线性收敛速度的结果。

好的,关于近端牛顿方法,我们就提这么多。

矩阵论/数值线性代数基础:浮点数运算

从这一部分开始,我们即进入了10-725这一门课的最后一部分。这一部分我们会关注一些更贴近人工智能,理论更不完备,思想也更加深刻的一些方法。不过在介绍它们之前,我们先打一些矩阵论/数值线性代数的基础,虽然并非之后的所有方法都需要依赖这一部分的工具,但不可否认的是,很多对于方法好坏和适用场景的分析都是依赖了矩阵论的工具,因此了解这一部分内容,对于优化这个学科而言自然是百利无害。

首先我们要关注的是我们之前一直在介绍,但是从来没有量化过的一个概念:时间复杂度计算。注意,这不是之前给出的一些迭代步数的分析结果(虽然它们也是时间复杂度分析的很重要的部分),而是与矩阵计算相关的一些时间复杂度。比方说牛顿法虽然具有局部上的二次收敛速度,但是它同时需要依赖海塞矩阵的计算,这就意味着即使迭代步数上牛顿法显著优于梯度法,但每一步迭代所花费的时间也急剧增大了。这一部分分析的主要方面也即是“每一步迭代所花费的时间”。

那么我们开始吧。首先给出的是浮点数运算(Floating Point Operation)的概念。

Definition 1: Floating point operation 定义一次浮点数运算为算法中的一次加减乘除运算。

一般来说我们用flop作为浮点数运算的单位,浮点数运算是非常重要的影响算法运行时间的部分。比方说计算

a + b

对于

a, b \in \mathbb R^{n}

就会进行

n

次运算。即

n

次加法,而计算

a^T b

则需要

2n

次运算,即计算完每一个

a_i b_i

n

次乘法)之后,还需要把它们都加起来(

n

次加法),当然有的人会说底层实现的时候,可以通过一开始赋值一个

a_1b_1

,这样加法就只需要

n-1

次了。emmmm这当然没有问题,不过浮点数运算次数的计算我们更多关心量级,不会太苛求于具体的次数。比方说对于

a+ b

a^Tb

,如果写成大O表达式,其实都是

O(n)

的时间复杂度(不过对于学计算机科学的同学来说,优化这些量级前的常数可能就是他们毕生的工作233

注意这只是涵盖了计算时间的一部分,比方说赋值就不会算作一次运算,但是这不代表它在计算机内运行的时候不消耗时间。但是对于分析来说这足够了。

接下来我们来看一些例子。

Example 1: Matrix-vector product 设

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

,

b \in \mathbb R^{n}

,考虑计算

Ab

对于一般情况来说,我们可以考虑设

A = \begin{bmatrix}a_1^T \\ \vdots \\ a_m^T \end{bmatrix}

,那么会有

Ab = \begin{bmatrix}a_1^T b \\ \vdots \\ a_m^T b\end{bmatrix}

,而对于每一个

a_i^T b

都需要

2n

次浮点数运算,所以对于

Ab

,需要

2mn

次浮点数运算(注意把它们都赋值到矩阵的对应位置是不计入浮点数运算的,也不在我们这里的分析范围内)。

但是这只是一般的情况,如果矩阵具备一些特征,那么浮点数运算的次数就会大大减少。比方说如果矩阵只有

s

个非零元素,那么注意,如果某个元素是

0

,计算机是不会计算的,因此对于这种情况,浮点数运算次数就是

2s

。同样的思路,如果矩阵是一个

k

带矩阵(k-banded),那么一个很重要的特征就是它每一行只有

k

个非零元素,所以浮点数运算次数就是

2nk

有一个稍微复杂点的情况是当矩阵可以分解为

A = \sum_{i = 1}^r u_i v_i^T

的浮点数运算次数,我们可以看一下这个时候会有

Ab = \sum_{i = 1}^r u_i v_i^Tb

,所以每一个

v_i^T b

需要

2n

次浮点数运算,一共是

2nr

次运算。每一个常数计算好之后,乘上

u_i

需要

m

次乘法(因为

u_i

是向量),一共就是

mr

次运算,最后把每一个

u_i v_i^T b

算好加在一起,就需要

mr

次加法。所以一共就是

2r(m + n)

次浮点数运算。这个结果说明了一般来说,低秩矩阵的运算次数会显著的变少。因为这里的

r

其实就是矩阵做分解之后所对应的秩。

最后还有一种情况就是矩阵为置换矩阵,这种情况因为实际上,

Ab

的作用就是对

b

的各行做了一个交换,不涉及到加减乘除的运算,所以这个运算次数是0。

Example 2: Matrix-matrix product 设

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

B \in \mathbb R^{n \times p}

,考虑计算

AB

还是一样,可以考虑设

B = \begin{bmatrix}b_1 & \ldots & b_p\end{bmatrix}

,所以可以得到

AB = \begin{bmatrix}Ab_1, \ldots, Ab_p\end{bmatrix}

,那么一般情况下,求解

Ab_i

需要

2mn

次浮点运算,那么一共有

p

列,所以对应的是

2mnp

次运算。这也是为什么我们一般说矩阵乘法的时间复杂度是

O(n^3)

,这其实就是

m = n = p

的一个特殊情况。

有一个情况就是

A

只有

s

个元素非零的情况,那么因为这个时候,每一个

Ab_i

需要

2s

次浮点运算,所以一共就是

2sp

次浮点数运算。

Example 3: Matrix-matrix-vector product 设

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

B \in \mathbb R^{n \times p}

c \in \mathbb R^p

,考虑计算

ABc

矩阵运算虽然没有交换律,但是具有结合律,因此按什么顺序计算就是一个很重要的问题。假如说我们按照

A(Bc)

的顺序计算,那么计算

Bc

需要

2np

次,然后得到一个计算好的新的向量

Bc

再和

A

相乘,因为

A

m

行,

Bc

n

行,那么相乘需要

2mn

次,所以一共是

2mn + 2np

次。

但是如果按照

(AB)c

的方式计算,那么先计算

AB

需要

2mnp

次,然后得到一个新的矩阵再和

c

相乘就需要

2mp

次,所以一共需要

2mnp+2mp

次,很明显这样计算的时间复杂度就大大增加了。

Example 4: Solving linear systems 设

A \in \mathbb R^{n \times n}

,考虑求解线性方程组

Ax = b

对于不同结构的

A

,我们会分析出非常不同的结果。

如果

A

对角阵

\operatorname{diag}(a_1, \ldots, a_n)

,那么容易得到解

x = (b_1/a_1, \cdots, b_n / a_n)^T

,所以需要

n

次除法运算,因此这种情况下求解方程组的时间复杂度就是

O(n)

如果

A

上三角阵或下三角阵(比如说上三角阵的意思就是

A_{ij} = 0, i < j

,当然了

A_{ij}

就表示矩阵的第

i

行第

j

列元素 ),那么这个时候是可以通过追赶法(当然这个方法严格来说是用来求解三对角矩阵的)来进行求解的,并且可以显式的写出解为

x_1 = b_1 / A_{11}
x_2 = (b_2 - A_{21} x_1) / A_{22}
\vdots
x_n = (b_n - A_{n, n- 1}x_{n - 1} - \cdots - A_{n1}x_1) / A_{nn}

计算每一个

x_i

的浮点数运算次数当然没问题,但是也没啥必要,因为容易看出每一个

x_i

大约需要

2i

步运算,所以一共可以算出大约是

n^2

次的浮点数运算次数。所以对于这种情况,我们的解方程组的时间复杂度大约是

O(n^2)

还有一个有趣的情况是

A

正交矩阵,这种情况下我们有

A^T = A^{-1}

,所以可以得到

x = A^T b

,那么这种情况就化归为了Example 1,所以这个浮点数运算次数是

2n^2

,依然是

O(n^2)

的时间复杂度。

当然了,最多的一种情况就是矩阵没有任何特殊结构,这种情况求解

x

,常见的情况就是高斯消元。可以思考一下,从把一般的矩阵消成上三角矩阵,就需要大约

\sum_{i = 1}^n 2i^2

次浮点数运算,因此对于这种情况,求解线性方程组就需要

O(n^3)

的时间复杂度。所以很多时候会说求解一个矩阵的逆需要

O(n^3)

的时间复杂度,是因为这个问题和求解

Ax = I

是等价的。

这个问题的一个最直接的应用就是在矩阵分解上,很多时候可以通过矩阵分解使得线性方程组的求解变得简单。具体来说,设

A = A_1 A_2 \ldots, A_k

,那么这个时候,可以得到

x = A_k^{-1} \ldots, A_2^{-1}A_1^{-1} b

,当然了求解

x

不可能是通过先求解

k

个矩阵的逆矩阵,再乘起来的这种操作,这是编程一开始最容易犯的错误。正确的做法应该是先通过解方程组来求解

y = A_1^{-1} b

,也就是求解

A_1 y = b

,再求解

z = A_2^{-1}y

,以此类推。这样做比较好的原因在于,如果每一个小的矩阵都有特殊的结构,那么每一步小的解方程组的运算就相比较之前来说会快很多。

当然了有的人会说分解本身就会需要很长时间。这当然是没错的,矩阵分解本身也是

O(n^3)

的时间复杂度的。但是如果说我们做了一次分解,可以给多次计算使用,那么就会大大降低时间复杂度。举个例子,假如说我们希望计算

AX = B

(注意这里的

X, B

都是矩阵了),那么就可以把它拆成一系列线性方程组

A x_i = b_i

那么实际上在解的时候,我们可以做一次分解,然后每一个方程组都使用分解后的

A

来求解,这样的话如果有

t

个这样的线性方程组,那么时间复杂度也只是

n^3 + tn^2

而已,但如果不做分解,时间复杂度就是

tn^3

,差距立现。

比方说任何一个矩阵都具备一个奇异值分解

A = U \Sigma V^T

,这里会发现

U, V^T

都是正交矩阵(一个是列正交,一个是行正交,注意区分),而

\Sigma

是对角矩阵,因此它们都算有特殊结构的矩阵,在求解的时候就会有比较大的加速。

比较常见的分解就是QR分解Cholesky分解,我们简单介绍一下。

Definition 2: QR Decomposition 设

A \in \mathbb R^{m \times n}, m \ge n

,那么存在矩阵分解

A = QR

,且

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

是列正交的,

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

是上三角阵。计算此分解需要

2mn^2 - \frac {n^3} 3

次浮点数运算。 Definition 3: Cholesky Decomposition 对任意的

A \in \mathbb S_{++}^n

,都存在一个Cholesky分解

A = LL^T

L \in \mathbb R^{n \times n}

是一个下三角矩阵。计算此分解需要

\frac {n^3} 3

次浮点数运算。

我们挑一个比较好算的Cholesky分解来看一下它的浮点数运算次数是如何计算的(当然了对于QR分解,准确的计算结果可能有些困难,但得到

2mn^2

的这个量级是不困难的,抓住主要矛盾就好)。

对于Cholesky分解,只需要注意到它的更新公式

l_{jj} = (a_{jj} - \sum_{k = 1}^{j - 1} l_{jk}^2)^{\frac 12}
l_{ij} = (a_{ij} - \sum_{k = 1}^{j - 1} l_{ik} l_{jk}) / l_{jj}, i = j+ 1, \cdots, n

(这个公式怎么推是数值分析的内容,这里就不多说了)那么注意,以

j

作为支点,我们能看出来,固定

j

之后,每一个

i

j

n

共有

n - j+ 1

项,每一项的浮点数运算次数大约为

2j

次,所以可以给出一个大致的值

\sum_{j = 1}^n 2 j (n - j + 1) = 2n \cdot \frac {n^2 } 2 - 2 \cdot \frac {n^3} 3 = \frac{n^3} 3

这就是它的浮点数运算次数。

接下来我们看一个实际的例子

Example 5 多元线性回归的解的表达式为

\hat \beta = (X^T X)^{-1} X^T y

,设

X \in \mathbb {R}^{n \times p}

,给出计算它的一种方案并给出它的浮点数运算次数。

首先我们要计算

X^T y

,这需要

2np

次,然后我们需要计算

X^TX

(显然我们不可能直接计算矩阵的逆),这需要

np^2

次(有的人可能会说这是

2np^2

次,但是注意到它是一个对称矩阵,所以只需要算一个矩阵的一半),接下来为了解方程组

(X^TX) \beta = X^T y

,我们需要对

X^TX

做矩阵分解,使用Cholesky分解需要

\frac {p^3} 3

次,最后再解一下方程组需要

2p^2

次(注意这是因为如果

X^T X = LL^T

,那么一个再解方程组其实就是两个

p^2

的操作)。所以加在一起,去除掉一些次要矛盾,就可以得到浮点数运算次数大致为

np^2 + \frac {p^3} 3

,如果

p

很小的话,那么

\frac {p^3} 3

其实也可以忽略不计。

当然了,如果分解使用QR分解的话,那么这个运算次数就大约为

3np^2

,实际上是不如Cholesky分解划算的。

矩阵论/数值线性代数基础:矩阵灵敏性分析

这里的灵敏性分析虽然也叫sensitivity analysis,但是和优化问题的灵敏性分析不太一样(优化问题的灵敏性分析会放到后面单独说),我们这里是研究解线性方程组的灵敏性,是数值分析的重要的内容。

问题的大背景是这样的:对于解线性方程组

A x = b

如果对

A, b

施加一些扰动,对于解会有什么影响?更具体点来说,就是现在,设

F \in \mathbb R^{n \times n}

f \in \mathbb R^n

,考虑问题

(A + \epsilon F)x(\epsilon) = (b + \epsilon f)

那么对应的新的解

x(\epsilon)

相比较

x

来说,会差距多少?这可以由下面这个性质所解答。

Proposition 2 记

x = x(0)

,那么有

\frac{\| x(\epsilon) - x \|_2^2}{\|x \|_2} \le \kappa (A) (\rho_A + \rho_b) + O(\epsilon^2)

,其中

\kappa(A)

是矩阵

A

的条件数

\frac{\sigma_1}{\sigma_n}

,且

\rho_A = |\epsilon| \|F\|_{2} / \|A \|_{2}

\rho_b = |\epsilon| \|f \|_2 / \|b \|_2

为相对误差。

我们证明一下这个结论。思路在于既然我们的解

x

是一个关于

\epsilon

的函数,那么不如就考虑它的一个变化率,毕竟我们的不等式左边也确实很像变化率。

既然需要得到变化率,我们就考虑对等式两边对

\epsilon

求微分,这样的话就会有

A x'(\epsilon) + F x(\epsilon) + \epsilon F x'(\epsilon) = f

那么这样的话,设

\epsilon = 0

就可以得到

x'(0) = A^{-1} (f - Fx)

,那么这样的话就可以根据Taylor展开得到

x(\epsilon) = x + \epsilon A^{-1} (f - Fx) + O(\epsilon^2)

所以这样的话就可以得到

\frac{\|x (\epsilon) - x \|_2}{\|x \|_2} \le |\epsilon| \|A^{-1} \|_2 (\frac{\|f\|_2}{\|x\|_2} + \frac{\|Fx\|_2}{\|x\|_2} ) + O(\epsilon^2)

这是通过一步三角不等式得到的。

最后我我们再做一步,就是两边先乘上一个

\|A\|_2

,再除以一个

\|A\|_2

,但区别在于第一次是为了构造出

\kappa(A)

这个量,第二次则是希望利用

\|A\|_2 \|x\|_2 \ge \| b\|_2

来修改分母。这些都做完之后,再利用矩阵范数的定义

\|F \|_2 = \max_{x} \frac{\|F x \|_2}{\|x \|_2}

,右边再取一个max就可以得到结论了。

这个性质告诉我们,如果线性方程组系数的条件数过大,那么解的变化率就会越大,因此在解线性方程组的时候也会具有更多的不稳定性。

一个比较常见的例子就是,有些时候在求解线性方程组的时候,会更多的希望对矩阵做QR分解而不是Cholesky分解。我们在《数值优化》第B节(数值优化(B)——二次规划(上):Schur补方法,零空间法,激活集方法)的零空间方法中提到过这个现象。

为什么呢?这也是因为做Cholesky分解对应的矩阵和QR分解对应的矩阵不同,而一般来说做Cholesky分解对应的矩阵条件数也会更大,所以稳定性会更差一些。当然Cholesky分解本身还有一个不稳定的问题,这个可以参考《数值优化》第C节(数值优化(C)——二次规划(下):内点法;现代优化:罚项法,ALM,ADMM;习题课)的习题课的部分。

好的,关于矩阵论的补充知识,我们就说这么多。

小结

这一节在算法层面上没有什么太多的干货,主要是对近端牛顿方法做了一个简介,更多的则是在做一些矩阵论知识的补充,它们也是优化的重要基础知识之一,也是很多时候用来分析复杂方法的工具。相信看完这一节之后,再回头看《数值优化》中对算法复杂度的分析,就会更加明晰和得心应手一些。

下一节我们会回归初心,继续以算法为主流介绍相关的内容。

本文参与 腾讯云自媒体分享计划,分享自微信公众号。
原始发表:2020-10-15,如有侵权请联系 cloudcommunity@tencent.com 删除

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

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

本文参与 腾讯云自媒体分享计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 目录
  • Source
  • 近端牛顿方法
  • 矩阵论/数值线性代数基础:浮点数运算
  • 矩阵论/数值线性代数基础:矩阵灵敏性分析
  • 小结
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档