图优化理论框架

Graph Optimization

1、非线性最小二乘问题

SLAM问题通常可以表述为非线性最小二乘问题
F ( x ) = 1 2 ∑ ( x i , x j ) ∈ C e i j ( x i , x j , z i j ) T Ω i j e i j ( x i , x j , z i j ) x ∗ = arg ⁡ min ⁡ x   F ( x ) \begin{align*} F(\mathbf{x}) &=\frac{1}{2}\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in{C}}\mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})^{T}\boldsymbol\Omega_{ij}\mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})\tag{1}\\ \mathbf{x}^{\ast} &=\underset{\mathbf{x}}{\arg\min}\,F(\mathbf{x})\tag{2} \end{align*} F(x)x=21(xi,xj)Ceij(xi,xj,zij)TΩijeij(xi,xj,zij)=xargminF(x)(1)(2)
其中 x = [ x 1 T , x 2 T , ⋯   , x n T ] T ∈ R N \mathbf{x}=[\mathbf{x}_{1}^{T},\mathbf{x}_{2}^{T},\cdots,\mathbf{x}_{n}^{T}]^{T}\in\mathbb{R}^{N} x=[x1T,x2T,,xnT]TRN为全部参数组成的向量, x k , k = 1 , ⋯   , n \mathbf{x}_{k},k=1,\cdots,n xk,k=1,,n表示一个参数块。 C C C是全部参与求和的参数块组合, e i j ( x i , x j , z i j ) ∈ R M i j \mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})\in\mathbb{R}^{M_{ij}} eij(xi,xj,zij)RMij称为误差函数, Ω i j ∈ R M i j × M i j \boldsymbol\Omega_{ij}\in\mathbb{R}^{M_{ij}\times{M_{ij}}} ΩijRMij×Mij称为信息矩阵,为对称正定矩阵。注意,这里为方便推导,假定误差函数的自变量是两个参数块,但在实际情况下,误差函数自变量可以仅包含一个参数块,也可以包含更多的参数块,此时该问题的求解方式可以很方便地通过拓展下面的推导得到。

非线性最小二乘问题可以通过图来构建,其中顶点表示参数快,边表示误差函数,此时称之为图优化问题

在这里插入图片描述

2、求解非线性最小二乘问题

略写误差函数中的观测量,并将误差函数视作全部参数的函数,即采用下面的记法
e i j ( x i , x j , z i j ) = e i j ( x i , x j ) = e i j ( x ) \mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})=\mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j})=\mathbf{e}_{ij}(\mathbf{x}) eij(xi,xj,zij)=eij(xi,xj)=eij(x)
通过泰勒展开得到误差函数的一阶近似
e i j ( x i + Δ x i , x j + Δ x j ) = e i j ( x + Δ x ) ≈ e i j + J i j Δ x (3) \mathbf{e}_{ij}(\mathbf{x}_{i}+\Delta\mathbf{x}_{i},\mathbf{x}_{j}+\Delta\mathbf{x}_{j}) =\mathbf{e}_{ij}(\mathbf{x}+\Delta\mathbf{x}) \approx\mathbf{e}_{ij}+\mathbf{J}_{ij}\Delta\mathbf{x}\tag{3} eij(xi+Δxi,xj+Δxj)=eij(x+Δx)eij+JijΔx(3)
其中 e i j = e i j ( x ) \mathbf{e}_{ij}=\mathbf{e}_{ij}(\mathbf{x}) eij=eij(x) J i j = lim ⁡ Δ x → 0 e i j ( x + Δ x ) − e i j ( x ) Δ x ∈ R M i j × N \mathbf{J}_{ij}=\lim_{\Delta{x}\rightarrow\mathbf{0}}\frac{\mathbf{e}_{ij}(\mathbf{x}+\Delta\mathbf{x})-\mathbf{e}_{ij}(\mathbf{x})}{\Delta\mathbf{x}}\in\mathbb{R}^{M_{ij}\times{N}} Jij=limΔx0Δxeij(x+Δx)eij(x)RMij×N

