牛顿法
梯度法仅仅依赖函数值和梯度的信息,如果函数 f ( x ) f(x) f(x)充分光滑,则可以利用二阶导数信息构造下降方向 d k d^k dk。牛顿类算法就是利用二阶导数信息来构造迭代格式的算法。由于利用信息变多,牛顿法的实际表现要好于梯度法,但是它对函数 f ( x ) f(x) f(x)的要求也相应变高。
一、经典牛顿法
对于二次函数
f
(
x
)
f(x)
f(x),考虑
f
(
x
)
f(x)
f(x)在点
x
k
x^k
xk处的二阶泰勒展开
f
(
x
)
=
f
(
x
k
)
+
∇
f
(
x
k
)
T
(
x
−
x
k
)
+
(
x
−
x
k
)
T
1
2
∇
2
f
(
x
k
)
(
x
−
x
k
)
+
o
(
∣
∣
x
−
x
k
∣
∣
2
)
(1)
f(x)=f(x^k)+\nabla f(x^k)^T(x-x^k)+(x-x^k)^T\frac{1}{2}\nabla^2f(x^k)(x-x^k) + o(||x-x^k||^2)\tag{1}
f(x)=f(xk)+∇f(xk)T(x−xk)+(x−xk)T21∇2f(xk)(x−xk)+o(∣∣x−xk∣∣2)(1)
对于(1)式忽略高阶项,极小化,根据一阶极小化的必要性条件,可以得到
∇
f
(
x
k
)
+
∇
2
f
(
x
k
)
(
x
−
x
k
)
=
0
\nabla f(x^k)+\nabla^2f(x^k)(x-x^k)=0
∇f(xk)+∇2f(xk)(x−xk)=0
若
∇
2
f
(
x
k
)
\nabla^2f(x^k)
∇2f(xk)可逆,于是有
x
=
x
k
−
(
∇
2
f
(
x
k
)
)
−
1
∇
f
(
x
k
)
x=x^k-(\nabla^2f(x^k))^{-1}\nabla f(x^k)
x=xk−(∇2f(xk))−1∇f(xk)
从而可得经典牛顿法的迭代格式
x
k
+
1
=
x
k
−
(
∇
2
f
(
x
k
)
)
−
1
∇
f
(
x
k
)
(2)
x^{k+1}=x^k-(\nabla^2f(x^k))^{-1}\nabla f(x^k) \tag{2}
xk+1=xk−(∇2f(xk))−1∇f(xk)(2)
注意在经典迭代格式中,步长
α
k
\alpha_k
αk恒为1,即可以不额外考虑步长的选取,我们也称步长为1的牛顿法为经典牛顿法。
我们这里不加证明给出经典牛顿法的局部收敛性
(1)如果初始点离
x
∗
x^*
x∗足够近,那么牛顿法产生的点列
{
x
k
}
\{x^k\}
{xk}收敛于
x
∗
x^*
x∗;
(2)
{
x
k
}
\{x^k\}
{xk}收敛到点
x
∗
x^*
x∗的速度是二次的
(3)
{
∣
∣
∇
f
(
x
k
)
∣
∣
}
\{||\nabla f(x^k)||\}
{∣∣∇f(xk)∣∣}二次收敛到0
上述性质表明经典牛顿法是收敛速度很快的算法,但它的收敛是有条件的:第一初始点
x
0
x^0
x0必须距离问题的解充分近,即牛顿法只有局部收敛性,当
x
0
x^0
x0距离问题的解较远时,牛顿算法在多数情况下会失效;第二,海塞矩阵
∇
2
f
(
x
∗
)
\nabla^2f(x^*)
∇2f(x∗)需要为正定矩阵;牛顿法适用于优化问题的高精度求解,但是它没有全局收敛性,因此在实际应用中,人们常会使用梯度类算法先求得较低精度的解,然后调用牛顿法来获得高精度的解。
二、修正牛顿法
尽管上节给出了经典牛顿法的迭代格式和收敛性质,在实际应用中这个迭代格式几乎时不能使用的,经典牛顿法有如下缺陷:
(1)每一步迭代需要求解一个
n
n
n维线性方程组,这导致在高维问题中计算量较大,海塞矩阵
∇
2
f
(
x
k
)
\nabla^2f(x^k)
∇2f(xk)既不容易计算又不容易存储;
(2)当
∇
2
f
(
x
k
)
\nabla^2f(x^k)
∇2f(xk)不正定时,牛顿方程给出的解的性质通常比较差
(3)当迭代点距离最优值较远时,直接选取步长
a
l
p
h
a
k
=
1
alpha_k=1
alphak=1会使得迭代机器不稳定,在有些情况下迭代点列会发散。
为了克服这些缺点,我们必须对经典牛顿法做出某些修正或者变形,使其成为真正可以使用的算法。这里介绍带线搜索的修正牛顿法,其基本思想是对牛顿方程中的海塞矩阵
∇
2
f
(
x
k
)
\nabla^2f(x^k)
∇2f(xk)进行修正,使其变成正定矩阵;同时引入线搜索以改善算法的稳定性;其算法迭代步骤如下
带线搜索的牛正牛顿法:
(1)给定初始点
x
0
x^0
x0;
(2)for k=0,1,2,…do
(3) 确定矩阵
E
k
E^k
Ek使得矩阵
B
k
=
∇
2
f
(
x
k
)
+
E
k
B^k=\nabla^2f(x^k)+E^k
Bk=∇2f(xk)+Ek正定且条件数小
(4) 求解修正的牛顿方程
B
k
(
d
k
)
=
−
∇
f
(
x
k
)
B^k(d^k)=-\nabla f(x^k)
Bk(dk)=−∇f(xk),得到下降方向
d
k
d^k
dk
(5) 使用任意一种线搜索准则确定步长
α
k
\alpha_k
αk
(6) 更行
x
k
+
1
=
x
k
+
α
k
d
k
x^{k+1}=x^k+\alpha_kd^k
xk+1=xk+αkdk
(7)endfor
一个直接的取法是取 E k = τ k I E^k=\tau_kI Ek=τkI,即取单位矩阵的常数倍,根据矩阵理论可知,当 τ k \tau_k τk充分大时,总可以保证 B k B^k Bk是正定矩阵,然而 τ k \tau_k τk不宜取得过大,这是因为 τ k \tau_k τk趋于无穷时, d k d^k dk的方向会接近负梯度方向。
由于经典牛顿法和修正牛顿法均需要对海塞矩阵求逆,因此牛顿法在工程优化设计中很少直接被利用,常常采用的是通过迭代的方式逐步计算处接近牛顿法下降方向的算法,这种方法称为拟牛顿法。由于工程应用很少,本节就不给出编程实现了。