Optimizing Note
本文档介绍了一些常用的优化方法,包括线性最小二乘及一些非线性的优化方法,如梯度下降法,牛顿法,高斯牛顿法,列文伯格-马尔夸特法,迭代非线性最小二乘法等。
1. Linear Optimazation
设有一线性系统
y
−
A
x
=
e
∼
N
(
0
,
Σ
)
\boldsymbol{y} - \boldsymbol{A}\boldsymbol{x} = \boldsymbol{e} \sim \mathcal{N}(\boldsymbol{0},\boldsymbol{\Sigma})
y−Ax=e∼N(0,Σ) 其中
e
e
e 表示误差,服从正太分布,
Σ
\boldsymbol{\Sigma}
Σ 为协方差矩阵,表示权重。
构造优化方程,
x
=
arg min
e
T
Σ
−
1
e
\boldsymbol{x}=\argmin \boldsymbol{e}^T\Sigma^{-1}\boldsymbol{e}
x=argmineTΣ−1e解得,
x
=
(
A
T
Σ
−
1
A
)
−
1
A
T
Σ
−
1
y
\boldsymbol{x} = (\boldsymbol{A}^T\Sigma^{-1}\boldsymbol{A})^{-1}\boldsymbol{A}^T\Sigma^{-1}\boldsymbol{y}
x=(ATΣ−1A)−1ATΣ−1y
2. NonLinear Optimazation
假设有优化问题:
min
x
F
(
x
)
=
1
2
∥
f
(
x
)
∥
2
2
\min_{x} F(\boldsymbol{x}) = \frac{1}{2}\|f(\boldsymbol{x})\|_2^2
xminF(x)=21∥f(x)∥22其中,
x
∈
R
n
x \in \mathbb{R}^{n}
x∈Rn,
f
f
f 是任意函数
f
(
x
)
:
R
n
↦
R
f(\boldsymbol{x}):\mathbb{R}^n \mapsto \mathbb{R}
f(x):Rn↦R。这里系数
1
2
\frac{1}{2}
21 不影响结果。
要求解该问题,只需求找到其极值点,即求解:
d
F
d
x
=
0
\frac{\mathrm{d}F}{\mathrm{d}\boldsymbol{x}}=0
dxdF=0
在大多数情况下,这个方程的求解起来较为困难,尝用迭代的方式求解。常常可以根据如下的流程求解:
- 给定初值 x 0 \boldsymbol{x_0} x0;
- 对第 k k k 次迭代,寻找一增量 Δ x k \Delta \boldsymbol{x}_k Δxk, 使得 ∥ f ( x + Δ x ) ∥ 2 2 \|f(\boldsymbol{x}+\Delta \boldsymbol{x})\|^2_2 ∥f(x+Δx)∥22 达到最小值;
- 若 Δ x k \Delta \boldsymbol{x}_k Δxk 足够小,则停止迭代;
- 否则,令 x k + 1 = x k + Δ x \boldsymbol{x}_{k+1}=\boldsymbol{x}_k+\Delta \boldsymbol{x} xk+1=xk+Δx ,返回第二步
以下是几种不同的寻找 Δ x \Delta{\boldsymbol{x}} Δx 的方法。
2.1. Newton Algorithm (梯度法)
对
F
(
x
)
F(\boldsymbol{x})
F(x) 进行泰勒展开:
F
(
x
k
+
Δ
x
k
)
≈
F
(
x
k
)
+
J
(
x
k
)
T
Δ
x
k
+
1
2
Δ
x
k
T
H
(
x
k
)
Δ
x
k
F(\boldsymbol{x}_k+\Delta{\boldsymbol{x}_k}) \approx F(\boldsymbol{x}_k)+\boldsymbol{J}(\boldsymbol{x}_k)^T\Delta{\boldsymbol{x}_k}+\frac{1}{2}\Delta{\boldsymbol{x}_k}^T \boldsymbol{H}(\boldsymbol{x}_k)\Delta{\boldsymbol{x}_k}
F(xk+Δxk)≈F(xk)+J(xk)TΔxk+21ΔxkTH(xk)Δxk其中,
J
(
x
)
\boldsymbol{J}(\boldsymbol{x})
J(x) 为雅可比(Jacobian)矩阵,是
F
(
x
)
F(\boldsymbol{x})
F(x)关于
x
x
x 的一阶导数;
H
(
x
)
\boldsymbol{H}(\boldsymbol{x})
H(x) 是二阶导数,称为海塞(Hessian)矩阵。
此时,优化问题可写成
Δ
x
=
arg min
Δ
x
(
F
(
x
)
+
J
(
x
)
T
Δ
x
+
1
2
Δ
x
T
H
(
x
)
Δ
x
)
\Delta{\boldsymbol{x}}=\argmin_{\Delta \boldsymbol{x}}(F(\boldsymbol{x})+\boldsymbol{J}(\boldsymbol{x})^T\Delta\boldsymbol{x}+\frac{1}{2}\Delta\boldsymbol{x}^T\boldsymbol{H}(\boldsymbol{x})\Delta\boldsymbol{x})
Δx=Δxargmin(F(x)+J(x)TΔx+21ΔxTH(x)Δx)在这个优化问题中,
x
\boldsymbol{x}
x 是已知量,求解时是对
Δ
x
\Delta\boldsymbol{x}
Δx求导。
2.1.1. 最速下降法(梯度下降法)
若保留一阶项,可让增量取梯度的反方向,可保证函数下降,有,
Δ
x
=
−
λ
J
(
x
)
\Delta\boldsymbol{x}=-\lambda\boldsymbol{J}(\boldsymbol{x})
Δx=−λJ(x)其中,
λ
\lambda
λ 为步长,在机器学习中称为学习率。
这种方法简单,但迭代次数较多。
2.1.2. 牛顿法(Newton Algorithm)
若保留二阶项,则对
Δ
x
\Delta\boldsymbol{x}
Δx 求导,有,
Δ
x
=
−
J
H
−
1
\Delta \boldsymbol{x}=-\boldsymbol{J}\boldsymbol{H}^{-1}
Δx=−JH−1
该方法称为牛顿法。这种方法相对最速下降法迭代次数较少,但矩阵
H
\boldsymbol{H}
H的计算往往伴随着巨大的计算量。
2.2. Gauss-Newton Equation (or Normal Equation)
高斯牛顿法直接对
f
(
x
)
f(\boldsymbol{x})
f(x) 进行泰勒展开,得如下优化方程:
Δ
x
=
arg min
Δ
x
1
2
∥
f
(
x
)
+
J
(
x
)
T
Δ
x
∥
2
2
\Delta{\boldsymbol{x}}=\argmin_{\Delta \boldsymbol{x}}\frac{1}{2}\|f(\boldsymbol{x})+\boldsymbol{J}(\boldsymbol{x})^T\Delta\boldsymbol{x}\|_2^2
Δx=Δxargmin21∥f(x)+J(x)TΔx∥22其中
J
(
x
)
\boldsymbol{J}(\boldsymbol{x})
J(x) 是
f
(
x
)
f(\boldsymbol{x})
f(x) 的导数,不是
F
(
x
)
F(\boldsymbol{x})
F(x) 的导数。
求上式关于
Δ
x
\Delta\boldsymbol{x}
Δx 的导数,令其等于零,则有
J
(
x
)
J
T
(
x
)
Δ
x
=
−
J
(
x
)
f
(
x
)
\boldsymbol{J}(\boldsymbol{x})\boldsymbol{J}^T(\boldsymbol{x})\Delta \boldsymbol{x}=-\boldsymbol{J}(x)f(\boldsymbol{x})
J(x)JT(x)Δx=−J(x)f(x)该方程称高斯牛顿方程(Gauss-Newton Equation) 或 正规方程(Normal Equation)。
求解得,
Δ
x
=
−
(
J
J
T
)
−
1
J
f
\Delta \boldsymbol{x}=-(\boldsymbol{J}\boldsymbol{J}^T)^{-1}\boldsymbol{J}f
Δx=−(JJT)−1Jf该方法对
J
J
T
\boldsymbol{J}\boldsymbol{J}^T
JJT 有一定的要求,如果
J
J
T
\boldsymbol{J}\boldsymbol{J}^T
JJT 为奇异矩阵或者病态矩阵,则可能导致算法不收敛。
2.3. Levenberg-Marquardt Algorithm
列文伯格-马夸尔特法引入了信赖区域(Trust Region) 的概念,构造了如下的优化方程:
Δ
x
=
arg min
Δ
x
1
2
∥
f
(
x
)
+
J
(
x
)
T
Δ
x
∥
2
2
,
s
.
t
.
∥
D
Δ
x
∥
≤
μ
\Delta{\boldsymbol{x}}=\argmin_{\Delta \boldsymbol{x}}\frac{1}{2}\|f(\boldsymbol{x})+\boldsymbol{J}(\boldsymbol{x})^T\Delta\boldsymbol{x}\|_2^2 ,\quad s.t.\|\boldsymbol{D}\Delta \boldsymbol{x}\| \leq \mu
Δx=Δxargmin21∥f(x)+J(x)TΔx∥22,s.t.∥DΔx∥≤μ其中,
μ
\mu
μ 为信赖区域的半径,
D
\boldsymbol{D}
D为系数矩阵,是一个非负数对角阵。
该方程可用拉格朗日乘数法求解,构造如下拉格朗日函数,
L
(
Δ
x
,
λ
)
=
1
2
∥
f
(
x
)
+
J
(
x
)
T
Δ
x
∥
2
2
+
λ
2
(
∥
D
Δ
x
∥
2
2
−
μ
)
\mathcal{L}(\Delta\boldsymbol{x},\lambda)=\frac{1}{2}\|f(\boldsymbol{x})+\boldsymbol{J}(\boldsymbol{x})^T\Delta\boldsymbol{x}\|_2^2+\frac{\lambda}{2}(\|\boldsymbol{D}\Delta\boldsymbol{x}\|_2^2-\mu)
L(Δx,λ)=21∥f(x)+J(x)TΔx∥22+2λ(∥DΔx∥22−μ)其中,
λ
\lambda
λ 为拉格朗日乘子。对上式求导,可得
Δ
x
=
−
(
J
J
T
+
λ
D
T
D
)
−
1
J
f
\Delta \boldsymbol{x} = -(\boldsymbol{J}\boldsymbol{J}^T+\lambda\boldsymbol{D}^T\boldsymbol{D})^{-1}\boldsymbol{J}f
Δx=−(JJT+λDTD)−1Jf若取
D
=
I
\boldsymbol{D}=\boldsymbol{I}
D=I,则有
Δ
x
=
−
(
J
J
T
+
λ
I
)
−
1
J
f
\Delta \boldsymbol{x} = -(\boldsymbol{J}\boldsymbol{J}^T+\lambda\boldsymbol{I})^{-1}\boldsymbol{J}f
Δx=−(JJT+λI)−1Jf可以看到,当
λ
\lambda
λ 较大时,列文伯格-马夸尔特法接近一阶梯度下降法,当
λ
\lambda
λ 较小时,列文伯格-马夸尔特法接近高斯牛顿法。
除此之外,列文伯格-马夸尔特法还对迭代流程进行了优化,引入了一个指标
ρ
\rho
ρ 来衡量局部的近似程度。定义如下,
ρ
=
f
(
x
+
Δ
x
)
−
f
(
x
)
J
(
x
)
T
Δ
x
\rho = \frac{f(\boldsymbol{x+\Delta x})-f(\boldsymbol{x})}{\boldsymbol{J}(\boldsymbol{x})^T \Delta \boldsymbol{x}}
ρ=J(x)TΔxf(x+Δx)−f(x)可以看到,分子为函数实际下降的值,分母为近似下降的值,
ρ
\rho
ρ 越接近 1 ,说明近似的越好;若
ρ
\rho
ρ 太大,则应缩小范围;若
ρ
\rho
ρ 太小,则应放大范围。列文伯格-马夸尔特法的具体流程如下:
- 给定初值 x 0 \boldsymbol{x_0} x0,系数矩阵 D \boldsymbol{D} D 和信赖区域半径 μ \mu μ;
- 对第 k k k 次迭代,寻找满足优化方程的增量 Δ x k \Delta \boldsymbol{x}_k Δxk
- 若 Δ x k \Delta \boldsymbol{x}_k Δxk 足够小,则停止迭代;
- 否则,计算 ρ \rho ρ
- 若 ρ > 3 4 \rho > \frac{3}{4} ρ>43, 则设置 μ = 2 μ \mu=2\mu μ=2μ,进入步骤 2
- 若 ρ < 1 4 \rho < \frac{1}{4} ρ<41, 则设置 μ = 1 2 μ \mu=\frac{1}{2}\mu μ=21μ,进入步骤 2
- 若 ρ \rho ρ 大小合适,则令 x k + 1 = x k + Δ x k \boldsymbol{x}_{k+1}=\boldsymbol{x}_k+\Delta \boldsymbol{x}_k xk+1=xk+Δxk ,返回步骤 2
一般而言,列文伯格-马夸尔特法要比高斯牛顿法更加稳定,受“病态”问题影响较小。
2.4. Iterative Nonlinear Least Square Method(迭代非线性最小二乘法)
迭代最小二乘法在测绘领域中使用的非常广泛。在测绘中,我们把 ( x + Δ x ) (\boldsymbol{x}+\Delta \boldsymbol{x}) (x+Δx) 看成是参数的真实值,把 x \boldsymbol{x} x 看成是参数的近似值,则 Δ x \Delta \boldsymbol{x} Δx 为真实值与观测值之差,即参数的估计误差。一般情况下,我们无法直接获得想要的参数 x \boldsymbol{x} x,而只能观测到参数 x \boldsymbol{x} x 的函数 f ( x ) f(\boldsymbol{x}) f(x) ,我们称之为观测值。
容易知道,如果我们参数估计的够准,即
Δ
x
→
0
\Delta x \to 0
Δx→0,则有,
f
(
x
+
Δ
x
)
→
f
(
x
)
+
J
(
x
)
T
Δ
x
f(\boldsymbol{x}+\Delta \boldsymbol{x}) \to f(\boldsymbol{x})+\boldsymbol{J}(x)^T\Delta\boldsymbol{x}
f(x+Δx)→f(x)+J(x)TΔx,因此可以构造如下的优化方程:
Δ
x
=
arg min
Δ
x
1
2
∥
f
(
x
)
+
J
(
x
)
T
Δ
x
−
f
(
x
+
Δ
x
)
∥
2
2
\Delta{\boldsymbol{x}}=\argmin_{\Delta \boldsymbol{x}}\frac{1}{2}\|f(\boldsymbol{x})+\boldsymbol{J}(\boldsymbol{x})^T\Delta\boldsymbol{x}-f(\boldsymbol{x}+\Delta\boldsymbol{x})\|_2^2
Δx=Δxargmin21∥f(x)+J(x)TΔx−f(x+Δx)∥22求解可得,
Δ
x
=
(
J
T
J
)
−
1
J
T
L
\Delta \boldsymbol{x} = (\boldsymbol{J}^T\boldsymbol{J})^{-1}\boldsymbol{J}^T\boldsymbol{L}
Δx=(JTJ)−1JTL其中
L
=
f
(
x
+
Δ
x
)
−
f
(
x
)
\boldsymbol{L} = f(\boldsymbol{x}+\Delta\boldsymbol{x}) - f(\boldsymbol{x})
L=f(x+Δx)−f(x)