考虑到参数 x \mathbf{x} x可能位于非欧式空间 D o m ( x ) \mathrm{Dom}(\mathbf{x}) Dom(x),而摄动量 Δ x \Delta\mathbf{x} Δx位于欧式空间 R N \mathbb{R}^{N} RN,故 ( 7 ) (7) (7)中的向量加法很可能导致 x + Δ x ∉ D o m ( x ) \mathbf{x}+\Delta\mathbf{x}\notin\mathrm{Dom}(\mathbf{x}) x+Δx/Dom(x)。为解决该问题,可采用广义加法 ⊕ : D o m ( x ) × R N → D o m ( x ) \oplus:\mathrm{Dom}(\mathbf{x})\times\mathbb{R}^{N}\rightarrow\mathrm{Dom}(\mathbf{x}) :Dom(x)×RNDom(x) ( 7 ) (7) (7)改写为
e i j ( x i ⊕ Δ x i , x j ⊕ Δ x j ) = e i j ( x ⊕ Δ x ) ≈ e i j + J i j Δ x \mathbf{e}_{ij}(\mathbf{x}_{i}\oplus\Delta\mathbf{x}_{i},\mathbf{x}_{j}\oplus\Delta\mathbf{x}_{j}) =\mathbf{e}_{ij}(\mathbf{x}\oplus\Delta\mathbf{x}) \approx\mathbf{e}_{ij}+\mathbf{J}_{ij}\Delta\mathbf{x} eij(xiΔxi,xjΔxj)=eij(xΔx)eij+JijΔx
其中 J i j = lim ⁡ Δ x → 0 e i j ( x ⊕ Δ x ) − e i j ( x ) Δ x ∈ R M i j × N \mathbf{J}_{ij}=\lim_{\Delta{x}\rightarrow\mathbf{0}}\frac{\mathbf{e}_{ij}(\mathbf{x}\oplus\Delta\mathbf{x})-\mathbf{e}_{ij}(\mathbf{x})}{\Delta\mathbf{x}}\in\mathbb{R}^{M_{ij}\times{N}} Jij=limΔx0Δxeij(xΔx)eij(x)RMij×N

F i j ( x i , x j , z i j ) = 1 2 e i j ( x i , x j , z i j ) T Ω i j e i j ( x i , x j , z i j ) F_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})=\frac{1}{2}\mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})^{T}\boldsymbol\Omega_{ij}\mathbf{e}_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij}) Fij(xi,xj,zij)=21eij(xi,xj,zij)TΩijeij(xi,xj,zij),同样略去误差观测量,并将 F i j F_{ij} Fij视作全部参数的函数,即采用下面的记法
F i j ( x i , x j , z i j ) = F i j ( x i , x j ) = F i j ( x ) F_{ij}(\mathbf{x}_{i},\mathbf{x}_{j},\mathbf{z}_{ij})=F_{ij}(\mathbf{x}_{i},\mathbf{x}_{j})=F_{ij}(\mathbf{x}) Fij(xi,xj,zij)=Fij(xi,xj)=Fij(x)
计算 F i j ( x + Δ x ) F_{ij}(\mathbf{x}+\Delta\mathbf{x}) Fij(x+Δx)
F i j ( x ⊕ Δ x ) = 1 2 e i j ( x ⊕ Δ x ) T Ω i j e i j ( x ⊕ Δ x ) = 1 2 ( e i j + J i j Δ x ) T Ω i j ( e i j + J i j Δ x ) = 1 2 ( e i j T Ω i j e i j ⏟ 2 c i j + 2 e i j T Ω i j J i j ⏟ b i j T Δ x + Δ x T J i j T Ω i j J i j ⏟ H i j Δ x ) = c i j + b i j T Δ x + 1 2 Δ x T H i j Δ x \begin{align*} F_{ij}(\mathbf{x}\oplus\Delta\mathbf{x}) &=\frac{1}{2}\mathbf{e}_{ij}(\mathbf{x}\oplus\Delta\mathbf{x})^{T}\boldsymbol\Omega_{ij}\mathbf{e}_{ij}(\mathbf{x}\oplus\Delta\mathbf{x})\\ &=\frac{1}{2}(\mathbf{e}_{ij}+\mathbf{J}_{ij}\Delta\mathbf{x})^{T}\boldsymbol\Omega_{ij}(\mathbf{e}_{ij}+\mathbf{J}_{ij}\Delta\mathbf{x})\\ &=\frac{1}{2}\left(\underbrace{\mathbf{e}_{ij}^{T}\boldsymbol\Omega_{ij}\mathbf{e}_{ij}}_{2c_{ij}}+2\underbrace{\mathbf{e}_{ij}^{T}\boldsymbol\Omega_{ij}\mathbf{J}_{ij}}_{\mathbf{b}_{ij}^{T}}\Delta\mathbf{x}+\Delta\mathbf{x}^{T}\underbrace{\mathbf{J}_{ij}^{T}\boldsymbol\Omega_{ij}\mathbf{J}_{ij}}_{\mathbf{H}_{ij}}\Delta\mathbf{x}\right)\\ &=c_{ij}+\mathbf{b}_{ij}^{T}\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}_{ij}\Delta\mathbf{x}\tag{4} \end{align*} Fij(xΔx)=21eij(xΔx)TΩijeij(xΔx)=21(eij+JijΔx)TΩij(eij+JijΔx)=21 2cij eijTΩijeij+2bijT eijTΩijJijΔx+ΔxTHij JijTΩijJijΔx =cij+bijTΔx+21ΔxTHijΔx(4)
( 4 ) (4) (4)代入 ( 1 ) (1) (1)计算 F ( x + Δ x ) F(\mathbf{x}+\Delta\mathbf{x}) F(x+Δx)
F ( x ⊕ Δ x ) = ∑ ( x i , x j ) ∈ C F i j ( x ⊕ Δ x ) = ∑ ( x i , x j ) ∈ C c i j + b i j Δ x + 1 2 Δ x T H i j Δ x = c + b T Δ x + 1 2 Δ x T H Δ x \begin{align*} F(\mathbf{x}\oplus\Delta\mathbf{x}) &=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}F_{ij}(\mathbf{x}\oplus\Delta\mathbf{x})\\ &=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}c_{ij}+\mathbf{b}_{ij}\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}_{ij}\Delta\mathbf{x}\\ &=c+\mathbf{b}^{T}\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}\Delta\mathbf{x}\tag{5} \end{align*} F(xΔx)=(xi,xj)CFij(xΔx)=(xi,xj)Ccij+bijΔx+21ΔxTHijΔx=c+bTΔx+21ΔxTHΔx(5)

