统计学习方法 牛顿法和拟牛顿法
学习李航的《统计学习方法》时,关于牛顿法和拟牛顿法的笔记。
牛顿法(Newton method)和拟牛顿法(quasi-Newton method)时求解无约束优化问题的常用方法,有收敛速度快的优点。牛顿法时迭代算法,每一步需要求解目标函数的 Hession 矩阵的逆矩阵,计算较为复杂;拟牛顿法通过正定矩阵近似 Hession 矩阵或其逆矩阵,简化了计算过程。
牛顿法
牛顿法的推导:考虑无约束最优化问题:
min
x
∈
R
n
f
(
x
)
\min\limits_{x\in \R^n} f(x)
x∈Rnminf(x)
设
x
∗
x^\ast
x∗ 是目标函数的极小点。
假设
f
(
x
)
f(x)
f(x) 具有二姐连续偏导数,若第
k
k
k 次迭代的值为
x
(
k
)
x^{(k)}
x(k) ,则对
f
(
x
)
f(x)
f(x) 在
x
(
k
)
x^{(k)}
x(k) 附近二阶泰勒展开:
f
(
x
)
=
f
(
x
(
k
)
)
+
g
k
T
(
x
−
x
(
k
)
)
+
1
2
(
x
−
x
(
k
)
)
H
k
(
x
−
x
(
k
)
)
f(x)=f(x^{(k)})+g_k^\mathrm{T}(x-x^{(k)})+\frac{1}{2}(x-x^{(k)})H_k(x-x^{(k)})
f(x)=f(x(k))+gkT(x−x(k))+21(x−x(k))Hk(x−x(k))
其中
g
k
g_k
gk 为
f
(
x
)
f(x)
f(x) 在
x
(
k
)
x^{(k)}
x(k) 处的梯度向量,
H
k
H_k
Hk 为
f
(
x
)
f(x)
f(x) 的 Hession 矩阵在
x
(
k
)
x^{(k)}
x(k) 处的值:
g
(
x
)
=
∇
f
(
x
)
=
[
∂
f
∂
x
i
]
n
×
1
,
H
(
x
)
=
[
∂
2
f
∂
x
i
∂
x
j
]
n
×
n
g(x)=\nabla f(x)=\left[\frac{\partial f}{\partial x_i}\right]_{n\times 1},\quad H(x)=\left[ \frac{\partial^2 f}{\partial x_i\partial x_j} \right]_{n\times n}
g(x)=∇f(x)=[∂xi∂f]n×1,H(x)=[∂xi∂xj∂2f]n×n
函数
f
(
x
)
f(x)
f(x) 有极值的必要条件是在极值点处一阶导数(即梯度向量)为
0
0
0 。特别地,当
H
(
x
)
H(x)
H(x) 在极值点处是正定矩阵时,函数
f
(
x
)
f(x)
f(x) 的极值为极小值。
牛顿法利用极小点的必要条件,每次从上一次迭代得到的极小点
x
(
k
)
x^{(k)}
x(k) 开始,得到当前目标函数的极小点
x
(
k
+
1
)
x^{(k+1)}
x(k+1) ,即:
∇
f
(
x
(
k
+
1
)
)
=
g
k
+
H
k
(
x
(
k
+
1
)
−
x
(
k
)
)
=
0
\nabla f(x^{(k+1)})=g_k+H_k(x^{(k+1)}-x^{(k)})=0
∇f(x(k+1))=gk+Hk(x(k+1)−x(k))=0
解得:
x
(
k
+
1
)
=
x
(
k
)
−
H
k
−
1
g
k
x^{(k+1)}=x^{(k)}-H_k^{-1}g_k
x(k+1)=x(k)−Hk−1gk
假设
p
k
p_k
pk 为一
n
×
1
n\times 1
n×1 的列向量,且满足:
H
k
p
k
=
−
g
k
H_kp_k=-g_k
Hkpk=−gk
则迭代公式可写为:
x
(
k
+
1
)
=
x
(
k
)
+
p
k
x^{(k+1)}=x^{(k)}+p_k
x(k+1)=x(k)+pk
算法:牛顿法
- 输入:目标函数 f ( x ) f(x) f(x) ,梯度 g ( x ) = ∇ f ( x ) g(x)=\nabla f(x) g(x)=∇f(x) ,Hession 矩阵 H ( x ) H(x) H(x) ,精度 ε \varepsilon ε ;
- 输出: f ( x ) f(x) f(x) 的极小值点 x ∗ x^\ast x∗ ;
- 取初始点 x ( 0 ) x^{(0)} x(0) ,置 k = 0 k=0 k=0 ;
- 计算 g k = g ( x ( k ) ) g_k=g(x^{(k)}) gk=g(x(k)) ;
- 若 ∥ g k ∥ < ε \|g_k\|\lt \varepsilon ∥gk∥<ε ,则停止计算,返回近似解 x ∗ = x ( k ) x^\ast=x^{(k)} x∗=x(k) ;
- 计算 H k = H ( x ( k ) ) H_k=H(x^{(k)}) Hk=H(x(k)) ,并求 p k p_k pk :
H k p k = − g k H_kp_k=-g_k Hkpk=−gk
- 更新 x ( k + 1 ) = x ( k ) + p k x^{(k+1)}=x^{(k)}+p_k x(k+1)=x(k)+pk , k = k + 1 k=k+1 k=k+1 ,跳转至 2;
拟牛顿法
牛顿法中需要计算 H k − 1 H_k^{-1} Hk−1 ,比较耗时,我们考虑使用一个与 Hession 矩阵具有类似性质的 n n n 阶矩阵 G k = G ( x ( k ) ) G_k=G(x^{(k)}) Gk=G(x(k)) 来近似代替 H k − 1 H_k^{-1} Hk−1 。
条件一:首先,对
x
(
k
)
x^{(k)}
x(k) 处的二阶泰勒展开进行求导,取
x
=
(
k
+
1
)
x=^{(k+1)}
x=(k+1) ,得:
g
k
+
1
−
g
k
=
H
k
(
x
(
k
+
1
)
−
x
(
k
)
)
g_{k+1}-g_k=H_k(x^{(k+1)}-x^{(k)})
gk+1−gk=Hk(x(k+1)−x(k))
记
y
k
=
g
k
+
1
−
g
k
y_k=g_{k+1}-g_k
yk=gk+1−gk ,
δ
k
=
x
(
k
+
1
)
−
x
(
k
)
\delta_{k}=x^{(k+1)}-x^{(k)}
δk=x(k+1)−x(k) ,则满足:
y
k
=
H
k
δ
k
y_k=H_k\delta_k
yk=Hkδk
或:
H
k
−
1
y
k
=
δ
k
H_k^{-1}y_k=\delta_k
Hk−1yk=δk
该条件称为拟牛顿条件。
条件二:如果
H
k
H_k
Hk 是正定的(
H
k
−
1
H_k^{-1}
Hk−1 也是正定的),那么可以保证牛顿法搜索方向
p
k
p_k
pk 是下降方向,因为:
x
(
k
+
1
)
=
x
(
k
)
+
λ
p
k
=
x
(
k
)
−
λ
H
k
−
1
g
k
x^{(k+1)}=x^{(k)}+\lambda p_k=x^{(k)}-\lambda H_k^{-1}g_k
x(k+1)=x(k)+λpk=x(k)−λHk−1gk
带入前面的泰勒展开式,为:
f
(
x
(
k
+
1
)
)
=
f
(
x
k
)
−
λ
g
k
T
H
k
−
1
g
k
+
1
2
λ
2
g
k
T
H
k
−
T
H
k
H
k
−
1
g
k
≈
f
(
x
k
)
−
λ
g
k
T
H
k
−
1
g
k
f(x^{(k+1)})=f(x^{k})-\lambda g_k^{\mathrm{T}}H_k^{-1}g_k+\frac{1}{2}\lambda^2g_k^\mathrm{T}H_k^{-\mathrm{T}}H_kH_{k}^{-1}g_k\approx f(x^{k})-\lambda g_k^{\mathrm{T}}H_k^{-1}g_k
f(x(k+1))=f(xk)−λgkTHk−1gk+21λ2gkTHk−THkHk−1gk≈f(xk)−λgkTHk−1gk
当
λ
\lambda
λ 为一个足够小的正数时,由于
H
k
H_k
Hk 正定,即
λ
g
k
T
H
k
−
1
g
k
>
0
\lambda g_k^{\mathrm{T}}H_k^{-1}g_k \gt 0
λgkTHk−1gk>0 ,因此有
f
(
x
(
k
+
1
)
)
<
f
(
x
(
k
)
)
f(x^{(k+1)})\lt f(x^{(k)})
f(x(k+1))<f(x(k)) ,也就是说
p
k
p_k
pk 时下降方向。
综合前面两个条件,我们在选择近似矩阵
G
k
G_k
Gk 时,首先要求每次迭代时
G
k
G_k
Gk 是正定的,其次
G
k
G_k
Gk 需要满足以下的拟牛顿条件:
G
k
+
1
y
k
=
δ
k
G_{k+1}y_{k}=\delta_k
Gk+1yk=δk
同时,按照拟牛顿条件,每次迭代可以选择更新矩阵
G
k
+
1
G_{k+1}
Gk+1 :
G
k
+
1
=
G
k
+
Δ
G
k
G_{k+1}=G_{k}+\Delta G_{k}
Gk+1=Gk+ΔGk
有多种选择近似矩阵的方法:
DFP 算法
DFP (Davidon-Fletcher-Powell)算法对
G
k
G_{k}
Gk 的更新为:
G
k
+
1
=
G
k
+
P
k
+
Q
k
G_{k+1}=G_{k}+P_k+Q_k
Gk+1=Gk+Pk+Qk
其中
P
k
P_k
Pk 和
Q
k
Q_k
Qk 是待定矩阵,此时:
G
k
+
1
y
k
=
G
k
y
k
+
P
k
y
k
+
Q
k
y
k
=
δ
k
G_{k+1}y_{k}=G_ky_k+P_ky_k+Q_ky_k=\delta_k
Gk+1yk=Gkyk+Pkyk+Qkyk=δk
为了满足拟牛顿条件,可以令:
P
k
y
k
=
δ
k
,
Q
k
y
k
=
−
G
k
y
k
P_ky_k=\delta_k,\quad Q_ky_k=-G_ky_k
Pkyk=δk,Qkyk=−Gkyk
书上给了一个例子:
P
k
=
δ
k
δ
k
T
δ
k
T
y
k
,
Q
k
=
−
G
k
y
k
y
k
T
G
k
y
k
T
G
k
y
k
P_k=\frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k},\quad Q_k=-\frac{G_ky_ky_k^{\mathrm{T}}G_k}{y_k^{\mathrm{T}}G_ky_k}
Pk=δkTykδkδkT,Qk=−ykTGkykGkykykTGk
可以证明(咱也不知道怎么证明),如果
G
0
G_0
G0 是正定的,则迭代过程中的每个矩阵
G
k
G_k
Gk 都是正定的。
算法:DFP 算法
- 输入:目标函数 f ( x ) f(x) f(x) ,梯度 g ( x ) = ∇ f ( x ) g(x)=\nabla f(x) g(x)=∇f(x) ,精度 ε \varepsilon ε ;
- 输出: f ( x ) f(x) f(x) 的极小点 x ∗ x^\ast x∗ ;
- 选定初始点 x ( 0 ) x^{(0)} x(0) ,取 G 0 G_0 G0 为正定对称矩阵,置 k = 0 k=0 k=0 ;
- 计算 g k = g ( x ( k ) ) g_k=g(x^{(k)}) gk=g(x(k)) ,若 ∥ g ( x ( k ) ) ∥ < ε \|g(x^{(k)})\| \lt \varepsilon ∥g(x(k))∥<ε ,则停止计算,返回近似解 x ∗ = x ( k ) x^\ast=x^{(k)} x∗=x(k) ,否则继续;
- 计算 p k = − G k g k p_k=-G_kg_k pk=−Gkgk ;
- 一维搜索:求 λ k \lambda_k λk 使得:
f ( x ( k ) + λ k p k ) = min λ ≥ 0 f ( x ( k ) + λ p k ) f(x^{(k)}+\lambda_kp_k)=\min_{\lambda \geq 0}f(x^{(k)}+\lambda p_k) f(x(k)+λkpk)=λ≥0minf(x(k)+λpk)
- 置 x ( k + 1 ) = x ( k ) + λ k p k x^{(k+1)}=x^{(k)}+\lambda_kp_k x(k+1)=x(k)+λkpk ;
- 计算 G k + 1 G_{k+1} Gk+1 ,置 k = k + 1 k=k+1 k=k+1 ,转 2;
G k + 1 = G k + δ k δ k T δ k T y k − G k y k y k T G k y k T G k y k G_{k+1}=G_{k}+\frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k}{y_k^{\mathrm{T}}G_ky_k} Gk+1=Gk+δkTykδkδkT−ykTGkykGkykykTGk
BFGS 算法
BFGS(Broyden-Fletcher-Goldfarb-Shanno)算法是最流行的拟牛顿算法。此时考虑使用
B
k
B_k
Bk 近似
H
k
H_k
Hk ,更新方法类似:
B
k
+
1
=
B
k
+
P
k
+
Q
k
B_{k+1}=B_{k}+P_{k}+Q_{k}
Bk+1=Bk+Pk+Qk
为了满足拟牛顿条件:
B
k
+
1
δ
k
=
B
k
δ
k
+
P
k
δ
k
+
Q
k
δ
k
=
y
k
B_{k+1}\delta_k=B_k\delta_k+P_k\delta_k+Q_k\delta_k=y_k
Bk+1δk=Bkδk+Pkδk+Qkδk=yk
取:
P
k
δ
k
=
y
k
,
Q
k
δ
k
=
−
B
k
δ
k
P_k\delta_k=y_k,\quad Q_k\delta_k=-B_k\delta_k
Pkδk=yk,Qkδk=−Bkδk
同样可令:
P
k
=
y
k
y
k
T
y
k
T
δ
k
,
Q
k
=
−
B
k
δ
k
δ
k
T
B
k
δ
k
T
B
k
δ
k
P_k=\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k},\quad Q_k=-\frac{B_k\delta_k\delta_k^{\mathrm{T}}B_k}{\delta_k^{\mathrm{T}}B_k\delta_k}
Pk=ykTδkykykT,Qk=−δkTBkδkBkδkδkTBk
算法:BFGS 算法
- 输入:目标函数 f ( x ) f(x) f(x) ,梯度 g ( x ) = ∇ f ( x ) g(x)=\nabla f(x) g(x)=∇f(x) ,精度 ε \varepsilon ε ;
- 输出: f ( x ) f(x) f(x) 的极小点 x ∗ x^\ast x∗ ;
- 选定初始点 x ( 0 ) x^{(0)} x(0) ,取 B 0 B_0 B0 为正定对称矩阵,置 k = 0 k=0 k=0 ;
- 计算 g k = g ( x ( k ) ) g_k=g(x^{(k)}) gk=g(x(k)) ,若 ∥ g ( x ( k ) ) ∥ < ε \|g(x^{(k)})\| \lt \varepsilon ∥g(x(k))∥<ε ,则停止计算,返回近似解 x ∗ = x ( k ) x^\ast=x^{(k)} x∗=x(k) ,否则继续;
- 由 B k p k = − g k B_kp_k=-g_k Bkpk=−gk 计算 p k p_k pk;
- 一维搜索:求 λ k \lambda_k λk 使得:
f ( x ( k ) + λ k p k ) = min λ ≥ 0 f ( x ( k ) + λ p k ) f(x^{(k)}+\lambda_kp_k)=\min_{\lambda \geq 0}f(x^{(k)}+\lambda p_k) f(x(k)+λkpk)=λ≥0minf(x(k)+λpk)
- 置 x ( k + 1 ) = x ( k ) + λ k p k x^{(k+1)}=x^{(k)}+\lambda_kp_k x(k+1)=x(k)+λkpk ;
- 计算 B k + 1 B_{k+1} Bk+1 ,置 k = k + 1 k=k+1 k=k+1 ,转 2;
B k + 1 = B k + y k y k T y k T δ k − B k δ k δ k T B k δ k T B k δ k B_{k+1}=B_{k}+\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}-\frac{B_k\delta_k\delta_k^{\mathrm{T}}B_k}{\delta_k^{\mathrm{T}}B_k\delta_k} Bk+1=Bk+ykTδkykykT−δkTBkδkBkδkδkTBk
Broyden 类算法
根据 BFGS 算法中
B
k
B_k
Bk 矩阵的迭代公式,我们可以得到对应的
G
k
G_k
Gk 矩阵。有
G
k
=
B
k
−
1
G_k=B_k^{-1}
Gk=Bk−1 ,
G
k
+
1
−
1
=
B
k
+
1
−
1
G_{k+1}^{-1}=B_{k+1}^{-1}
Gk+1−1=Bk+1−1 ;有:
G
k
+
1
−
1
=
G
k
−
1
+
y
k
y
k
T
y
k
T
δ
k
−
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
T
G
k
−
1
δ
k
⇒
G
k
+
1
=
(
G
k
−
1
+
y
k
y
k
T
y
k
T
δ
k
+
−
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
T
G
k
−
1
δ
k
)
−
1
\begin{aligned} G_{k+1}^{-1}=&\, G_{k}^{-1}+\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}-\frac{G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k} \\ \Rightarrow G_{k+1}=&\, \left( G_{k}^{-1}+\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}+\frac{-G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k} \right)^{-1} \\ \end{aligned}
Gk+1−1=⇒Gk+1=Gk−1+ykTδkykykT−δkTGk−1δkGk−1δkδkTGk−1(Gk−1+ykTδkykykT+δkTGk−1δk−Gk−1δkδkTGk−1)−1
需要用到 Shermax-Morrison 公式:
(
A
+
u
v
T
)
−
1
=
A
−
1
−
A
−
1
u
v
T
A
−
1
1
+
v
T
A
−
1
u
(A+uv^\mathrm{T})^{-1}=A^{-1}-\frac{A^{-1}uv^\mathrm{T}A^{-1}}{1+v^{\mathrm{T}}A^{-1}u}
(A+uvT)−1=A−1−1+vTA−1uA−1uvTA−1
其中
u
u
u 和
v
v
v 为
n
n
n 维向量,
A
A
A 为可逆矩阵;记:
A
=
G
k
−
1
+
y
k
y
k
T
y
k
T
δ
k
,
u
=
−
G
k
−
1
δ
k
,
v
T
=
δ
k
T
G
k
−
1
δ
k
T
G
k
−
1
δ
k
A=G_k^{-1}+\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k},\quad u=-G_k^{-1}\delta_k,\quad v^{\mathrm{T}}=\frac{\delta_k^{\mathrm{T}}G_k^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k}
A=Gk−1+ykTδkykykT,u=−Gk−1δk,vT=δkTGk−1δkδkTGk−1
则:
G
k
+
1
=
(
A
+
−
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
T
G
k
−
1
δ
k
)
−
1
=
A
−
1
+
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
T
G
k
−
1
δ
k
A
−
1
1
−
δ
k
T
G
k
−
1
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
=
A
−
1
+
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
A
−
1
δ
k
T
G
k
−
1
δ
k
−
δ
k
T
G
k
−
1
A
−
1
G
k
−
1
δ
k
\begin{aligned} G_{k+1} =&\, \left(A+\frac{-G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k}\right)^{-1} \\ =&\, A^{-1}+\frac{A^{-1}\frac{G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k}A^{-1}}{1-\frac{\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}G_k^{-1}\delta_{k}} {\delta_k^{\mathrm{T}}G_k^{-1}\delta_k}} \\ =&\, A^{-1}+\frac{A^{-1}G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k-\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}G_k^{-1}\delta_{k}} \end{aligned}
Gk+1===(A+δkTGk−1δk−Gk−1δkδkTGk−1)−1A−1+1−δkTGk−1δkδkTGk−1A−1Gk−1δkA−1δkTGk−1δkGk−1δkδkTGk−1A−1A−1+δkTGk−1δk−δkTGk−1A−1Gk−1δkA−1Gk−1δkδkTGk−1A−1
对
A
A
A 再次使用公式,此时
u
=
y
k
u=y_k
u=yk ,
v
T
=
y
k
T
y
k
T
δ
k
v^{\mathrm{T}}=\frac{y_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}
vT=ykTδkykT ,得到:
A
−
1
=
(
G
k
−
1
+
y
k
y
k
T
y
k
T
δ
k
)
−
1
=
G
k
−
G
k
y
k
y
k
T
y
k
T
δ
k
G
k
1
+
y
k
T
G
k
y
k
y
k
T
δ
k
=
G
k
−
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
A^{-1}=\left(G_k^{-1}+\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right)^{-1} =G_{k}-\frac{G_k\frac{y_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}G_k} {1+\frac{y_k^{\mathrm{T}}G_ky_k}{y_k^{\mathrm{T}}\delta_k}} =G_{k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k}
A−1=(Gk−1+ykTδkykykT)−1=Gk−1+ykTδkykTGkykGkykTδkykykTGk=Gk−ykTδk+ykTGkykGkykykTGk
将
A
−
1
A^{-1}
A−1 逐步代入
G
k
+
1
G_{k+1}
Gk+1 ,得到:
G
k
+
1
=
A
−
1
+
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
A
−
1
δ
k
T
G
k
−
1
δ
k
−
δ
k
T
G
k
−
1
A
−
1
G
k
−
1
δ
k
=
A
−
1
+
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
A
−
1
δ
k
T
G
k
−
1
δ
k
−
δ
k
T
G
k
−
1
(
G
k
−
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
G
k
−
1
δ
k
=
A
−
1
+
A
−
1
G
k
−
1
δ
k
δ
k
T
G
k
−
1
A
−
1
δ
k
T
y
k
y
k
T
δ
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
=
A
−
1
+
(
G
k
−
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
(
G
k
−
1
δ
k
δ
k
T
G
k
−
1
δ
k
T
y
k
y
k
T
δ
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
(
G
k
−
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
=
A
−
1
+
(
I
−
G
k
y
k
y
k
T
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
(
δ
k
δ
k
T
δ
k
T
y
k
y
k
T
δ
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
(
I
−
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
=
A
−
1
+
(
(
δ
k
δ
k
T
)
(
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
δ
k
T
y
k
y
k
T
δ
k
)
−
(
δ
k
δ
k
T
y
k
y
k
T
G
k
δ
k
T
y
k
y
k
T
δ
k
)
−
(
G
k
y
k
y
k
T
δ
k
δ
k
T
δ
k
T
y
k
y
k
T
δ
k
)
+
(
G
k
y
k
y
k
T
δ
k
δ
k
T
y
k
y
k
T
G
k
(
δ
k
T
y
k
y
k
T
δ
k
)
(
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
)
=
G
k
−
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
+
(
(
δ
k
δ
k
T
)
(
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
δ
k
T
y
k
y
k
T
δ
k
)
−
(
δ
k
δ
k
T
y
k
y
k
T
G
k
δ
k
T
y
k
y
k
T
δ
k
)
−
(
G
k
y
k
y
k
T
δ
k
δ
k
T
δ
k
T
y
k
y
k
T
δ
k
)
+
(
G
k
y
k
y
k
T
G
k
y
k
T
δ
k
+
y
k
T
G
k
y
k
)
=
G
k
+
(
(
δ
k
δ
k
T
)
(
y
k
T
δ
k
)
δ
k
T
y
k
y
k
T
δ
k
)
+
(
(
δ
k
δ
k
T
)
(
y
k
T
G
k
y
k
)
δ
k
T
y
k
y
k
T
δ
k
)
−
(
δ
k
y
k
T
G
k
y
k
T
δ
k
)
−
(
G
k
y
k
δ
k
T
δ
k
T
y
k
)
=
G
k
−
(
G
k
y
k
δ
k
T
δ
k
T
y
k
)
−
(
δ
k
y
k
T
G
k
y
k
T
δ
k
)
+
(
δ
k
(
y
k
T
G
k
y
k
)
δ
k
T
δ
k
T
y
k
y
k
T
δ
k
)
+
(
δ
k
δ
k
T
δ
k
T
y
k
)
=
G
k
(
I
−
y
k
δ
k
T
δ
k
T
y
k
)
−
(
δ
k
y
k
T
G
k
y
k
T
δ
k
)
(
I
−
y
k
δ
k
T
δ
k
T
y
k
)
+
(
δ
k
δ
k
T
δ
k
T
y
k
)
=
(
I
−
δ
k
y
k
T
y
k
T
δ
k
)
G
k
(
I
−
y
k
δ
k
T
δ
k
T
y
k
)
+
(
δ
k
δ
k
T
δ
k
T
y
k
)
=
(
I
−
δ
k
y
k
T
y
k
T
δ
k
)
G
k
(
I
−
δ
k
y
k
T
y
k
T
δ
k
)
T
+
(
δ
k
δ
k
T
δ
k
T
y
k
)
\begin{aligned} G_{k+1} =&\, A^{-1}+\frac{A^{-1}G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k-\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}G_k^{-1}\delta_{k}} \\ =&\, A^{-1}+\frac{A^{-1}G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}}{\delta_k^{\mathrm{T}}G_k^{-1}\delta_k-\delta_k^{\mathrm{T}}G_k^{-1} \left( G_{k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) G_k^{-1}\delta_{k}} \\ =&\, A^{-1}+\frac{A^{-1}G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}A^{-1}} {\frac{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k}{y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k}} \\ =&\, A^{-1}+ \left( G_{k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) \left( \frac{G_k^{-1}\delta_k\delta_k^{\mathrm{T}}G_k^{-1}} {\frac{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k}{y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k}} \right) \left( G_{k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) \\ =&\, A^{-1}+ \left( I-\frac{G_ky_ky_k^{\mathrm{T}}} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) \left( \frac{\delta_k\delta_k^{\mathrm{T}}} {\frac{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k}{y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k}} \right) \left( I-\frac{y_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) \\ =&\, A^{-1}+ \left( \frac{(\delta_k\delta_k^{\mathrm{T}})(y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k)} {\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{\delta_k\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}G_k}{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{G_ky_ky_k^{\mathrm{T}}\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) +\left( \frac{G_ky_ky_k^{\mathrm{T}}\delta_k\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}G_k} {(\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k)(y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k)} \right) \\ =&\, G_{k}-\frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \\ +&\,\left( \frac{(\delta_k\delta_k^{\mathrm{T}})(y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k)} {\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{\delta_k\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}G_k}{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{G_ky_ky_k^{\mathrm{T}}\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) +\left( \frac{G_ky_ky_k^{\mathrm{T}}G_k} {y_k^{\mathrm{T}}\delta_k+y_k^{\mathrm{T}}G_ky_k} \right) \\ =&\, G_{k}+ \left( \frac{(\delta_k\delta_k^{\mathrm{T}})(y_k^{\mathrm{T}}\delta_k)} {\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) +\left( \frac{(\delta_k\delta_k^{\mathrm{T}})(y_k^{\mathrm{T}}G_ky_k)} {\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{\delta_ky_k^{\mathrm{T}}G_k}{y_k^{\mathrm{T}}\delta_k} \right) -\left( \frac{G_ky_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) \\ =&\, G_{k} -\left( \frac{G_ky_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) -\left( \frac{\delta_ky_k^{\mathrm{T}}G_k}{y_k^{\mathrm{T}}\delta_k} \right) +\left( \frac{\delta_k(y_k^{\mathrm{T}}G_ky_k)\delta_k^{\mathrm{T}}} {\delta_k^{\mathrm{T}}y_ky_k^{\mathrm{T}}\delta_k} \right) +\left( \frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) \\ =&\, G_{k}\left( I- \frac{y_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) -\left( \frac{\delta_ky_k^{\mathrm{T}}G_k}{y_k^{\mathrm{T}}\delta_k} \right) \left( I-\frac{y_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) +\left( \frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) \\ =&\, \left( I- \frac{\delta_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right) G_k \left( I-\frac{y_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) +\left( \frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) \\ =&\, \left( I- \frac{\delta_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right) G_k \left( I- \frac{\delta_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right)^\mathrm{T} +\left( \frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k} \right) \\ \end{aligned}
Gk+1=======+=====A−1+δkTGk−1δk−δkTGk−1A−1Gk−1δkA−1Gk−1δkδkTGk−1A−1A−1+δkTGk−1δk−δkTGk−1(Gk−ykTδk+ykTGkykGkykykTGk)Gk−1δkA−1Gk−1δkδkTGk−1A−1A−1+ykTδk+ykTGkykδkTykykTδkA−1Gk−1δkδkTGk−1A−1A−1+(Gk−ykTδk+ykTGkykGkykykTGk)
ykTδk+ykTGkykδkTykykTδkGk−1δkδkTGk−1
(Gk−ykTδk+ykTGkykGkykykTGk)A−1+(I−ykTδk+ykTGkykGkykykT)
ykTδk+ykTGkykδkTykykTδkδkδkT
(I−ykTδk+ykTGkykykykTGk)A−1+(δkTykykTδk(δkδkT)(ykTδk+ykTGkyk))−(δkTykykTδkδkδkTykykTGk)−(δkTykykTδkGkykykTδkδkT)+((δkTykykTδk)(ykTδk+ykTGkyk)GkykykTδkδkTykykTGk)Gk−ykTδk+ykTGkykGkykykTGk(δkTykykTδk(δkδkT)(ykTδk+ykTGkyk))−(δkTykykTδkδkδkTykykTGk)−(δkTykykTδkGkykykTδkδkT)+(ykTδk+ykTGkykGkykykTGk)Gk+(δkTykykTδk(δkδkT)(ykTδk))+(δkTykykTδk(δkδkT)(ykTGkyk))−(ykTδkδkykTGk)−(δkTykGkykδkT)Gk−(δkTykGkykδkT)−(ykTδkδkykTGk)+(δkTykykTδkδk(ykTGkyk)δkT)+(δkTykδkδkT)Gk(I−δkTykykδkT)−(ykTδkδkykTGk)(I−δkTykykδkT)+(δkTykδkδkT)(I−ykTδkδkykT)Gk(I−δkTykykδkT)+(δkTykδkδkT)(I−ykTδkδkykT)Gk(I−ykTδkδkykT)T+(δkTykδkδkT)
因此得到 BFGS 算法的
G
k
G_k
Gk 的迭代公式:
G
k
+
1
BFGS
=
(
I
−
δ
k
y
k
T
y
k
T
δ
k
)
G
k
BFGS
(
I
−
δ
k
y
k
T
y
k
T
δ
k
)
T
+
δ
k
δ
k
T
δ
k
T
y
k
G_{k+1}^{\text{BFGS}}= \left( I- \frac{\delta_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right) G_k^{\text{BFGS}} \left( I- \frac{\delta_ky_k^{\mathrm{T}}}{y_k^{\mathrm{T}}\delta_k}\right)^\mathrm{T} +\frac{\delta_k\delta_k^{\mathrm{T}}}{\delta_k^{\mathrm{T}}y_k}
Gk+1BFGS=(I−ykTδkδkykT)GkBFGS(I−ykTδkδkykT)T+δkTykδkδkT
由 DFP 算法得到的
G
k
G_k
Gk 记作
G
k
DFP
G_k^{\text{DFP}}
GkDFP ,可知二者的线性组合也满足拟牛顿条件式,而且是正定的:
G
k
+
1
=
α
G
k
DFP
+
(
1
−
α
)
G
k
BFGS
,
0
≤
α
≤
1
G_{k+1}=\alpha G_k^{\text{DFP}}+(1-\alpha)G_k^{\text{BFGS}},\quad 0\leq\alpha\leq 1
Gk+1=αGkDFP+(1−α)GkBFGS,0≤α≤1
这样就得到了一类拟牛顿算法,称为 Broyden 类算法。