最优化问题—非线性规划(二)
在之前的文章最优化问题—非线性规划(一)里面,我们主要关注了对于非线性规划的一般形式和最优条件,在最后我们介绍了关于凸规划的相关定义和其最优条件。下面,我们要介绍的是关于对于非线性规划问题的一般的求解思路。
1. 非线性规划的求解思路
1.1 迭代法以及其基本思想
对于非线性规划问题的求解,一般情况下,我们都是采用迭代法进行求解,这种方法的思路很简单,首先给定一个初始的点 x 0 ∈ R n x_0∈R^n x0∈Rn,按照某种迭代法则产生一个点 x k x_k xk,使得当 { x k } \{x_k\} {xk}是又穷点列是,其最后一个点是最优化模型问题的最优解:当 { x k } \{x_k\} {xk}是无穷点列的时候,它有一个极限点,其极限点是最优化模型的最优解。
举一个具体的例子来说,在最经典的BP神经网络中,我们随机初始化网络参数,正是通过这种迭代的思路在最小化误差的基础上来迭代生成最优的参数。
这里,我们在进一步来说明,在上述的思路中,我们提到了,迭代法需要根据某一种方法来产生一个点 x K x_K xK,那么我们应该如何选择这种迭代方法或者说如何评价这种方法的优劣呢?
迭代点列能够稳定的接近局部极小点 x ∗ x^* x∗的邻域,然后能够迅速收敛到 x ∗ x^* x∗,当给定的某种收敛规则满足之后,迭代终止。一般的,我们迭代产生的点的集合 { x k } \{x_k\} {xk}的极限为一个局部极小点,我们以此来作为一个迭代算法优劣的评价标准。
1.2 迭代法的核心步骤
在了解了基本思想之后,我们将目光聚焦在这种方法的核心步骤上。
首先,设
x
k
∈
R
n
x_k∈R^n
xk∈Rn是某种迭代算法的第k轮迭代点,
x
k
+
1
∈
R
n
x_{k+1}∈R^n
xk+1∈Rn是第k+1轮的迭代点,令
△
x
k
=
x
k
+
1
−
x
k
△x_k=x_{k+1}-x_k
△xk=xk+1−xk,也就是
x
k
+
1
=
x
k
+
△
x
k
x_{k+1}=x_k+△x_k
xk+1=xk+△xk。这里根据
x
k
和
x
k
+
1
x_k和x_{k+1}
xk和xk+1我们可以知道,
△
x
k
△x_k
△xk也是一个向量的形式,并且,这个向量表示以
x
k
x_k
xk为起点,
x
k
+
1
x_{k+1}
xk+1为重点的向量表示。
则进一步,我们去这个向量的方向,也就是令 p k = x k + 1 ∣ ∣ x k + 1 ∣ ∣ p_k=\frac{x_{k+1}}{||x_{k+1}||} pk=∣∣xk+1∣∣xk+1表示在 △ x k △x_k △xk方向上的单位向量。我们最后在取一个关于这个单位向量的长度系数 λ k λ_k λk,最后,对于 △ x k = λ k p k △x_k=λ_kp_k △xk=λkpk来进行最终的表示。
最后,我们表达式也变成了:
x
k
+
1
=
x
k
+
λ
k
p
k
x_{k+1}=x_k+λ_kp_k
xk+1=xk+λkpk
我们将 λ k λ_k λk称为步长,将 p k p_k pk称为第k轮迭代的搜索方向。由此可以知道,当我们选择不同的步长和不同的确定方向之后,这也就构成了不同的迭代方法。
小结一下,这里我们确定了如何进行更新的核心步骤,在这个步骤中,我们发现了一个新的问题,那就是如何选择步长 λ k λ_k λk和如何确定更新的方向 p k p_k pk。关于这两个步骤的求解,我们将在下面以及之后的文章中进行具体的介绍。
1.3 迭代法的基本流程
在了解了迭代法的核心步骤以及基本思路之后,我们对迭代法的流程进行描述:
- 确定搜索方向 p k p_k pk:即依照一定的准则,构造f在 x k x_k xk点处的下降方向作为搜索方向。
- 确定步长λ,是目标函数值有某种意义的下降。
- 对当前搜索的点进行更新,也及时 x k + 1 = x k + λ k p k x_{k+1}=x_{k}+λ_kp_k xk+1=xk+λkpk,如果 x k + 1 x_{k+1} xk+1满足某种终止条件,则停止迭代,得到近似最优解 x k + 1 x_{k+1} xk+1,否则,重复上面的步骤。
算法框架如下图所示:
这这一段的描述中,我们除了对于步长和方向的求解方法不确定之外,还引入了另外的一个问题,那就是如何控制迭代的停止?也就是应该如何定义停止条件?
小结一下,通过该段的介绍,能够知道迭代法的基本流程,根据上一段和本段的,描述,我们也知道的了下面需要解的问题包括步长,方向和迭代停止条件。
1.4 方向问题
在这一节中,我们首先要解决的是方向问题,这里我们需要先了解到下降方向 和 可行下降方向。
下降方向: 是指对于目标函数 f : R n − > R 1 , x − ∈ R n f:R^n->R^1,x^-∈R^n f:Rn−>R1,x−∈Rn,向量p∈ R n R^n Rn(p≠0),如果存在δ>0,任意的λ∈(0,δ),都有 f ( x − + λ p ) < f ( x − ) f(x^-+λp)<f(x^-) f(x−+λp)<f(x−),则称向量p为函数在 x − x^- x−的下降方向。同时,凡是满足这种迭代性质的最优化方法都可以成为下降方法。
通俗一点的来说,我们所获得的方向p,要能够保证x的邻域内的其他的点的函数值是要小于x的函数值,这种我们能够保证在不断的对点x进行迭代更新的时候,函数值是不断下降的。
可行下降方向: 在有约束的规划问题中,我们首先要保证的是x在其可行域内进行迭代,也就是说,x在通过步长和方向进行转移之后新生成的点也必须保证在可行域内。这种下降的方向称为可行下降方向。
在上篇文章中,我们着重回顾了梯度这个概念,梯度表示了在当前点x的最快的是的函数值增大的方向,同理,当我们把梯度值取到相反的方向的时候,我们能够获得也是最快下降的一个方向。所以这里我们就给出了使用负梯度作为下降方向的方法。当然,这是一种最为常见的方法,有兴趣的读者可以可以查阅求他的关于下降反向如何确定的方法。
1.5 步长因子
在解决了关于下降方向如何确定之后,我们来关注上面的第二个问题,如何确定步长因子?这里,我先回顾目标函数的形式为
f
(
x
)
f(x)
f(x),在进行迭代之后,函数值变成了
f
(
x
+
λ
p
)
f(x+λp)
f(x+λp),步长因子是一个标量,p代表一个反向,我们利用上一段的方法以及可以进行确定了。换一个角度来讲,我们之前是将输入的点x作为f的变量,同理对于位置量λ而言,其也可以是f的一个变量。我们用公式来表示就是:
f
(
x
k
+
λ
k
p
k
)
=
−
m
i
n
λ
f
(
x
k
+
λ
p
k
)
f(x_k+λ_kp_k)=-min_{λ}f(x_k+λp_k)
f(xk+λkpk)=−minλf(xk+λpk)
对于第k+1轮迭代的函数值,显然我们已经能够知道当前的点
x
k
x_k
xk,以及在
x
k
x_k
xk处的下降方向
p
k
p_k
pk,λ就顺理成章的称为了我们在第k+1轮要确定的参数。同时,我们希望尽量的减少迭代的次数,这样可以降低算法的时间复杂度。那么我们就希望每一次能够尽量的去降低目标函数值。这也就是上面公式成立的意义。
在确定了最优步长和下降方向之后,就可以给出一个使用最优步长和下降方向进行更新的定理:
定理:设目标函数f(x)具有连续的一阶偏导数, x k + 1 x_{k+1} xk+1则可以有下面的规则进行产生:
f
(
x
k
+
λ
k
p
k
)
=
m
i
n
λ
f
(
x
k
+
λ
p
k
)
f(x_k+λ_kp_k)=min_{λ}f(x_k+λp_k)
f(xk+λkpk)=minλf(xk+λpk)
x
k
+
1
=
x
k
+
λ
k
p
k
x_{k+1}=x_k+λ_kp_k
xk+1=xk+λkpk
则有,
▽
f
(
x
k
+
1
)
T
p
k
=
0
▽f(x_{k+1})^Tp_k=0
▽f(xk+1)Tpk=0
前两个公式想必大家能够明白其中的含义了,我们来看最优的这个公式,在这个公式中,其说明了一维的最优搜索步长
λ
k
λ_k
λk,所对应的点
x
k
+
1
x_{k+1}
xk+1的目标函数的梯度
▽
f
(
x
k
+
1
)
T
▽f(x_{k+1})^T
▽f(xk+1)T与下降方向
p
p
p正交。
这个部分很好说明,我们以步长λ作为变量,f作为目标函数,其最优的步长
λ
k
λ_k
λk肯定对应着函数f关于λ的导数为0(极值点)。那么则有
∂
f
∂
λ
k
=
0
\frac{∂f}{∂λ_k}=0
∂λk∂f=0,则将f的公式展开就有
▽
f
(
x
k
+
1
)
T
p
k
=
0
▽f(x_{k+1})^Tp_k=0
▽f(xk+1)Tpk=0
对于最优步长的搜索,一般是采用线性搜索的方式来确定最优的λ,对于线性搜索的相关思想和方法,我们将在下一篇文章中进行具体的介绍。
1.6 小结
通过上几个小结的描述,我们知道了如何确定方向,如何确定步长(具体算法在下一篇文章中说明)。显然这种利用梯度作为方向和线性搜索最优步长的方法只能算作是一类求解非线性规划问题的方法。其不能适应与所有的非线性最优化问题,最简单的原因在于我们无法保证目标函数是有梯度值的。
对于这种方法而言,其迭代映射的方法我们假设为
T
:
R
n
−
>
R
n
T:R^n->R^n
T:Rn−>Rn,则下一个迭代的点为:
x
k
+
1
=
T
(
x
k
)
x_{k+1}=T(x_k)
xk+1=T(xk)
而我们最终迭代的目的是希望能够找到一个使得目标函数最小的极值点
x
∗
x^*
x∗,并且迭代最终能够停在极值点上面,也就是说
x
∗
=
T
(
x
∗
)
x^*=T(x^*)
x∗=T(x∗)。显然,这个最优的极值点就是一个不动点。对于这种迭代查找的方式,我们只需要这种迭代映射是关于某一种范数是收缩映射。也就是||T(x)-T(y)||≤β||x-y||,任意的x,y∈
R
n
R^n
Rn,β<1。这意味这这种映射T有唯一的不动点,从任意的初始点x开始,迭代之后一定能够收敛域
x
∗
x^*
x∗点。
上面这段话主要是从数学的角度出发来进行阐述的,从理解的层面,我们想要找到一个极值点 x ∗ x^* x∗,我们需要的是这个极值点是能够收敛到的,如果我们当前的迭代算法无法收敛,那么也就找不到最优解。而通过某一种收缩的迭代映射T,我们能够保证这种迭代是能够收敛到极值点的。
举一个例子来说:
2 终止条件与收敛速度
在整个上一节的介绍中,我们给出了非线性规划的问题使用迭代法的求解思路,包括如何确定步长,如何确定下降方向等等。在本节中,我们重点讨论关于如何停止迭代的方法和算法收敛的速度。
2.1 终止条件
- 最自然的停止准则: ∣ ∣ ▽ f ( x k ) ∣ ∣ < ε ||▽f(x_k)||<ε ∣∣▽f(xk)∣∣<ε,其中 ε ε ε为实现定义的阈值,这表明当梯度趋于0,并且迭代序列 { x k ∥ \{x_k\| {xk∥收敛到某一个点。
- 相邻两次迭代点的关系: ∣ ∣ x k + 1 − x k ∣ ∣ < ε ||x_{k+1}-x_k||<ε ∣∣xk+1−xk∣∣<ε,或者 ∣ ∣ x k + 1 − x k ∣ ∣ ∣ ∣ x k ∣ ∣ < ε \frac{||x_{k+1}-x_k||}{||x_k||}<ε ∣∣xk∣∣∣∣xk+1−xk∣∣<ε。这个停止条件将注意力关注再来“收敛”上面,当不存在梯度的时候,可以采用这种方法。
- 相邻两次迭代函数值的关系: f ( x k ) − f ( x k + 1 ) < ε f(x_k)-f(x_{k+1})<ε f(xk)−f(xk+1)<ε,或者 f ( x k ) − f ( x k + 1 ) f ( x k ) < ε \frac{f(x_k)-f(x_{k+1})}{f(x_k)}<ε f(xk)f(xk)−f(xk+1)<ε。
- 同时使用迭代值和迭代函数作为收敛的条件,当 ∣ ∣ x k ∣ ∣ > ε ||x_k||>ε ∣∣xk∣∣>ε和| ∣ f ( x k ) ∣ > ε |f(x_k)|>ε ∣f(xk)∣>ε时,采用 ∣ ∣ x k + 1 − x k ∣ ∣ ∣ ∣ x k ∣ ∣ ≤ ε ′ , ∣ f ( x k + 1 ) − f ( x k ) ∣ ∣ f ( x k ) ∣ ≤ ε ′ \frac{||x_{k+1}-x_k||}{||x_k||}≤ε',\frac{|f(x_{k+1})-f(x_k)|}{|f(x_k)|}≤ε' ∣∣xk∣∣∣∣xk+1−xk∣∣≤ε′,∣f(xk)∣∣f(xk+1)−f(xk)∣≤ε′,否则采用 ∣ ∣ x k + 1 − x k ∣ ∣ < ε ′ , ∣ f ( x k ) − f ( x k + 1 ) ∣ ≤ ε ′ ||x_{k+1}-x_{k}||<ε',|f(x_k)-f(x_{k+1})|≤ε' ∣∣xk+1−xk∣∣<ε′,∣f(xk)−f(xk+1)∣≤ε′。
2.2 收敛速度
收敛速度是算法的局部刻画,其主要度量了优化算法的有效性,使得迭代序列 { x k } \{x_k\} {xk}依照某一个范数收敛到 x ∗ x^* x∗。
Q-α收敛速度:
l
i
m
k
−
>
∞
∣
∣
x
k
−
x
∗
∣
∣
=
0
lim_{k->∞}||x_k-x^*||=0
limk−>∞∣∣xk−x∗∣∣=0,如果存在实数α≥1,和迭代次数K无关的正数β>0,使得:
l
i
m
k
−
>
∞
∣
∣
x
k
+
1
−
x
∗
∣
∣
∣
∣
x
k
−
x
∗
∣
∣
α
lim_{k->∞}\frac{||x_{k+1}-x^*||}{||x_k-x^*||^α}
limk−>∞∣∣xk−x∗∣∣α∣∣xk+1−x∗∣∣
这称为算法产生的点列具有
Q
−
α
Q-α
Q−α阶收敛速度。
- 当α=1,β∈(0,1),迭代点列 { x k } \{x_k\} {xk}叫做具有Q-线性收敛。
- 当α∈(1,2),β>0或者α=1,β=0的时候, { x k } \{x_k\} {xk}叫做具有Q-超线性收敛。
- 当α=2, { x k } \{x_k\} {xk}叫做具有Q-二阶收敛。
相关定理:如果迭代点列 { x k } \{x_k\} {xk}具有Q-超线性收敛到 x ∗ x^* x∗,则 l i m k − > ∞ ∣ ∣ x k + 1 − x ∗ ∣ ∣ ∣ ∣ x k − x ∗ ∣ ∣ lim_{k->∞}\frac{||x_{k+1}-x^*||}{||x_k-x^*||} limk−>∞∣∣xk−x∗∣∣∣∣xk+1−x∗∣∣,反之一般不成立。
2.3 如何选择一个好的停止条件?
对于具有一阶的导数信息,其收敛速度不太快的算法,例如共轭梯度方法,可以采用第一种的终止条件: ∣ ∣ ▽ f ( x k ) ∣ ∣ < ε ||▽f(x_k)||<ε ∣∣▽f(xk)∣∣<ε。但是由于临界点点可能是鞍点,因此建议可以结合上述的第一种和第四种方法综合来进行使用,综合使用的时候,一般第一种方法中的 ε = 1 0 − 4 ε=10^{-4} ε=10−4,第二种方法中的 ε = ε ′ = 1 0 − 5 ε=ε'=10^{-5} ε=ε′=10−5。
最后,我们对终止条件简单的扩展的一下,这里我们先介绍二次终止性的概念:
一个算法用于其求解正定二次函数的无约束极小值的时候,有限的迭代步骤可以到达最优解,则该算法具有二次终止性。
大多数的优化算法本身是基于二次函数的模型导出,因此在极小值点的附件可以采用二次函数进行很好的近似,如果算法具有二次终止性,则会很快的进行收敛。
2.4 共轭方向
在上一个小节里面,我们主要给出的是一般的算法的终止条件,并在最后强调的二次函数的优势和二次终止性的概念。
进一步,我们知道了二次终止性的概念,对于二次终止性而言,二次终止性=共轭方向+精确的一维搜索。一维搜索的过程,我们在下一篇文章中具体介绍。这里我们首先来介绍共轭方向的问题。
共轭方向: 设 A n ∗ n A_{n*n} An∗n为对称的正定矩阵,若 p 1 , p 2 ∈ R n , p 1 ≠ 0 , p 2 ≠ 0 p_1,p_2∈R^n,p_1≠0,p_2≠0 p1,p2∈Rn,p1=0,p2=0,满足 p 1 T A p 2 = 0 p_1^TAp_2=0 p1TAp2=0,则称 p 1 , p 2 p_1,p_2 p1,p2关于矩阵A共轭,如果m个不同的非零向量 p 1 , p 2 , . . . , p m ∈ R n , m < n p_1,p_2,...,p_m∈R^n,m<n p1,p2,...,pm∈Rn,m<n,满足 p i T A p j = 0 , i f i ≠ j p_i^TAp_j=0,if i≠j piTApj=0,ifi=j,则称这m个向量为共轭向量组,或者A共轭,如果A=I(单位矩阵),则共轭等价于正交。
换句话来讲,当两个向量 p i , p j p_i,p_j pi,pj关于对称正定矩阵A的乘积结果 p 1 T A p j = 0 p_1^TAp_j=0 p1TApj=0,则说明 p i 和 p j p_i和p_j pi和pjl两个向量代表的方向是共轭的方向。
2.5 一般的共轭方向法(GCD)
在介绍完共轭方向的概念之后,我们下面给出的是一般的使用共轭方向的最优化的方法。
- 给定初始点 x 0 , ε > 0 , k = 0 x_0,ε>0,k=0 x0,ε>0,k=0,计算 ▽ f ( x ) ∣ x 0 = g ( x 0 ) ▽f(x)|x_0=g(x_0) ▽f(x)∣x0=g(x0),再计算方向 p 0 p_0 p0,使得 p 0 T g 0 < 0 p_0^Tg_0<0 p0Tg0<0;
- 如果 ∣ ∣ g k ∣ ∣ ≤ ε ||g_k||≤ε ∣∣gk∣∣≤ε,则迭代终止。
- 计算步长因子 λ k λ_k λk,令 x k + 1 = x k + λ k p k x_{k+1}=x_k+λ_kp_k xk+1=xk+λkpk。
- 通过某种共轭方向算法计算 p k + 1 p_{k+1} pk+1,使得 p k + 1 T G p j = 0 , j = 0 , 1 , , , , , k p_{k+1}^TGp_{j}=0,j=0,1,,,,,k pk+1TGpj=0,j=0,1,,,,,k。
- 令k=k+1,转移到第二步。
总结一下,共轭方向法实际上针对的是梯度下降的方向来进行改变的。同时,在这种方式中,如果函数正定的,则GCD方法最多经过n步精确线性搜索终止。同时还可以证明的是,如果多个向量关于正定矩阵A两两是共轭的,则这些向量线性无关。
3 总结
在整个的这一个部分中,我们重点介绍了对于非线性规划算法的一般的求解思路和求解框架,并介绍了下降方向步长,终止条件,收敛速度等概念,其中也给出量常见的梯度迭代的算法和共轭方向的迭代算法。
3.1 回顾步长问题
在上文的描述中,我们有一个问题的始终没有得到具体的解决,那就是步长因子的确定。对于最优步长而言,其求解的算法有很多,但是其一般不具有解析解。此时,我们选择的是采用近似求解的方式。其求解过程如下:
- 首先,确定包含最小值的搜索区间。
- 其次,利用某些区间的分割技术或者插值技术来迭代减少区间,一直到最后的区间宽度小于预先设定的精度,就可以获得近似解。
具体而言,这种搜索近似解的方法可以分成试探法和插值法。这些方法我们将在下一篇文章中进行介绍。
4 参考
- 哈工大—组合优化与凸优化