其中 b = ∑ ( x i , x j ) ∈ C b i j ∈ R n \mathbf{b}=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}\mathbf{b}_{ij}\in\mathbb{R}^{n} b=(xi,xj)CbijRn H = ∑ ( x i , x j ) ∈ C H i j ∈ R N × N \mathbf{H}=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}\mathbf{H}_{ij}\in\mathbb{R}^{N\times{N}} H=(xi,xj)CHijRN×N

实际上, b \mathbf{b} b等于 F F F x \mathbf{x} x点处的梯度,而采用 H \mathbf{H} H近似 F F F x \mathbf{x} x点处的海森矩阵
b = ∑ ( x i , x j ) ∈ C [ ∇ e i j ( x ) Ω i j e i j ( x ) ] = ∇ F ( x ) H = ∑ ( x i , x j ) ∈ C [ ∇ e i j ( x ) Ω i j ∇ e i j ( x ) T ] → ∇ 2 F ( x ) = ∑ ( x i , x j ) ∈ C [ ⋯ + ∇ e i j ( x ) Ω i j ∇ e i j ( x ) T ] \begin{align*} \mathbf{b}&=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}\left[\nabla\mathbf{e}_{ij}(\mathbf{x})\boldsymbol{\Omega}_{ij}\mathbf{e}_{ij}(\mathbf{x})\right]=\nabla{F}(\mathbf{x})\\ \mathbf{H}&=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}\left[\nabla\mathbf{e}_{ij}(\mathbf{x})\boldsymbol{\Omega}_{ij}\nabla\mathbf{e}_{ij}(\mathbf{x})^{T}\right]\rightarrow\nabla^{2}F(\mathbf{x})=\sum_{(\mathbf{x}_{i},\mathbf{x}_{j})\in\mathcal{C}}\left[\cdots+\nabla\mathbf{e}_{ij}(\mathbf{x})\boldsymbol{\Omega}_{ij}\nabla\mathbf{e}_{ij}(\mathbf{x})^{T}\right] \end{align*} bH=(xi,xj)C[eij(x)Ωijeij(x)]=F(x)=(xi,xj)C[eij(x)Ωijeij(x)T]2F(x)=(xi,xj)C[+eij(x)Ωijeij(x)T]
2.1、梯度下降法

梯度下降法通过取增量方向为目标函数的负梯度方向来保证目标函数的一阶近似下降

保留 ( 5 ) (5) (5)至一阶项
F ( x ⊕ Δ x ) = c + b T Δ x F(\mathbf{x}\oplus\Delta\mathbf{x})=c+\mathbf{b}^{T}\Delta\mathbf{x} F(xΔx)=c+bTΔx
此时取
Δ x = − λ b (6) \Delta\mathbf{x}=-\lambda\mathbf{b}\tag{6} Δx=λb(6)
即可保证目标函数的一阶近似下降,其中 λ > 0 \lambda>0 λ>0称为步长

