格物致知之卡尔曼滤波(Kalman Filtering)——(7)扩展卡尔曼滤波(EKF)一阶滤波
卡尔曼滤波器是为线性系统设计的算法,不能用于非线性系统中。但是事实上在实际工程中多数系统都是非线性的,所以如果卡尔曼滤波器不能用在非线性系统中的话,那么它的应用范围就非常有限了。扩展卡尔曼滤波(EKF)就是为了打破卡尔曼滤波器只能用于线性系统的局限性而设计出来的,一般解决非线性系统线性化是通过泰勒级数展开,省略高阶项得到一个近似的线性系统。这种处理非线性系统的方法被称为扩展卡尔曼滤波(EKF)。
一、状态空间模型
二、泰勒级数展开公式
泰勒级数展开是将一个在
x
=
x
0
x=x_{0}
x=x0 处具有
n
n
n 阶导数的函数
f
(
x
)
f(x)
f(x), 利用关于
(
x
−
x
0
)
\left(x-x_{0}\right)
(x−x0) 的
n
n
n 次多项式逼近函数 值的方法。
若函数
f
(
x
)
f(x)
f(x) 在包含
x
0
x_{0}
x0 的某个闭区间
[
a
,
b
]
[a, b]
[a,b] 上具有
n
n
n 阶导数,且在开区间
(
a
,
b
)
(a, b)
(a,b) 上具有
(
n
+
1
)
(n+1)
(n+1) 阶导数,则 对闭区间
[
a
,
b
]
[a, b]
[a,b] 上的任意一点
x
x
x ,都有:
f
(
x
)
=
f
(
x
0
)
0
!
+
f
′
(
x
0
)
1
!
(
x
−
x
0
)
+
…
+
f
(
n
)
(
x
0
)
n
!
(
x
−
x
0
)
n
+
R
n
(
x
)
f(x)=\frac{f\left(x_{0}\right)}{0 !}+\frac{f^{\prime}\left(x_{0}\right)}{1 !}\left(x-x_{0}\right)+\ldots+\frac{f^{(n)}\left(x_{0}\right)}{n !}\left(x-x_{0}\right)^{n}+R_{n}(x)
f(x)=0!f(x0)+1!f′(x0)(x−x0)+…+n!f(n)(x0)(x−x0)n+Rn(x)
其中
f
(
n
)
(
x
0
)
f^{(n)}\left(x_{0}\right)
f(n)(x0) 表示函数
f
(
x
)
f(x)
f(x) 在
x
=
x
0
x=x_{0}
x=x0 处的
n
n
n 阶导数,等式右边成为泰勒展开式,剩余的
R
n
(
x
)
R_{n}(x)
Rn(x) 是泰勒展 开式的余项,是
(
x
−
x
0
)
n
\left(x-x_{0}\right)^{n}
(x−x0)n 的高阶无穷小。
三、麦克劳林展开公式
函数的麦克劳林展开指上面泰勒公式中
x
0
x_{0}
x0 取0的情况,即是泰勒公式的特殊形式,若
f
(
x
)
f(x)
f(x) 在
x
=
0
x=0
x=0 处n阶连续可导,则下式成立:
f
(
x
)
=
f
(
0
)
+
f
′
(
0
)
1
!
x
+
f
′
′
(
0
)
2
!
x
2
+
f
′
′
′
(
0
)
3
!
x
3
+
⋯
+
f
(
n
)
(
0
)
n
!
x
n
+
R
n
(
x
)
f(x)=f(0)+\frac{f^{\prime}(0)}{1 !} x+\frac{f^{\prime \prime}(0)}{2 !} x^{2}+\frac{f^{\prime \prime \prime}(0)}{3 !} x^{3}+\cdots+\frac{f^{(n)}(0)}{n !} x^{n}+R_{n}(x)
f(x)=f(0)+1!f′(0)x+2!f′′(0)x2+3!f′′′(0)x3+⋯+n!f(n)(0)xn+Rn(x)
其中
f
(
n
)
(
x
)
f^{(n)}(x)
f(n)(x) 表示
f
(
x
)
f(x)
f(x) 的n阶导数。
麦克劳林展开的形象描述图:
实际应用中,泰勒公式需要截断,只取有限项,一个函数的有限项的泰勒级数叫做泰勒展开式。泰勒公式的余项可以用于估算这种近似的误差。 泰勒展开式的重要性体现在以下三个方面:
-
幂级数的求导和积分可以逐项进行,因此求和函数相对比较容易。
-
一个解析函数可被延伸为一个定义在复平面上的一个开片上的解析函数,并使得复分析这种手法可行。
-
泰勒级数可以用来近似计算函数的值。
三、多维向量函数泰勒级数一维展开
当变量是多维向量时,一维的泰勒展开就需要做拓展,具体形式如下:
f
(
x
)
=
f
(
x
k
)
+
[
∇
f
(
x
k
)
]
T
(
x
−
x
k
)
+
o
n
f(\mathbf{x})=f\left(\mathbf{x}_{k}\right)+\left[\nabla f\left(\mathbf{x}_{k}\right)\right]^{T}\left(\mathbf{x}-\mathbf{x}_{k}\right)+o^{n}
f(x)=f(xk)+[∇f(xk)]T(x−xk)+on
其中,
[
∇
f
(
x
k
)
]
T
=
J
(
x
k
)
\left[\nabla f\left(\mathbf{x}_{k}\right)\right]^{T}=\mathbf{J}\left(\mathbf{x}_{k}\right)
[∇f(xk)]T=J(xk) 表示雅克比矩阵,
o
n
\mathbf{o}^{n}
on 表示高阶无穷小。
[
∇
f
(
x
)
]
T
=
J
(
x
)
=
[
∂
f
(
x
)
1
∂
x
1
…
∂
f
(
x
)
1
∂
x
n
⋮
⋱
⋮
∂
f
(
x
)
m
∂
x
1
…
∂
f
(
x
)
m
∂
x
n
]
[\nabla f(\mathbf{x})]^{T}=\mathbf{J}(\mathbf{x})=\left[\begin{array}{ccc} \frac{\partial f(\mathbf{x})_{1}}{\partial \mathbf{x}_{1}} & \ldots & \frac{\partial f(\mathbf{x})_{1}}{\partial \mathbf{x}_{n}} \\ \vdots & \ddots & \vdots \\ \frac{\partial f(\mathbf{x})_{m}}{\partial \mathbf{x}_{1}} & \ldots & \frac{\partial f(\mathbf{x})_{m}}{\partial \mathbf{x}_{n}} \end{array}\right]
[∇f(x)]T=J(x)=⎣⎢⎢⎡∂x1∂f(x)1⋮∂x1∂f(x)m…⋱…∂xn∂f(x)1⋮∂xn∂f(x)m⎦⎥⎥⎤
这里,
f
(
x
k
)
f\left(\mathbf{x}_{k}\right)
f(xk) 为
m
m
m 维,
x
k
\mathbf{x}_{k}
xk 状态向量为
n
n
n 维,
∂
2
f
(
x
k
)
∂
x
m
∂
x
n
=
∂
f
(
x
k
)
T
∂
x
m
∂
f
(
x
k
)
∂
x
n
\frac{\partial^{2} f\left(\mathbf{x}_{k}\right)}{\partial x_{m} \partial x_{n}}=\frac{\partial f\left(\mathbf{x}_{k}\right)^{T}}{\partial x_{m}} \frac{\partial f\left(\mathbf{x}_{k}\right)}{\partial x_{n}}
∂xm∂xn∂2f(xk)=∂xm∂f(xk)T∂xn∂f(xk) 。
一般来说,EKF在对非线性函数做泰勒展开时,只取到一阶导和二阶导,而由于二阶导的计算复杂性,更多的实际应用只取到一阶导,同样也能有较好的结果。取一阶导时,状态转移方程和观测方程 就近似为线性方程,高斯分布的变量经过线性变换之后仍然是高斯分布,这样就能够延用标准卡尔曼滤波的框架。
五、 一阶泰勒展式取近似
非线性状态方程
{
X
k
=
f
(
x
k
−
1
,
u
k
−
1
,
w
k
−
1
)
Z
k
=
h
(
x
k
,
v
k
)
w
k
−
1
∼
N
(
0
,
Q
)
,
v
k
−
1
∼
N
(
0
,
P
)
\begin{aligned} \begin{cases} X_k=f(x_{k-1},u_{k-1},w_{k-1}) \\ Z_k=h(x_k,v_k) \end{cases} \\w_{k-1}\sim N(0,Q) ,v_{k-1}\sim N(0,P) \end{aligned}
{Xk=f(xk−1,uk−1,wk−1)Zk=h(xk,vk)wk−1∼N(0,Q),vk−1∼N(0,P)
虽然系统有误差但是其误差任然服从正态分布,由于系统有误差,也不能在真实值处进行线性化,所以
f
(
x
k
)
f({x}_{k})
f(xk)在
x
^
k
−
1
\hat{x}_{k-1}
x^k−1(
k
−
1
k-1
k−1时刻的后验估计即上一次)处进行线性化操作。 1
用泰勒展开进行线性化操作:
X
k
=
f
(
x
k
−
1
,
u
k
−
1
,
w
k
−
1
)
+
A
(
x
k
−
x
^
k
−
1
)
+
w
w
k
−
1
\begin{aligned} &X_k=f(x_{k-1},u_{k-1},w_{k-1})+A(x_k-\hat{x}_{k-1})+ww_{k-1}\\ \end{aligned}
Xk=f(xk−1,uk−1,wk−1)+A(xk−x^k−1)+wwk−1
假设误差
w
k
−
1
=
0
w_{k-1}=0
wk−1=0则:
X
k
=
f
(
x
k
−
1
,
u
k
−
1
,
0
)
+
A
(
x
k
−
x
^
k
−
1
)
+
w
w
k
−
1
A
=
∂
f
∂
x
∣
x
^
k
−
1
,
u
k
−
1
w
=
∂
f
∂
w
∣
x
^
k
−
1
,
u
k
−
1
\begin{aligned} \\&X_k=f(x_{k-1},u_{k-1},0)+A(x_k-\hat{x}_{k-1})+ww_{k-1} \\&A=\frac{\partial f}{\partial x}\mid_{ \hat{x}_{k-1},u_{k-1}} \\&w=\frac{\partial f}{\partial w}\mid_{ \hat{x}_{k-1},u_{k-1}} \end{aligned}
Xk=f(xk−1,uk−1,0)+A(xk−x^k−1)+wwk−1A=∂x∂f∣x^k−1,uk−1w=∂w∂f∣x^k−1,uk−1
令
f
(
x
k
−
1
,
u
k
−
1
,
0
)
=
x
~
k
f(x_{k-1},u_{k-1},0)=\tilde{x}_k
f(xk−1,uk−1,0)=x~k,则
Z
k
Z_k
Zk在
x
~
k
\tilde{x}_k
x~k处线性化:
Z
k
=
h
(
x
~
k
,
v
k
)
+
H
(
x
k
−
x
~
k
)
+
v
v
k
\begin{aligned} \\&Z_k=h(\tilde{x}_k,v_k)+H(x_k-\tilde{x}_k)+vv_k \end{aligned}
Zk=h(x~k,vk)+H(xk−x~k)+vvk
假设误差
v
k
=
0
v_{k}=0
vk=0则:
Z
k
=
h
(
x
~
k
,
0
)
+
H
(
x
k
−
x
~
k
)
+
v
v
k
H
=
∂
h
∂
x
∣
x
~
k
v
=
∂
h
∂
v
∣
x
~
k
\begin{aligned} \\&Z_k=h(\tilde{x}_k,0)+H(x_k-\tilde{x}_k)+vv_k \\&H=\frac{\partial h}{\partial x}\mid_{ \tilde{x}_{k}} \\&v=\frac{\partial h}{\partial v}\mid_{ \tilde{x}_{k}} \end{aligned}
Zk=h(x~k,0)+H(xk−x~k)+vvkH=∂x∂h∣x~kv=∂v∂h∣x~k
令
h
(
x
~
k
,
0
)
=
z
~
k
h(\tilde{x}_k,0)=\tilde{z}_k
h(x~k,0)=z~k,则:
{
X
k
=
x
~
k
+
A
(
x
k
−
x
^
k
−
1
)
+
w
w
k
−
1
Z
k
=
z
~
k
+
H
(
x
k
−
x
~
k
)
+
v
v
k
\begin{aligned} \\\begin{cases} X_k=\tilde{x}_k+A(x_k-\hat{x}_{k-1})+ww_{k-1} \\Z_k=\tilde{z}_k+H(x_k-\tilde{x}_k)+vv_k \end{cases} \end{aligned}
{Xk=x~k+A(xk−x^k−1)+wwk−1Zk=z~k+H(xk−x~k)+vvk
由于
P
(
w
)
∼
N
(
0
,
Q
)
P(w)\sim N(0,Q)
P(w)∼N(0,Q),则:
{
P
(
w
w
k
−
1
)
∼
N
(
0
,
W
Q
W
T
)
P
(
v
v
k
)
∼
N
(
0
,
V
P
V
T
)
\begin{aligned} \begin{cases} P(ww_{k-1})\sim N(0,WQW^T) \\P(vv_k)\sim N(0,VPV^T) \end{cases} \end{aligned}
{P(wwk−1)∼N(0,WQWT)P(vvk)∼N(0,VPVT)
六、kalman与扩展kalman 公式对比
方法 | kalman 公式 | 扩展kalman 公式 |
---|---|---|
先验估计 | x ^ k − = A x ^ k − 1 + B u k − 1 \hat{x}_{k}^-=A \hat{x}_{k-1}+B u_{k-1} x^k−=Ax^k−1+Buk−1 | X ^ k − = f ( x k − 1 , u k − 1 , 0 ) \hat{X}_k^-=f(x_{k-1},u_{k-1},0) X^k−=f(xk−1,uk−1,0) |
先验误差协方差矩阵 | P k − = A P k − 1 A T + Q P_{k}^{-}=A P_{k-1} A^{T}+Q Pk−=APk−1AT+Q | P k − = A P k − 1 A T + W Q W T P_{k}^{-}=A P_{k-1} A^{T}+WQW^T Pk−=APk−1AT+WQWT |
卡尔曼增益 | K k = P k − H T H P k − H T + R K_{k}=\frac{P_{k}^{-} H^{T}}{H P_{k}^{-} H^{T}+R} Kk=HPk−HT+RPk−HT | K k = P k − H T H P k − H T + V R V T K_{k}=\frac{P_{k}^{-} H^{T}}{H P_{k}^{-} H^{T}+VRV^T} Kk=HPk−HT+VRVTPk−HT |
后验估计 | x ^ k = x ^ k − + K k ( z k − H x ^ k − ) \hat{x}_{k}=\hat{x}_{k}^{-}+K_{k}\left(z_{k}-H \hat{x}_{k}^{-}\right) x^k=x^k−+Kk(zk−Hx^k−) | x ^ k = x ^ k − + K k ( z k − h ( x ~ k , 0 ) ) \hat{x}_{k}=\hat{x}_{k}^{-}+K_{k}\left(z_{k}-h(\tilde{x}_k,0)\right) x^k=x^k−+Kk(zk−h(x~k,0)) |
更新误差协方差矩阵 | P k = ( I − K k H ) P k − P_{k}=\left(I-K_{k} H\right) P_{k}^{-} Pk=(I−KkH)Pk− | P k = ( I − K k H ) P k − P_{k}=\left(I-K_{k} H\right) P_{k}^{-} Pk=(I−KkH)Pk− |
参考文献
【卡尔曼滤波器】6_扩展卡尔曼滤波器_Extended Kalman Filter
线性化
在扩展卡尔曼滤波中,我们并不用前一个时刻的先验值 X k - X_k^- Xk-(卡尔曼滤波器未经过修正的预测值)作为参考点,而是用前一个时刻的估计值作为参考点做线性化。这是因为相对于先验值,前一个时刻的估计值更加贴近于真实值,将估计值作为线性化参考点可以得到一个更加贴近于实际的线性化系统模型。 ↩︎