第二章: 无约束优化基础
文章目录
在无约束优化中, 我们试图在对变量 x x x上不加任何约束时, 寻求目标函数 f f f的最小值点, 即
min x f ( x ) , \min_xf(x), xminf(x),
其中 x ∈ R n , f : R n → R x\in\mathbb{R}^n, f:\mathbb{R}^n\to\mathbb{R} x∈Rn,f:Rn→R为一光滑函数. 基于目的, 我们需要了解:
- 解是什么? 判定条件如何?
- 根据判定条件, 有何高效的算法?
而这也是本章的架构.
1. 关于解
1.1 解的定义
- 全局极小点: x ∗ x^* x∗为全局极小点, 若 f ( x ∗ ) ≤ f ( x ) , ∀ x f(x^*)\le f(x),\quad\forall x f(x∗)≤f(x),∀x.
- (弱)局部极小点: x ∗ x^* x∗为局部极小点, 若 ∃ x ∗ \exists x^* ∃x∗的邻域 N \mathcal{N} N, 使得 f ( x ∗ ) ≤ f ( x ) , ∀ x ∈ N f(x^*)\le f(x),\quad\forall x\in\mathcal{N} f(x∗)≤f(x),∀x∈N.
- 严格(或强)局部极小点: x ∗ x^* x∗为严格局部极小点, 若 ∃ x ∗ \exists x^* ∃x∗的邻域 N \mathcal{N} N, 使得 f ( x ∗ ) < f ( x ) , ∀ x ∈ N ∖ { x ∗ } f(x^*)<f(x),\quad\forall x\in\mathcal{N}\setminus\{x^*\} f(x∗)<f(x),∀x∈N∖{x∗}.
- 孤立局部极小点: x ∗ x^* x∗为孤立局部极小点, 若 ∃ x ∗ \exists x^* ∃x∗的邻域 N \mathcal{N} N, 使得 x ∗ x^* x∗为 N \mathcal{N} N中唯一的局部极小点.
例 f ( x ) = x 4 cos ( 1 / x ) + 2 x 4 , f ( 0 ) = 0 f(x)=x^4\cos(1/x)+2x^4, f(0)=0 f(x)=x4cos(1/x)+2x4,f(0)=0. 显然有严格局部极小点 x ∗ = 0 x^*=0 x∗=0, 但它附近还有无穷多不同于它的局部极小点 x j x_j xj, 且 x j → 0 , j → ∞ x_j\to0,j\to\infty xj→0,j→∞.
从上例可以看出, 严格局部极小点不一定是孤立局部极小点, 但反之是对的.
一般来说, 我们不苛求一个算法能够找到全局极小点, 因为往往在迭代的过程中陷入局部极小 (尤其是那种到处是坑的函数, 可见’墨西哥草帽’). 因此我们往往追求和识别局部极小. 当然不可否认, 在特殊情形下, 全局极小点是显而易见的, 例如凸优化.
1.2 局部极小的判别
判别局部极小的主要数学工具是Taylor定理 (或Taylor展开).
定理1 (Taylor定理) 设
f
:
R
n
→
R
f:\mathbb{R}^n\to\mathbb{R}
f:Rn→R连续可微, 有一点
p
∈
R
n
p\in\mathbb{R}^n
p∈Rn. 于是我们有
f
(
x
+
p
)
=
f
(
x
)
+
∇
f
(
x
+
t
p
)
T
p
,
t
∈
(
0
,
1
)
.
f(x+p)=f(x)+\nabla f(x+tp)^Tp, \quad t\in(0,1).
f(x+p)=f(x)+∇f(x+tp)Tp,t∈(0,1).
进一步, 若
f
f
f二阶连续可微, 则
∇
f
(
x
+
p
)
=
∇
f
(
x
)
+
∫
0
1
∇
2
f
(
x
+
t
p
)
p
 