梯度下降法求解非线性最小二乘问题的具体步骤如下:


  1. k = 0 k=0 k=0,给定初始值 x 0 \mathbf{x}_{0} x0 λ \lambda λ
  2. k k k达到最大迭代次数,则停止迭代;否则根据 x k \mathbf{x}_{k} xk求出当前的和 b k \mathbf{b}_{k} bk
  3. Δ x k = − λ b k \Delta\mathbf{x}_{k}=-\lambda\mathbf{b}_{k} Δxk=λbk
  4. 如果 Δ x k \Delta\mathbf{x}_{k} Δxk足够小,则停止迭代;否则令 x k + 1 = x k ⊕ Δ x k \mathbf{x}_{k+1}=\mathbf{x}_{k}\oplus\Delta\mathbf{x}_{k} xk+1=xkΔxk k = k + 1 k=k+1 k=k+1,返回第2步

2.2、高斯牛顿法

高斯牛顿法采用 H \mathbf{H} H近似 F ( x ) F(\mathbf{x}) F(x) x \mathbf{x} x点处的海森矩阵,避免了海森矩阵的计算,又可以通过求解二阶近似的最小值实现更快的收敛

因为 H \mathbf{H} H矩阵半正定,故根据 ( 8 ) (8) (8) Δ x \Delta\mathbf{x} Δx的函数 F ( x ⊕ Δ x ) F(\mathbf{x}\oplus\Delta\mathbf{x}) F(xΔx)为凸函数,因此其最小值在 ∂ F ( x ⊕ Δ x ) ∂ Δ x = 0 \frac{\partial{F(\mathbf{x}\oplus\Delta\mathbf{x})}}{\partial\Delta\mathbf{x}}=\mathbf{0} ΔxF(xΔx)=0处取得,即要求 Δ x \Delta\mathbf{x} Δx满足下式时取得
H Δ x = − b (7) \mathbf{H}\Delta\mathbf{x}=-\mathbf{b}\tag{7} HΔx=b(7)
( 11 ) (11) (11)称为增量方程,求解增量方程是整个优化问题的核心。

高斯牛顿法求解非线性最小二乘问题的具体步骤如下:


  1. k = 0 k=0 k=0,给定初始值 x 0 \mathbf{x}_{0} x0
  2. k k k达到最大迭代次数,则停止迭代;否则根据 x k \mathbf{x}_{k} xk求出当前的 H k \mathbf{H}_{k} Hk b k \mathbf{b}_{k} bk
  3. 求解增量方程: H k Δ x k = − b k \mathbf{H}_{k}\Delta\mathbf{x}_{k}=-\mathbf{b}_{k} HkΔxk=bk:当 H \mathbf{H} H正定时,增量方程可以通过Cholesky分解高效求解;当 H \mathbf{H} H非正定时,可取 Δ x k = Δ x k − 1 \Delta\mathbf{x}_{k}=\Delta\mathbf{x}_{k-1} Δxk=Δxk1
  4. 如果 Δ x k \Delta\mathbf{x}_{k} Δxk足够小,则停止迭代;否则令 x k + 1 = x k ⊕ Δ x k \mathbf{x}_{k+1}=\mathbf{x}_{k}\oplus\Delta\mathbf{x}_{k} xk+1=xkΔxk k = k + 1 k=k+1 k=k+1,返回第2步

2.3、Levenberg-Marquardt算法

Levenberg-Marquardt算法,简称LM算法。考虑到高斯牛顿法中对海森矩阵的近似不一定准确,二阶泰勒展开式也只能才展开点附近有较好的效果,为结合梯度下降法和高斯牛顿法的优点,可以将增量方程改写为
( H + λ I ) Δ x = − b (8) (\mathbf{H}+\lambda\mathbf{I})\Delta\mathbf{x}=-\mathbf{b}\tag{8} (H+λI)Δx=b(8)
此方程称为LM算法的增量方程

此外,通过 ρ = F ( x ⊕ Δ x ) − F ( x ) ∇ F ( x ) T Δ x \rho=\frac{F(\mathbf{x}\oplus\Delta\mathbf{x})-F(\mathbf{x})}{\nabla{F(\mathbf{x})^{T}\Delta\mathbf{x}}} ρ=F(x)TΔxF(xΔx)F(x)反映一阶近似的准确度,如果一阶近似准确度高,则 λ \lambda λ值增大,更接近梯度下降法,如果一阶近似准确度低,即函数非线性的性质更加突出,则 λ \lambda λ值减小,更接近高斯牛顿法。

