按秩1法(详见博文《最优化方法Python计算:秩1拟牛顿法》)计算的修正矩阵
Q
k
+
1
=
Q
k
+
E
k
\boldsymbol{Q}_{k+1}=\boldsymbol{Q}_k+\boldsymbol{E}_k
Qk+1=Qk+Ek无法保证其正定性。这时,
d
k
+
1
=
−
Q
k
+
1
g
k
+
1
\boldsymbol{d}_{k+1}=-\boldsymbol{Q}_{k+1}\boldsymbol{g}_{k+1}
dk+1=−Qk+1gk+1可能不是
f
(
x
)
f(\boldsymbol{x})
f(x)在
x
k
+
1
\boldsymbol{x}_{k+1}
xk+1处的下降方向,致使算法失败。要摆脱秩1法的这一窘境,需另辟蹊径。很自然的想法是“秩2”修正:设
Q
k
\boldsymbol{Q}_k
Qk对称正定,令修正矩阵
E
k
=
Δ
x
k
Δ
x
k
⊤
Δ
x
k
⊤
Δ
g
k
−
Q
k
Δ
g
k
Δ
g
k
⊤
Q
k
Δ
g
k
⊤
Q
k
Δ
g
k
\boldsymbol{E}_k=\frac{\Delta\boldsymbol{x}_k\Delta\boldsymbol{x}_k^\top}{\Delta\boldsymbol{x}_k^\top\Delta\boldsymbol{g}_k}-\frac{\boldsymbol{Q}_k\Delta\boldsymbol{g}_k\Delta\boldsymbol{g}_k^\top\boldsymbol{Q}_k}{\Delta\boldsymbol{g}_k^\top\boldsymbol{Q}_k\Delta\boldsymbol{g}_k}
Ek=Δxk⊤ΔgkΔxkΔxk⊤−Δgk⊤QkΔgkQkΔgkΔgk⊤Qk
于是得到
Q
k
+
1
\boldsymbol{Q}_{k+1}
Qk+1的秩2修正公式
Q
k
+
1
=
Q
k
+
E
k
=
Q
k
+
Δ
x
k
Δ
x
k
⊤
Δ
x
k
⊤
Δ
g
k
−
Q
k
Δ
g
k
Δ
g
k
⊤
Q
k
Δ
g
k
⊤
Q
k
Δ
g
k
.
(
1
)
\boldsymbol{Q}_{k+1}=\boldsymbol{Q}_k+\boldsymbol{E}_k=\boldsymbol{Q}_k+\frac{\Delta\boldsymbol{x}_k\Delta\boldsymbol{x}_k^\top}{\Delta\boldsymbol{x}_k^\top\Delta\boldsymbol{g}_k}-\frac{\boldsymbol{Q}_k\Delta\boldsymbol{g}_k\Delta\boldsymbol{g}_k^\top\boldsymbol{Q}_k}{\Delta\boldsymbol{g}_k^\top\boldsymbol{Q}_k\Delta\boldsymbol{g}_k}.\quad\quad(1)
Qk+1=Qk+Ek=Qk+Δxk⊤ΔgkΔxkΔxk⊤−Δgk⊤QkΔgkQkΔgkΔgk⊤Qk.(1)
利用式(1)作为正定阵
Q
k
+
1
\boldsymbol{Q}_{k+1}
Qk+1修正公式的拟牛顿法是由Broyden,Fletcher,Goldfarb和Shanno在20世纪70年代各自独立提出来的,故常称为BFGS算法。可以证明:
定理1 (1)设
Q
k
\boldsymbol{Q}_k
Qk对称正定,由式(1)确定的
Q
k
+
1
\boldsymbol{Q}_{k+1}
Qk+1正定,当且仅当
Δ
x
k
⊤
Δ
g
k
>
0
\Delta\boldsymbol{x}_k^\top\Delta\boldsymbol{g}_k>0
Δxk⊤Δgk>0。
(2)设目标函数
f
(
x
)
f(\boldsymbol{x})
f(x),
x
∈
R
n
\boldsymbol{x}\in\text{R}^n
x∈Rn一阶连续可微,且有极小值点
x
0
\boldsymbol{x}_0
x0。则BFGS算法每次迭代均有
Δ
x
k
⊤
Δ
g
k
>
0
\Delta\boldsymbol{x}_k^\top\Delta\boldsymbol{g}_k>0
Δxk⊤Δgk>0。
BFGS算法是一个改进了的拟牛顿算法,读者可作为练习用Python实现BFGS算法。Python的scipy.optimize为用户提供了BFGS方法,只需要在调用minimize时将’BFGS’传递给method参数即可用BFGS方法计算目标函数的最优解。
例1 用scipy.optimize提供的BFGS方法计算Rosenbrock函数的最优解,给定初始点
x
1
=
(
100
100
)
\boldsymbol{x}_1=\begin{pmatrix}100\\100\end{pmatrix}
x1=(100100)。
解:下列代码完成本例计算。
import numpy as np #导入numpy
from scipy.optimize import rosen, minimize #导入rosen, minimize
x=np.array([100,100]) #设置初始点
res=minimize(rosen,x,method='BFGS') #计算最优解
print(res)
运行程序,输出
fun: 1.8831204186846363e-11
hess_inv: array([[0.49113161, 0.98272927],
[0.98272927, 1.97132641]])
jac: array([ 2.16488885e-06, -9.42470479e-07])
message: 'Optimization terminated successfully.'
nfev: 1488
nit: 385
njev: 496
status: 0
success: True
x: array([0.99999566, 0.99999131])
这意味着BFGS方法从初始点
x
1
=
(
100
100
)
\boldsymbol{x}_1=\begin{pmatrix}100\\100\end{pmatrix}
x1=(100100)起,迭代385次算得Rosenbrock函数的最优解
x
0
\boldsymbol{x}_0
x0的近似值为
(
0.99999566
0.99999131
)
\begin{pmatrix}0.99999566\\0.99999131\end{pmatrix}
(0.999995660.99999131)。虽然运行效率未必优于秩1算法(见博文《最优化方法Python计算:秩1拟牛顿法》中例),但根据定理1,算法运行的可靠性得到了保证。
写博不易,敬请支持:
如果阅读本文于您有所获,敬请点赞、评论、收藏,谢谢大家的支持!