d
t
,
\nabla f(x+p)=\nabla f(x)+\int_0^1\nabla^2f(x+tp)p\,\mathrm{d}t,
∇f(x+p)=∇f(x)+∫01∇2f(x+tp)pdt,
从而
f
(
x
+
p
)
=
f
(
x
)
+
∇
f
(
x
)
T
p
+
1
2
p
T
∇
2
f
(
x
+
t
p
)
p
,
t
∈
(
0
,
1
)
.
f(x+p)=f(x)+\nabla f(x)^Tp+\frac{1}{2}p^T\nabla^2f(x+tp)p,\quad t\in(0,1).
f(x+p)=f(x)+∇f(x)Tp+21pT∇2f(x+tp)p,t∈(0,1).
为什么不展开至更高的项呢? 对于多元函数, 展开至更高阶时会牵涉张量, 而对于曲面分类 (
f
f
f可以看成是高维直角坐标系上的曲面), 高等代数所研究的限于二次型, 也就是展开至二阶的情形. 至于为什么要说曲面分类, 我们知道高等代数根据标准型将二次型分为几大类, 对于函数曲面来说, 展开至二阶的多项式可以视作
f
f
f较好的局部逼近, 从而
- 若展开式二阶项的矩阵为正定或负定矩阵, 则说明局部 f f f近似为椭圆抛物面;
- 若展开式二阶项的矩阵为不定矩阵, 则说明局部 f f f近似为双曲抛物面 (即马鞍面);
- 若展开式二阶项的矩阵为半定矩阵, 则说明局部 f f f近似为抛物柱面.
因此, 如果 x ∗ x^* x∗是局部极小点, 那么它周围应当类似是椭圆抛物面, 它位于顶点. 这样, 我们有以下判定局部极小的必要条件和充分条件.
定理2 (一阶必要条件) 若
x
∗
x^*
x∗为局部极小点且
f
f
f在
x
∗
x^*
x∗的一个开邻域内连续可微, 则
∇
f
(
x
∗
)
=
0
\nabla f(x^*)=0
∇f(x∗)=0. 称满足
∇
f
(
x
)
=
0
\nabla f(x)=0
∇f(x)=0的点为稳定点.
定理3 (二阶必要条件) 若
x
∗
x^*
x∗为局部极小点且
∇
2
f
\nabla^2 f
∇2f在
x
∗
x^*
x∗的一个开邻域内存在且连续, 则
∇
f
(
x
∗
)
=
0
\nabla f(x^*)=0
∇f(x∗)=0且
∇
2
f
(
x
∗
)
\nabla^2 f(x^*)
∇2f(x∗)半正定.
定理4 (二阶充分条件) 设
∇
2
f
\nabla^2 f
∇2f在
x
∗
x^*
x∗的一个开邻域内连续且
∇
f
(
x
∗
)
=
0
\nabla f(x^*)=0
∇f(x∗)=0,
∇
2
f
(
x
∗
)
\nabla^2 f(x^*)
∇2f(x∗)正定, 则
x
∗
x^*
x∗是
f
f
f的一个严格局部极小点.
我们可以举出例子说明以上条件为何仅仅是必要或充分. 对于必要条件, 考虑 f ( x ) = x 3 f(x)=x^3 f(x)=x3, 在 x = 0 x=0 x=0处满足导数和二阶导数均为0, 但它不是局部极小点 (鞍点, saddle point). 对于充分条件, 考虑 f ( x ) = x 4 f(x)=x^4 f(x)=x4, x = 0 x=0 x=0是严格局部极小点, 但二阶导为0.
对于凸函数, 我们有
定理5 当
f
f
f是凸函数, 任一局部极小点
x
∗
x^*
x∗都是全局极小点. 若
f
f
f还可微, 则任一稳定点都是全局极小点.
若碰上了非光滑问题 (即目标函数不光滑), 现常用的手段是次梯度法 (subgradient methods) , 这里不做讨论. 特别的非光滑问题, 如
f
(
x
)
=
∥
r
(
x
)
∥
1
,
f
(
x
)
=
∥
r
(
x
)
∥
∞
f(x)=\Vert r(x)\Vert_1,\quad f(x)=\Vert r(x)\Vert_{\infty}
f(x)=∥r(x)∥1,f(x)=∥r(x)∥∞可以化成光滑的约束优化问题. 这在后面会提到.
2. 无约束优化算法综述
这里只讨论光滑函数情形. 算法流程可见第一章的博客. 对于无约束优化问题, 常用以下两种策略.
2.1 线搜索 (Line Search)
线搜索策略中, 算法在生成下一个迭代点前将先选取合适的搜索方向, 再在当前搜索方向上, 寻找具有更小函数值的点作为下一迭代点. 由此带来两个问题:
- 何为合适的搜索方向?
- 何为合适的搜索距离? 我们当然希望能够精确求解以下问题: 设 p k p_k pk为当前的搜索方向, min α > 0 f ( x k + α p k ) . \min_{\alpha>0}f(x_k+\alpha p_k). α>0minf(xk+αpk). 当精确求解消耗过大时, 我们应当考虑有限的求解方案. 它需要满足一定的下降等条件. 这一点我们放到下一章讨论.
2.1.1 线搜索的搜索方向
我们选取搜索方向的一项重要准则就是, 这一方向要能够为我们带来目标函数值的下降. 从微积分中我们知道, 函数的梯度为它上升最快的方向, 从而负梯度就是它下降最快的方向. 进一步地, 任一与负梯度夹角小于
π
2
\frac{\pi}{2}
2π的方向都是函数值下降方向. 对于这些能够带来下降的方向, 我们也并不是照单全收. 事实上, 根据下一章的Zoutendijk条件, 我们需要
cos
θ
k
≥
δ
>
0
\cos\theta_k\ge\delta>0
cosθk≥δ>0, 其中
θ
k
\theta_k
θk是当前搜索方向与负梯度的夹角, 即
θ
k
=
∠
(
p
k
,
−
∇
f
(
x
k
)
)
.
\theta_k=\angle(p_k,-\nabla f(x_k)).
θk=∠(pk,−∇f(xk)).
下面是一些常用的搜索方向.
- 最速下降方向. 也即负梯度方向. 这一方向的选取是最直观的, 它能够带来当前迭代点下最大的函数值下降 (当然未必是最好的选择, 正如贪心策略不一定最好一样). 最速下降法在问题后期可能会受制于锯齿效应 (Zigzag), 这在变量量纲不一致时尤其明显.
- 牛顿方向. 此方向源于
f
(
x
k
+
p
)
f(x_k+p)
f(xk+p)的二阶Taylor展开,
f
(
x
k
+
p
)
≈
f
k
+
p
T
∇
f
k
+
1
2
p
T
∇
2
f
k
p
≜
m
k
(
p
)
.
f(x_k+p)\approx f_k+p^T\nabla f_k+\frac{1}{2}p^T\nabla^2 f_kp\triangleq m_k(p).
f(xk+p)≈fk+pT∇fk+21pT∇2fkp≜mk(p). 假设
∇
2
f
k
\nabla^2 f_k
∇2fk是正定矩阵, 则牛顿方向为
p
k
N
=
arg
min
p
m
k
(
p
)
=
−
(
∇
2
f
k
)
−
1
∇
f
k
.
p_k^N=\arg\min_{p}m_k(p)=-(\nabla^2 f_k)^{-1}\nabla f_k.
pkN=argpminmk(p)=−(∇2fk)−1∇fk.当
f
f
f与二次模型
m
k
m_k
mk相差不大时, 牛顿方向是很好的选择. 下面验证它在
∇
2
f
k
\nabla^2 f_k
∇2fk是正定时是
f
f
f的下降方向.
∇
f
k
T
p
k
N
=
−
(
p
k
N
)
T
∇
2
f
k
p
k
N
<
0.
\nabla f_k^Tp_k^N=-(p_k^N)^T\nabla^2 f_kp_k^N<0.
∇fkTpkN=−(pkN)T∇2fkpkN<0.当
∇
2
f
k
\nabla^2 f_k
∇2fk不是正定矩阵时, 牛顿方向可能甚至无法定义 (如半定时). 即使可以确定 (如不定矩阵), 也无法保证牛顿方向仍然是下降方向. 此时应当对牛顿方向加修正, 使它成为下降方向同时保留
∇
2
f
k
\nabla^2 f_k
∇2fk中的二阶信息. 我们将这些放在下一章中讨论.
需要说明是, 牛顿方向在应用时具有天然的步长, 即 α k = 1 \alpha_k=1 αk=1. 这就省去了步长选取时的线搜索; 牛顿方向需要 f f f的二阶信息才可计算, 但事实上二阶信息往往难以获取. 一种替代方案是利用后面会提到的有限差分逼近, 此时只需要一阶信息. 另一种替代方案则是下面要提到的拟牛顿方向. - 拟牛顿 (Quasi-Newton)方向. 计算拟牛顿方向时我们利用的模型与计算牛顿方向类似, 只是将
m
k
(
p
)
m_k(p)
mk(p)中的
∇
2
f
k
\nabla^2 f_k
∇2fk替换成它的近似矩阵
B
k
B_k
Bk. 为了不在每步都重新计算近似矩阵, 拟牛顿法将在每步迭代对
B
k
B_k
Bk做一次更新. 更新所遵循的法则有
- 相邻两 B k B_k Bk的差不"大";
- B k B_k Bk是对称矩阵;
-
B
k
B_k
Bk满足割线方程 (The secant equation):
B
k
+
1
s
k
=
y
k
B_{k+1}s_k=y_k
Bk+1sk=yk, 其中
s
k
=
x
k
+
1
−
x
k
,
y
k
=
∇
f
k
+
1
−
∇
f
k
s_k=x_{k+1}-x_k,y_k=\nabla f_{k+1}-\nabla f_k
sk=xk+1−xk,yk=∇fk+1−∇fk.
其中第1条放在极小点 x ∗ x^* x∗的附近是无可厚非的, 我们自然可以要求在算法后期二阶信息的变动不大, 这一点可以用"特定的"矩阵范数度量; 第2条, 由于 B k B_k Bk是二阶信息的近似, 这一点提法也对; 第3条, 则要利用梯度在 p k p_k pk上的变化反映二阶信息. 根据Taylor定理, 有 ∇ f ( x + p ) = ∇ f ( x ) + ∇ 2 f ( x ) p + ∫ 0 1 [ ∇ 2 f ( x + t p ) − ∇ 2 f ( x ) ] p   d t . \nabla f(x+p)=\nabla f(x)+\nabla^2 f(x)p+\int_0^1[\nabla^2 f(x+tp)-\nabla^2 f(x)]p\,\mathrm{d}t. ∇f(x+p)=∇f(x)+∇2f(x)p+∫01[∇2f(x+tp)−∇2f(x)]pdt. 其中积分由二阶量的连续性, 就是 o ( ∥ p ∥ ) o(\Vert p\Vert) o(∥p∥), 再令 p = α k p k , x = x k p=\alpha_kp_k,x=x_k p=αkpk,x=xk, 我们就有 ∇ f k + 1 = ∇ f k + ∇ 2 f k ( x k + 1 − x k ) + o ( ∥ x k + 1 − x k ∥ ) . \nabla f_{k+1}=\nabla f_k+\nabla^2 f_k(x_{k+1}-x_k)+o(\Vert x_{k+1}-x_k\Vert). ∇fk+1=∇fk+∇2fk(xk+1−xk)+o(∥xk+1−xk∥).略去最后的高阶量, ∇ 2 f k ( x k + 1 − x k ) ≈ ∇ f k + 1 − ∇ k . \nabla^2f_k(x_{k+1}-x_k)\approx\nabla f_{k+1}-\nabla_k. ∇2fk(xk+1−xk)≈∇fk+1−∇k.再由第1条, 得到割线方程. 综上所述, 计算 B k + 1 B_{k+1} Bk+1实质上就是求解以下约束优化问题: min B ∥ B − B k ∥ s u b j e c t    t o B = B T , B s k = y k . \begin{array}{rll}\min_B & \Vert B-B_k\Vert & \\\mathrm{subject\,\,to} & B=B^T, & Bs_k=y_k.\end{array} minBsubjectto∥B−Bk∥B=BT,Bsk=yk.后面我们会推导上面问题的求解过程, 这里先给出更新公式. B k + 1 = ( I − ρ k y k s k T ) B k ( I − ρ k s k y k T ) + ρ k y k y k T , ρ k = 1 y k T s k . B_{k+1}=(I-\rho_ky_ks_k^T)B_k(I-\rho_ks_ky_k^T)+\rho_ky_ky_k^T,\quad\rho_k=\frac{1}{y_k^Ts_k}. Bk+1=(I−ρkykskT)Bk(I−ρkskykT)+ρkykykT,ρk=ykTsk1.这一公式称为Davidon-Fletcher-Powell(DFP)公式. 注意到若将上述问题改为对 B k − 1 B_{k}^{-1} Bk−1直接更新 (这样一来省去了求解线性方程组 B k p = − ∇ f k B_kp=-\nabla f_k Bkp=−∇fk的工作, 而只需要一步矩阵-向量乘积, 计算量从 O ( n 3 ) O(n^3) O(n3)降至 O ( n 2 ) O(n^2) O(n2), 当然也可通过直接更新 B k B_k Bk的Cholesky因子 L k L_k Lk, 此时计算量也在 O ( n 2 ) O(n^2) O(n2)), 我们就会得到Broyden-Fletcher-Goldfarb-Shanno(BFGS)公式.1这样一来, p k = − B k − 1 ∇ f k p_k=-B_k^{-1}\nabla f_k pk=−Bk−1∇fk. 在 B k B_k Bk正定时 (这只需要最初 B 0 B_0 B0是正定的), 我们也可以保证 p k p_k pk一定是下降方向.
- 共轭方向. 它来源于非线性共轭梯度法, 具有形式
p
k
=
−
∇
f
k
+
β
k
p
k
−
1
,
p_k=-\nabla f_k+\beta_kp_{k-1},
pk=−∇fk+βkpk−1,其中
β
k
\beta_k
βk保证了
p
k
p_k
pk和
p
k
−
1
p_{k-1}
pk−1是"共轭的"2. 从
p
k
p_k
pk的形式中我们得到了两个信息:
- p k p_k pk是负梯度方向和之前搜索方向的线性组合. 相对于最速下降方向, 共轭方向在二维平面上达到最小. 相当于最速下降急着下山, 一条道走到黑; 而共轭梯度法会环顾四周, 得到更佳的方向. 这个方向可能在当前不是最好的, 但最后表现可能会更佳. p k p_k pk是下降方向则是不言而喻的.
- 在计算 p k p_k pk时, 只需一次一阶信息 ∇ f k \nabla f_k ∇fk的获取、一个向量 p k − 1 p_{k-1} pk−1的储存以及一个系数 β k \beta_k βk的计算 (实际上后面可以看到, β k \beta_k βk的计算也至多牵涉一次矩阵-向量乘积). 这说明共轭梯度法很适合于大规模问题的求解.
注意我们暂时不牵涉有关全局收敛和收敛速度的讨论. 但给予一定条件下, 关于收敛速度我们有如下比较: 牛 顿 方 向 > 拟 牛 顿 法 > 共 轭 梯 度 法 ≈ 最 速 下 降 法 牛顿方向>拟牛顿法>共轭梯度法\approx最速下降法 牛顿方向>拟牛顿法>共轭梯度法≈最速下降法当然正如在第一章提到的, 算法没有绝对的好与坏. 收敛速度不是唯一判定算法适合与否的标准.
2.2 信赖域方法 (Trust-Region Methods)
信赖域方法使用的模型 (这里假定) 仍然是
m
k
(
p
)
m_k(p)
mk(p), 只是提出的问题有所不同. 信赖域考虑的是,
m
k
(
p
)
m_k(p)
mk(p)在
x
k
x_k
xk附近可能是
f
f
f的一个较好近似, 但离得远就不一定了. 因此需要在求解方向时, 限定求解范围在
x
k
x_k
xk的一个邻域内. 若这一步迭代带来了较大的函数值下降, 我们就有理由相信这个方向上
m
k
(
p
)
m_k(p)
mk(p)是
f
f
f更大范围内的好近似, 从而扩大 (expand) 邻域; 反之, 则缩小 (shrink) 邻域. 所谓的信赖域子问题为
min
p
m
k
(
p
)
,
∥
p
∥
≤
Δ
k
,
\min_p m_k(p),\quad \Vert p\Vert\le\Delta_k,
pminmk(p),∥p∥≤Δk,这里
B
k
B_k
Bk是
∇
2
f
k
\nabla^2 f_k
∇2fk或者其近似. 求解信赖域子问题是信赖域方法的关键步骤.
信赖域与线搜索的区别还在于, 线搜索先求搜索方向
p
k
p_k
pk, 再求步长
α
k
\alpha_k
αk, 而信赖域则是二者同时进行.
m
k
m_k
mk选取的一个特例是, 令
B
k
=
0
B_k=0
Bk=0, 此时
m
k
m_k
mk为
f
f
f的一阶近似. 我们有基于此模型的Cauchy点. 这对探讨不同搜索方向的收敛性有着重要作用. 信赖域可以和线搜索结合起来使用. 例如,
B
k
B_k
Bk就取二阶信息
∇
2
f
k
\nabla^2 f_k
∇2fk, 这时为信赖域牛顿方法 (Trust-Region Newton Methods) ;
B
k
B_k
Bk以拟牛顿方式更新, 则有信赖域拟牛顿方法 (Trust-Region Quasi-Newton Methods).
3. 缩放
算法的表现极大取决于问题构成. 其中之一就是尺度. 在无约束优化中, 我们称一个问题尺度病态, 若对
x
x
x朝某一方向 (相对于其他方向而言) 微小移动会导致
f
f
f的极大变化. 例如
f
(
x
)
=
1
0
9
x
1
2
+
x
2
2
f(x)=10^9x_1^2+x_2^2
f(x)=109x12+x22, 对
x
1
x_1
x1的影响极其敏感, 而对
x
2
x_2
x2的影响可能无动于衷.
尺度病态的函数可能来自于物理或化学反应, 其中不同的过程以不同速率进行. 例如一个化学反应用有四个过程正在进行. 对应于每个过程都有一个反应速度. 比如说, 四个过程的反应速度大致具有如下量级:
x
1
≈
1
0
−
10
,
x
2
≈
x
3
≈
1
,
x
4
≈
1
0
5
.
x_1\approx10^{-10},\quad x_2\approx x_3\approx1,\quad x_4\approx10^5.
x1≈10−10,x2≈x3≈1,x4≈105.此时我们应当引入新变量
z
z
z:
(
x
1
x
2
x
3
x
4
)
=
(
1
0
−
10
0
0
0
0
1
0
0
0
0
1
0
0
0
0
1
0
5
)
(
z
1
z
2
z
3
z
4
)
,
\begin{pmatrix}x_1\\x_2\\x_3\\x_4\end{pmatrix}=\begin{pmatrix}10^{-10} & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\0 & 0 & 0 & 10^5\end{pmatrix}\begin{pmatrix}z_1\\z_2\\z_3\\z_4\end{pmatrix},
⎝⎜⎜⎛x1x2x3x4⎠⎟⎟⎞=⎝⎜⎜⎛10−1000001000010000105⎠⎟⎟⎞⎝⎜⎜⎛z1z2z3z4⎠⎟⎟⎞,再求解新变量下的问题. 这种放缩技巧称为对角线尺度变换.
有些算法对于尺度尤其敏感, 例如最速下降法; 而有些却少有影响, 如牛顿法. 我们可以在算法中加入尺度不变性, 例如迭代策略以及收敛判别. 一般说来, 线搜索比信赖域法更容易保持尺度不变.
事实上通过Shermann-Morrison公式我们可以分别得到DFP与BFGS的两种版本. ↩︎
共轭方向对于二次函数的优化具有重要作用, 而在局部 f f f是可以用二次模型近似的. ↩︎
1236
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