LM算法求解非线性最小二乘问题的具体步骤如下:


  1. k = 0 k=0 k=0,给定初始值 x 0 \mathbf{x}_{0} x0 λ 0 \lambda_{0} λ0
  2. k k k达到最大迭代次数,则停止迭代;否则根据 x k \mathbf{x}_{k} xk求出当前的 H k , b k , ρ k \mathbf{H}_{k},\mathbf{b}_{k},\rho_{k} Hk,bk,ρk
  3. ρ k > 3 4 \rho_{k}>\frac{3}{4} ρk>43,则 λ k + 1 = 2 λ k \lambda_{k+1}=2\lambda_{k} λk+1=2λk;若 ρ k < 1 4 \rho_{k}<\frac{1}{4} ρk<41,则 λ k + 1 = 1 2 λ k \lambda_{k+1}=\frac{1}{2}\lambda_{k} λk+1=21λk;若 1 4 ⩽ ρ k ⩽ 3 4 \frac{1}{4}\leqslant\rho_{k}\leqslant\frac{3}{4} 41ρk43,则 λ k + 1 = λ k \lambda_{k+1}=\lambda_{k} λk+1=λk
  4. 求解 ( H k + λ k + 1 I ) Δ x k = − b k (\mathbf{H}_{k}+\lambda_{k+1}\mathbf{I})\Delta\mathbf{x}_{k}=-\mathbf{b}_{k} (Hk+λk+1I)Δxk=bk:因为 ( H k + λ k + 1 I ) (\mathbf{H}_{k}+\lambda_{k+1}\mathbf{I}) (Hk+λk+1I)恒为正定,故可以通过Cholesky分解高效求解
  5. 如果 Δ x k \Delta\mathbf{x}_{k} Δxk足够小,则停止迭代; x k + 1 = x k ⊕ Δ x k \mathbf{x}_{k+1}=\mathbf{x}_{k}\oplus\Delta\mathbf{x}_{k} xk+1=xkΔxk k = k + 1 k=k+1 k=k+1,返回第2步

3、核函数

注意到,当某条边的误差 e i j \mathbf{e}_{ij} eij很大时, F i j ( x ) = e i j T Ω i j e i j F_{ij}(\mathbf{x})=\mathbf{e}_{ij}^{T}\boldsymbol{\Omega}_{ij}\mathbf{e}_{ij} Fij(x)=eijTΩijeij会很大,其梯度 ∇ F i j ( x ) = J i j T Ω i j e i j \nabla{F}_{ij}(\mathbf{x})=\mathbf{J}_{ij}^{T}\boldsymbol{\Omega}_{ij}\mathbf{e}_{ij} Fij(x)=JijTΩijeij也会很大,而算法会根据梯度更大幅度地调整这条边所连接的节点的估计值,而掩盖这些节点与其他边的关系。

可以定义核函数 ρ ( ⋅ ) \rho(\cdot) ρ(),在非线性最小二乘问题中用 ρ ( F i j ( x ) ) \rho(F_{ij}(\mathbf{x})) ρ(Fij(x))代替 F i j ( x ) F_{ij}(\mathbf{x}) Fij(x)

