连续函数的最优化方法-LM法
1、介绍
列文伯格(1944)和马夸尔特(1963)先后对高斯牛顿法进行了改进，求解过程中引入了阻尼因子。
将公式 36 的无约束最小二乘问题转变为公式 44 有约束最小二乘问题，其中，
1
2
×
(
∥
D
Δ
X
k
∥
2
−
μ
)
⩽
0
\frac{1}{2} \times ( \begin{Vmatrix} D\Delta X_k \end{Vmatrix}^2 -\mu )\leqslant 0
21×(∥∥DΔXk∥∥2−μ)⩽0 是信赖区间对应的条件，
D
D
D 是系数矩阵(列文伯格方法和马夸尔特方法的主要区别就是系数矩阵
D
D
D 不同)，
μ
\mu
μ 是信赖 半径。
F
(
X
k
+
Δ
X
k
)
≈
1
2
∥
L
(
X
k
)
+
J
(
X
k
)
⏟
L
T
Δ
X
k
∥
2
,
s
.
t
.
(
1
2
×
(
∥
D
Δ
X
k
∥
2
−
μ
)
⩽
0
)
(
公
式
44
)
F(X_k+\Delta X_k)\approx \frac{1}{2} \begin{Vmatrix} \\ L(X_k)+\underbrace{J(X_k)}_{L}{^T} \Delta X_k \end{Vmatrix}^2 \qquad ,s.t.( \frac{1}{2} \times ( \begin{Vmatrix} D\Delta X_k \end{Vmatrix}^2 -\mu )\leqslant 0 ) \qquad (公式44)
F(Xk+ΔXk)≈21∥∥∥∥∥∥L(Xk)+L
J(Xk)TΔXk∥∥∥∥∥∥2,s.t.(21×(∥∥DΔXk∥∥2−μ)⩽0)(公式44)
2、数学原理
将公式 44 的约束条件引入，可以定义拉格朗日函数
L
a
g
(
X
)
Lag(X)
Lag(X)：
L
a
g
(
Δ
X
)
=
d
e
f
1
2
∥
L
(
X
k
)
+
J
(
X
k
)
⏟
L
T
Δ
X
k
∥
2
+
λ
k
×
1
2
×
(
∥
D
Δ
X
k
∥
2
−
μ
)
,
s
.
t
.
(
λ
k
⩾
0
)
(
公
式
45
)
Lag(\Delta X)\stackrel{\mathrm{def}}{=} \frac{1}{2} \begin{Vmatrix} \\ L(X_k)+\underbrace{J(X_k)}_{L}{^T} \Delta X_k \end{Vmatrix}^2 +{\lambda}_k \times \frac{1}{2} \times ( \begin{Vmatrix} D\Delta X_k \end{Vmatrix}^2 -\mu ) \qquad ,s.t.( {\lambda}_k \geqslant 0 ) \qquad (公式45)
Lag(ΔX)=def21∥∥∥∥∥∥L(Xk)+L
J(Xk)TΔXk∥∥∥∥∥∥2+λk×21×(∥∥DΔXk∥∥2−μ),s.t.(λk⩾0)(公式45)
公式45中
J
(
X
k
)
⏟
L
\underbrace{J(X_k)}_{L}
L
J(Xk)是残差函数
L
(
X
k
)
L(X_k)
L(Xk)的雅可比矩阵。当
∥
D
Δ
X
k
∥
2
⩾
0
\begin{Vmatrix} D\Delta X_k \end{Vmatrix}^2 \geqslant 0
∥∥DΔXk∥∥2⩾0的时候，在优化目标函数
F
(
X
)
F(X)
F(X)中引入阻尼因子
λ
k
{\lambda}_k
λk作为惩罚项。
该表达式中
L
a
g
(
Δ
X
k
)
Lag(\Delta X_k)
Lag(ΔXk)、
L
(
X
k
)
L(X_k)
L(Xk)是一个常数，
J
(
X
k
)
⏟
L
\underbrace{J(X_k)}_{L}
L
J(Xk)、
D
D
D是一个常数矩阵，
Δ
X
k
\Delta X_k
ΔXk是一个变量矩阵，即函数
L
a
g
(
Δ
X
k
)
Lag(\Delta X_k)
Lag(ΔXk)是以
Δ
X
k
\Delta X_k
ΔXk为自变量的二次函数。综上所述，当
L
a
g
(
Δ
X
k
)
Lag(\Delta X_k)
Lag(ΔXk)一阶导数为0的时候，函数
L
a
g
(
Δ
X
k
)
Lag(\Delta X_k)
Lag(ΔXk)取得极值，可推得：
L
a
g
Δ
X
k
′
(
Δ
X
k
)
=
0
(
公
式
46
)
Lag'_{\Delta X_k}(\Delta X_k)=0 \qquad (公式46)
LagΔXk′(ΔXk)=0(公式46)
由公式 45、公式 46 可推得：
0
=
λ
k
D
T
D
Δ
X
k
+
L
(
X
k
)
J
(
X
k
)
⏟
L
+
J
(
X
k
)
⏟
L
J
(
X
k
)
⏟
L
T
Δ
X
k
(
公
式
47
)
0= {\lambda}_k D^TD \Delta X_k + L(X_k) \underbrace{J(X_k)}_{L} + \underbrace{J(X_k)}_{L} \underbrace{J(X_k)}_{L}{^T} \Delta X_k \qquad (公式47)
0=λkDTDΔXk+L(Xk)L
J(Xk)+L
J(Xk)L
J(Xk)TΔXk(公式47)
即：
0
=
λ
k
D
T
D
Δ
X
k
+
∑
i
=
1
m
(
L
i
(
X
k
)
J
(
X
k
)
⏟
L
i
+
J
(
X
k
)
⏟
L
i
J
(
X
k
)
⏟
L
i
T
Δ
X
k
)
(
公
式
48
)
0= {\lambda}_k D^TD \Delta X_k+ \sum_{i=1}^m ( L_i(X_k) \underbrace{J(X_k)}_{L_i} + \underbrace{J(X_k)}_{L_i} \underbrace{J(X_k)}_{L_i}{^T} \Delta X_k ) \qquad (公式48)
0=λkDTDΔXk+i=1∑m(Li(Xk)Li
J(Xk)+Li
J(Xk)Li
J(Xk)TΔXk)(公式48)
由公式 47 可推得:
0
=
L
(
X
k
)
J
(
X
k
)
⏟
L
+
(
J
(
X
k
)
⏟
L
J
(
X
k
)
⏟
L
T
+
λ
k
D
T
D
)
Δ
X
k
(
公
式
49
)
0= L(X_k) \underbrace{J(X_k)}_{L} + ( \underbrace{J(X_k)}_{L} \underbrace{J(X_k)}_{L}{^T} +{\lambda}_k D^TD ) \Delta X_k \qquad (公式49)
0=L(Xk)L
J(Xk)+(L
J(Xk)L
J(Xk)T+λkDTD)ΔXk(公式49)
设：由公式49结合公式26结构形式，可近似推得函数
F
(
X
)
F(X)
F(X)的黑塞矩阵
H
(
X
k
)
⏟
F
\underbrace{H(X_k)}_{F}
F
H(Xk)和雅克比矩阵
J
(
X
k
)
⏟
F
\underbrace{J(X_k)}_{F}
F
J(Xk)：
H
(
X
k
)
⏟
F
≈
d
e
f
J
(
X
k
)
⏟
L
J
(
X
k
)
⏟
L
T
(
公
式
50
)
\underbrace{H(X_k)}_{F} \stackrel{\mathrm{def}}{\approx} \underbrace{J(X_k)}_{L} \underbrace{J(X_k)}_{L}{^T} \qquad (公式50)
F
H(Xk)≈defL
J(Xk)L
J(Xk)T(公式50)
J
(
X
k
)
⏟
F
≈
d
e
f
L
(
X
k
)
J
(
X
k
)
⏟
L
(
公
式
51
)
\underbrace{J(X_k)}_{F} \stackrel{\mathrm{def}}{\approx} L(X_k) \underbrace{J(X_k)}_{L} \qquad (公式51)
F
J(Xk)≈defL(Xk)L
J(Xk)(公式51)
由公式 49、公式 50、公式 51 可推得：
Δ
X
k
=
−
(
H
(
X
k
)
⏟
F
+
λ
k
D
T
D
)
−
1
J
(
X
k
)
⏟
F
(
公
式
52
)
\Delta X_k=- ( {\underbrace{H(X_k)}_{F}} + {\lambda}_k D^TD )^{-1} \underbrace{J(X_k)}_{F} \qquad (公式52)
ΔXk=−(F
H(Xk)+λkDTD)−1F
J(Xk)(公式52)
由公式 52 可推得目标函数
F
(
X
)
F(X)
F(X)的最优化迭代公式：
X
k
+
1
=
d
e
f
X
k
−
(
H
(
X
k
)
⏟
F
+
λ
k
D
T
D
)
−
1
J
(
X
k
)
⏟
F
(
公
式
53
)
X_{k+1}\stackrel{\mathrm{def}}{=} X_{k} -( {\underbrace{H(X_k)}_{F}} + {\lambda}_k D^TD )^{-1} \underbrace{J(X_k)}_{F} \qquad (公式53)
Xk+1=defXk−(F
H(Xk)+λkDTD)−1F
J(Xk)(公式53)
当阻尼因子
λ
k
>
0
{\lambda}_k>0
λk>0的时候，保证
H
(
X
k
)
⏟
F
+
λ
k
D
T
D
{\underbrace{H(X_k)}_{F}} + {\lambda}_k D^TD
F
H(Xk)+λkDTD正定，迭代朝着下降方向进行；当
λ
k
{\lambda}_k
λk趋近于 0 的时候，退化为高斯牛顿算法；当
λ
k
{\lambda}_k
λk 趋近于正无穷大的时候，退化为最速下降算法；
3、阻尼因子更新策略
目标函数
F
(
X
)
F(X)
F(X)在
X
=
X
k
X=X_k
X=Xk处的不含皮亚诺余项的二阶泰勒公式如下：
F
(
X
k
+
Δ
X
k
)
≈
F
(
X
k
)
+
J
(
X
k
)
⏟
F
T
Δ
X
k
＋
1
2
Δ
X
k
T
H
(
X
k
)
⏟
F
Δ
X
k
(
公
式
54
)
F(X_k+\Delta X_k)\approx F(X_k)+ \underbrace{J(X_k)}_{F}{^T} \Delta X_k ＋\frac{1}{2} {\Delta X_k}^T \underbrace{H(X_k)}_{F} \Delta X_k \qquad (公式54)
F(Xk+ΔXk)≈F(Xk)+F
J(Xk)TΔXk＋21ΔXkTF
H(Xk)ΔXk(公式54)
根据公式 50、公式 51 的近似关系和公式 54 可推得：
U
(
Δ
X
k
)
=
d
e
f
F
(
X
k
+
Δ
X
k
)
≈
F
(
X
k
)
+
J
(
X
k
)
⏟
L
T
L
(
X
k
)
Δ
X
k
＋
1
2
Δ
X
k
T
J
(
X
k
)
⏟
L
T
J
(
X
k
)
⏟
L
Δ
X
k
(
公
式
55
)
U(\Delta X_k)\stackrel{\mathrm{def}}{=} F(X_k+\Delta X_k)\approx F(X_k)+ {\underbrace{J(X_k)}_{L}}^T L(X_k) \Delta X_k ＋\frac{1}{2} {\Delta X_k}^T {\underbrace{J(X_k)}_{L}}^T {\underbrace{J(X_k)}_{L}} \Delta X_k \qquad (公式55)
U(ΔXk)=defF(Xk+ΔXk)≈F(Xk)+L
J(Xk)TL(Xk)ΔXk＋21ΔXkTL
J(Xk)TL
J(Xk)ΔXk(公式55)
阻尼因子
λ
k
{\lambda}_k
λk的初始值
λ
0
{\lambda}_0
λ0的选取，
H
0
⏟
F
≈
d
e
f
J
(
X
0
)
⏟
L
J
(
X
0
)
⏟
L
T
\underbrace{H_0}_{F} \stackrel{\mathrm{def}}{\approx} \underbrace{J(X_0)}_{L} \underbrace{J(X_0)}_{L}{^T}
F
H0≈defL
J(X0)L
J(X0)T对角线上面的数
h
i
i
h_{ii}
hii有关。
λ
0
=
d
e
f
τ
×
m
a
x
{
h
i
i
}
,
s
.
t
.
(
τ
∈
[
1
0
−
8
,
1
]
)
(
公
式
56
)
{\lambda}_0\stackrel{\mathrm{def}}{=} \tau \times max\begin{Bmatrix} h_{ii} \end{Bmatrix} \qquad ,s.t.( \tau \in[10^{-8},1] ) \qquad (公式56)
λ0=defτ×max{hii},s.t.(τ∈[10−8,1])(公式56)
而阻尼因子
λ
k
{\lambda}_k
λk 的更新可以由增益比
β
k
{\beta}_k
βk 来定量分析：
β
k
=
d
e
f
F
(
X
k
)
−
F
(
X
k
+
Δ
X
k
)
U
(
0
)
−
U
(
Δ
X
k
)
(
公
式
57
)
{\beta}_k\stackrel{\mathrm{def}}{=} \frac{F(X_k)-F(X_k+\Delta X_k)}{U(0)-U(\Delta X_k)} \qquad (公式57)
βk=defU(0)−U(ΔXk)F(Xk)−F(Xk+ΔXk)(公式57)
其中公式57的分子为目标函数
F
(
X
k
)
F(X_k)
F(Xk) 的值在步长
Δ
X
k
\Delta X_k
ΔXk 下的实际变化，分母为目标函数
F
(
X
k
)
F(X_k)
F(Xk) 的值二阶近 似的变化，可得出：
如果
β
k
{\beta}_k
βk 的值越大，那么目标函数的二阶近似变化对实际变化效果越好，可以缩小
λ
k
{\lambda}_k
λk 的值接近高斯牛顿算法;
如果
β
k
{\beta}_k
βk 的值越小，那么目标函数的二阶近似变化对实际变化效果越差，可以增大
λ
k
{\lambda}_k
λk 的值接近最速下降算法;
阻尼因子
λ
k
{\lambda}_k
λk 的
Nielsen
\text {Nielsen}
Nielsen (1991)更新策略如下：
如
果
(
β
k
>
0
)
:
λ
k
=
d
e
f
λ
k
−
1
×
m
a
x
{
1
3
，
1
−
(
2
β
k
−
1
)
3
}
v
k
=
d
e
f
2
否
则
(
β
k
⩽
0
)
:
λ
k
=
d
e
f
λ
k
−
1
×
v
k
−
1
v
k
=
2
×
v
k
−
1
(
公
式
58
)
\begin{aligned} 如果(\beta_k>0): & \qquad \\ &\qquad {\lambda}_k\stackrel{\mathrm{def}}{=} {\lambda}_{k-1} \times max\begin{Bmatrix} \frac{1}{3}，1-(2\beta_k-1)^3 \end{Bmatrix}\\ &\qquad v_k\stackrel{\mathrm{def}}{=}2\\ 否则(\beta_k \leqslant 0): &\\ &\qquad {\lambda}_k\stackrel{\mathrm{def}}{=} {\lambda}_{k-1} \times v_{k-1} \\ &\qquad v_k=2 \times v_{k-1} \end{aligned} \qquad (公式58)
如果(βk>0):否则(βk⩽0):λk=defλk−1×max{31，1−(2βk−1)3}vk=def2λk=defλk−1×vk−1vk=2×vk−1(公式58)
4、列文伯格方法
在列文伯格方法中，系数矩阵 D D D 是单位矩阵 I I I ，信赖区域是一个球形状。
5、马夸尔特方法
在马夸尔特方法中，系数矩阵 D D D 是 H ( X k ) ⏟ F ≈ d e f J ( X k ) ⏟ L J ( X k ) ⏟ L T \underbrace{H(X_k)}_{F} \stackrel{\mathrm{def}}{\approx} \underbrace{J(X_k)}_{L} \underbrace{J(X_k)}_{L}{^T} F H(Xk)≈defL J(Xk)L J(Xk)T 对角元素的平方根，信赖区域是一个椭球形状。
1163
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