最常用的核函数是Huber核函数,其定义为
ρ ( x ) = { x 0 ⩽ x ⩽ δ 2 2 δ ( x − δ 2 ) x > δ 2 \rho(x)=\left\{ \begin{array}{ll} x&0\leqslant{x}\leqslant\delta^{2}\\ 2\delta(\sqrt{x}-\frac{\delta}{2})&{x}>\delta^{2} \end{array} \right. ρ(x)={x2δ(x 2δ)0xδ2x>δ2
Huber函数图像如下

Huber

Huber核函数是二阶可导的,其一阶导数和二阶导数为
ρ ′ ( x ) = { 1 0 ⩽ x ⩽ δ 2 δ x x > δ 2 ρ ′ ′ ( x ) = { 0 0 ⩽ x ⩽ δ 2 − 1 2 δ ( x ) 3 x > δ 2 \begin{align*} \rho^{\prime}(x)=\left\{ \begin{array}{ll} 1&0\leqslant{x}\leqslant\delta^{2}\\ \frac{\delta}{\sqrt{x}}&{x}>\delta^{2} \end{array} \right.\quad\quad\quad \rho^{\prime\prime}(x)=\left\{ \begin{array}{ll} 0&0\leqslant{x}\leqslant\delta^{2}\\ -\frac{1}{2}\frac{\delta}{(\sqrt{x})^{3}}&{x}>\delta^{2} \end{array} \right.\\ \end{align*} ρ(x)={1x δ0xδ2x>δ2ρ′′(x)={021(x )3δ0xδ2x>δ2
ρ ( F i j ( x ⊕ Δ x ) ) \rho(F_{ij}(\mathbf{x}\oplus\Delta\mathbf{x})) ρ(Fij(xΔx))做如下展开
ρ ( F i j ( x ⊕ Δ x ) ) = ρ ( F i j ( x ) ) + ∇ ρ ( F i j ( x ) ) T Δ x + 1 2 Δ x T ∇ 2 ρ ( F i j ( x ) ) Δ x \rho(F_{ij}(\mathbf{x}\oplus\Delta\mathbf{x}))=\rho(F_{ij}(\mathbf{x}))+\nabla\rho(F_{ij}(\mathbf{x}))^{T}\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\nabla^{2}\rho(F_{ij}(\mathbf{x}))\Delta\mathbf{x} ρ(Fij(xΔx))=ρ(Fij(x))+ρ(Fij(x))TΔx+21ΔxT2ρ(Fij(x))Δx
其中
∇ ρ ( F i j ( x ) ) = ∂ ρ ∂ F i j ∣ F i j ( x ) ∇ F i j ( x ) = ρ ′ ( F i j ( x ) ) b i j = J i j T ( ρ ′ ( F i j ( x ) ) Ω i j ) e i j ∇ 2 ρ ( F i j ( x ) ) = ∂ 2 ρ ∂ F i j 2 ∣ F i j ( x ) ∇ F i j ( x ) ∇ F i j ( x ) T + ∂ ρ ∂ F i j ∣ F i j ( x ) ∇ 2 F i j ( x ) = ρ ′ ′ ( F i j ( x ) ) b i j b i j T + ρ ′ ( F i j ( x ) ) H i j → ρ ′ ( F i j ( x ) ) H i j = J i j T ( ρ ′ ( F i j ( x ) ) Ω i j ) J i j \begin{align*} \nabla\rho(F_{ij}(\mathbf{x}))&=\left.\frac{\partial\rho}{\partial{F}_{ij}}\right|_{F_{ij}(\mathbf{x})}\nabla{F}_{ij}(\mathbf{x})=\rho^{\prime}(F_{ij}(\mathbf{x}))\mathbf{b}_{ij}=\mathbf{J}_{ij}^{T}\left(\rho^{\prime}(F_{ij}(\mathbf{x}))\boldsymbol{\Omega}_{ij}\right)\mathbf{e}_{ij}\\ \nabla^{2}\rho(F_{ij}(\mathbf{x}))&=\left.\frac{\partial^{2}\rho}{\partial{F}_{ij}^{2}}\right|_{F_{ij}(\mathbf{x})}\nabla{F}_{ij}(\mathbf{x})\nabla{F}_{ij}(\mathbf{x})^{T}+\left.\frac{\partial\rho}{\partial{F}_{ij}}\right|_{F_{ij}(\mathbf{x})}\nabla^{2}{F}_{ij}(\mathbf{x})=\rho^{\prime\prime}(F_{ij}(\mathbf{x}))\mathbf{b}_{ij}\mathbf{b}_{ij}^{T}+\rho^{\prime}(F_{ij}(\mathbf{x}))\mathbf{H}_{ij}\rightarrow\rho^{\prime}(F_{ij}(\mathbf{x}))\mathbf{H}_{ij}=\mathbf{J}_{ij}^{T}\left(\rho^{\prime}(F_{ij}(\mathbf{x}))\boldsymbol{\Omega}_{ij}\right)\mathbf{J}_{ij} \end{align*} ρ(Fij(x))2ρ(Fij(x))=Fijρ Fij(x)Fij(x)=ρ(Fij(x))bij=JijT(ρ(Fij(x))Ωij)eij=Fij22ρ Fij(x)Fij(x)Fij(x)T+Fijρ Fij(x)2Fij(x)=ρ′′(Fij(x))bijbijT+ρ(Fij(x))Hijρ(Fij(x))Hij=JijT(ρ(Fij(x))Ωij)Jij
上式第二行中,因为考虑到可能存在 ρ ′ ′ ( F i j ( x ) ) < 0 \rho^{\prime\prime}(F_{ij}(\mathbf{x}))<0 ρ′′(Fij(x))<0,故为确保 ∇ 2 ρ ( F i j ( x ) ) \nabla^{2}\rho(F_{ij}(\mathbf{x})) 2ρ(Fij(x))正定性与 H i j \mathbf{H}_{ij} Hij一致,舍去含 ρ ′ ′ ( F i j ( x ) ) \rho^{\prime\prime}(F_{ij}(\mathbf{x})) ρ′′(Fij(x))的项。由上面的分析知,对某条边应用核函数后,在优化过程中只需修改其对应的信息矩阵即可。

b i j ← ∇ ρ ( F i j ( x ) ) \mathbf{b}_{ij}\leftarrow\nabla\rho(F_{ij}(\mathbf{x})) bijρ(Fij(x)) H i j ← ∇ 2 ρ ( F i j ( x ) ) \mathbf{H}_{ij}\leftarrow\nabla^{2}\rho(F_{ij}(\mathbf{x})) Hij2ρ(Fij(x)),从而可以构建新的增量方程。

通过上述过程,可以解决因误差 e i j \mathbf{e}_{ij} eij过大而导致的问题。

附录

§1、标量函数的二阶泰勒展开

若一个多元标量函数 f : R N → R f:\mathbb{R}^{N}\rightarrow\mathbb{R} f:RNR a \mathbf{a} a点处二阶可微,则存在 a \mathbf{a} a的邻域 U ( x , δ ) = { x + Δ x ∣ ∥ Δ x ∥ < δ } U(\mathbf{x},\delta)=\{\mathbf{x}+\Delta\mathbf{x}|\|\Delta\mathbf{x}\|<\delta\} U(x,δ)={x+Δx∣∥Δx<δ},使得在该邻域中
f ( x + Δ x ) = f ( x ) + J f ( x ) Δ x + 1 2 Δ x T H f ( x ) Δ x + o ( ∥ Δ x ∥ 2 ) (A1) f(\mathbf{x}+\Delta\mathbf{x})=f(\mathbf{x})+\mathbf{J}f(\mathbf{x})\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}f(\mathbf{x})\Delta\mathbf{x}+\mathrm{o}(\|\Delta\mathbf{x}\|^{2})\tag{A1} f(x+Δx)=f(x)+Jf(x)Δx+21ΔxTHf(x)Δx+o(∥Δx2)(A1)
其中

J f ( x ) = ∇ f ( x ) T = lim ⁡ Δ x → 0 f ( x + Δ x ) − f ( x ) Δ x T ∈ R 1 × N \mathbf{J}f(\mathbf{x})=\nabla{f}(\mathbf{x})^{T}=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{f(\mathbf{x}+\Delta\mathbf{x})-f(\mathbf{x})}{\Delta\mathbf{x}^{T}}\in\mathbb{R}^{1\times N} Jf(x)=f(x)T=limΔx0ΔxTf(x+Δx)f(x)R1×N称为雅可比矩阵;

H f ( x ) = ∇ 2 f ( x ) = lim ⁡ Δ x → 0 ∇ f ( x + Δ x ) − ∇ f ( x ) Δ x T ∈ R N × N \mathbf{H}f(\mathbf{x})=\nabla^{2}{f}(\mathbf{x})=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{\nabla{f}(\mathbf{x}+\Delta\mathbf{x})-\nabla{f}(\mathbf{x})}{\Delta\mathbf{x}^{T}}\in\mathbb{R}^{N\times{N}} Hf(x)=2f(x)=limΔx0ΔxTf(x+Δx)f(x)RN×N称为海森矩阵。

( A 1 ) \mathrm{(A1)} (A1)称为标量函数 f f f的二阶泰勒展开。

下面根据上述理论推导当 f \mathbf{f} f的定义域为非欧式空间 D o m ( x ) \mathrm{Dom}(\mathbf{x}) Dom(x)时的二阶泰勒展开

h ( Δ x ) = f ( x ⊕ Δ x ) : R N → R h(\Delta\mathbf{x})=f(\mathbf{x}\oplus\Delta\mathbf{x}):\mathbb{R}^{N}\rightarrow\mathbb{R} h(Δx)=f(xΔx):RNR ⊕ : D o m ( x ) × R N → D o m ( x ) \oplus:\mathrm{Dom}(\mathbf{x})\times\mathbb{R}^{N}\rightarrow\mathrm{Dom}(\mathbf{x}) :Dom(x)×RNDom(x),若 h h h 0 \mathbf{0} 0点处二阶可微,则存在 0 \mathbf{0} 0的邻域 U ( 0 , δ ) = { Δ x ∣ ∥ Δ x ∥ < δ } U(\mathbf{0},\delta)=\{\Delta\mathbf{x}|\|\Delta\mathbf{x}\|<\delta\} U(0,δ)={Δx∣∥Δx<δ},使得在该邻域中
h ( Δ x ) = h ( 0 ) + J h ( 0 ) T Δ x + 1 2 Δ x T H h ( 0 ) Δ x + o ( ∥ Δ x ∥ 2 ) (A2) h(\Delta\mathbf{x})=h(\mathbf{0})+\mathbf{J}h(\mathbf{0})^{T}\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}h(\mathbf{0})\Delta\mathbf{x}+\mathrm{o}(\|\Delta\mathbf{x}\|^{2})\tag{A2} h(Δx)=h(0)+Jh(0)TΔx+21ΔxTHh(0)Δx+o(∥Δx2)(A2)
其中

J h ( 0 ) = ∇ h ( 0 ) T = lim ⁡ Δ x → 0 h ( Δ x ) − h ( 0 ) Δ x T = lim ⁡ Δ x → 0 f ( x ⊕ Δ x ) − f ( x ) Δ x T ∈ R 1 × N \mathbf{J}h(\mathbf{0})=\nabla{h}(\mathbf{0})^{T}=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{h(\Delta\mathbf{x})-h(\mathbf{0})}{\Delta\mathbf{x}^{T}}=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{f(\mathbf{x}\oplus\Delta\mathbf{x})-f(\mathbf{x})}{\Delta\mathbf{x}^{T}}\in\mathbb{R}^{1\times N} Jh(0)=h(0)T=limΔx0ΔxTh(Δx)h(0)=limΔx0ΔxTf(xΔx)f(x)R1×N,可以发现 J h ( 0 ) \mathbf{J}h(\mathbf{0}) Jh(0) ( A 1 ) \mathrm{(A1)} (A1) J f ( x ) \mathbf{J}f(\mathbf{x}) Jf(x)具有类似的结构,故记 J f ( x ) = ∇ f ( x ) T = J h ( 0 ) \mathbf{J}f(\mathbf{x})=\nabla{f}(\mathbf{x})^{T}=\mathbf{J}h(\mathbf{0}) Jf(x)=f(x)T=Jh(0),称 J f ( x ) \mathbf{J}f(\mathbf{x}) Jf(x)为对增量 Δ x \Delta\mathbf{x} Δx的雅可比矩阵;

H h ( 0 ) = ∇ 2 h ( 0 ) = lim ⁡ Δ x → 0 ∇ h ( Δ x ) − ∇ h ( 0 ) Δ x T = lim ⁡ Δ x → 0 ∇ f ( x ⊕ Δ x ) − ∇ f ( x ) Δ x T ∈ R N × N \mathbf{H}h(\mathbf{0})=\nabla^{2}h(\mathbf{0})=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{\nabla{h}(\Delta\mathbf{x})-\nabla{h}(\mathbf{0})}{\Delta\mathbf{x}^{T}}=\lim_{\Delta\mathbf{x}\rightarrow\mathbf{0}}\frac{\nabla{f}(\mathbf{x}\oplus\Delta\mathbf{x})-\nabla{f}(\mathbf{x})}{\Delta\mathbf{x}^{T}}\in\mathbb{R}^{N\times{N}} Hh(0)=2h(0)=limΔx0ΔxTh(Δx)h(0)=limΔx0ΔxTf(xΔx)f(x)RN×N,可以发现 H h ( 0 ) \mathbf{H}h(\mathbf{0}) Hh(0) ( A 1 ) \mathrm{(A1)} (A1) H f ( x ) \mathbf{H}f(\mathbf{x}) Hf(x)具有类似的结构,故记 H f ( x ) = ∇ 2 f ( x ) = H h ( 0 ) \mathbf{H}f(\mathbf{x})=\nabla^{2}f(\mathbf{x})=\mathbf{H}h(\mathbf{0}) Hf(x)=2f(x)=Hh(0),称 H f ( x ) \mathbf{H}f(\mathbf{x}) Hf(x)为对增量 Δ x \Delta\mathbf{x} Δx的海森矩阵。

综上, ( A 2 ) \mathrm{(A2)} (A2)也可以写成
f ( x ⊕ Δ x ) = f ( x ) + J f ( x ) Δ x + 1 2 Δ x T H f ( x ) Δ x + o ( ∥ Δ x ∥ 2 ) (A3) f(\mathbf{x}\oplus\Delta\mathbf{x})=f(\mathbf{x})+\mathbf{J}f(\mathbf{x})\Delta\mathbf{x}+\frac{1}{2}\Delta\mathbf{x}^{T}\mathbf{H}f(\mathbf{x})\Delta\mathbf{x}+\mathrm{o}(\|\Delta\mathbf{x}\|^{2})\tag{A3} f(xΔx)=f(x)+Jf(x)Δx+21ΔxTHf(x)Δx+o(∥Δx2)(A3)
( A 3 ) \mathrm{(A3)} (A3)称为 D o m ( x ) \mathrm{Dom}(\mathbf{x}) Dom(x)上的标量函数 f \mathbf{f} f的二阶泰勒展开。

为简便起见,常将上面的雅可比矩阵和海森矩阵分别用 J ( ⋅ ) \mathbf{J}(\cdot) J() H ( ⋅ ) \mathbf{H}(\cdot) H()代替

  • 5
    点赞
  • 8
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

小于小于大橙子

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值