1 原理
通过对比残差协方差的在线估计值及其理论值来自适应调整过程噪声方差矩阵
Q
d
\mathbf{Q}_{\mathrm{d}}
Qd以及测量噪声方差矩阵
R
\mathbf{R}
R,直到残差协方差的估计值和理论值相等。
假定
R
\mathbf{R}
R已知,残差协方差的实际值为
A
^
=
1
N
∑
j
=
i
−
N
+
1
i
r
j
r
j
T
(1)
\widehat{\mathbf{A}}=\frac{1}{N} \sum_{j=i-N+1}^{i} \mathbf{r}_{j} \mathbf{r}_{j}^{\mathrm{T}}\tag{1}
A
=N1j=i−N+1∑irjrjT(1)
与其相对应的理论值
[
H
(
t
i
)
P
(
t
i
−
)
H
T
(
t
i
)
+
R
(
t
i
)
]
\left[\mathbf{H}\left(t_{i}\right) \mathbf{P}\left(t_{i}^{-}\right) \mathbf{H}^{\mathrm{T}}\left(t_{i}\right)+\mathbf{R}\left(t_{i}\right)\right]
[H(ti)P(ti−)HT(ti)+R(ti)]由卡尔曼滤波器提供。如果
A
^
\widehat{\mathbf{A}}
A
超过了滤波器推导的理论值(这个超过可以指特征值,对角元或者范数意义上),那么过程噪声
Q
d
\mathbf{Q}_{\mathrm{d}}
Qd就应该增加。具体方法如下。由
A
(
t
i
)
=
H
P
(
t
i
−
)
H
T
+
R
=
H
[
Φ
P
(
t
i
−
1
+
)
Φ
T
+
G
d
Q
d
G
d
T
]
H
T
+
R
\mathbf{A}\left(t_{i}\right)=\mathbf{H} \mathbf{P}\left(t_{i}^{-}\right) \mathbf{H}^{\mathrm{T}}+\mathbf{R}=\mathbf{H}\left[\mathbf{\Phi} \mathbf{P}\left(t_{i-1}^{+}\right) \mathbf{\Phi}^{\mathrm{T}}+\mathbf{G}_{\mathrm{d}} \mathbf{Q}_{\mathrm{d}} \mathbf{G}_{\mathrm{d}}^{\mathrm{T}}\right] \mathbf{H}^{\mathrm{T}}+\mathbf{R}
A(ti)=HP(ti−)HT+R=H[ΦP(ti−1+)ΦT+GdQdGdT]HT+R可得
H
G
d
Q
d
G
d
T
H
T
=
A
(
t
i
)
−
H
Φ
P
(
t
i
−
1
+
)
Φ
T
H
T
−
R
(2)
\mathbf{H G}_{\mathrm{d}} \mathbf{Q}_{\mathrm{d}} \mathbf{G}_{\mathrm{d}}^{\mathrm{T}} \mathbf{H}^{\mathrm{T}}=\mathbf{A}\left(t_{i}\right)-\mathbf{H} \mathbf{\Phi} \mathbf{P}\left(t_{i-1}^{+}\right) \mathbf{\Phi}^{\mathrm{T}} \mathbf{H}^{\mathrm{T}}-\mathbf{R}\tag{2}
HGdQdGdTHT=A(ti)−HΦP(ti−1+)ΦTHT−R(2)
具体求解
Q
d
\mathbf{Q}_{\mathrm{d}}
Qd可能需要借助伪逆。
如果
Q
d
\mathbf{Q}_{\mathrm{d}}
Qd已知,则
R
\mathbf{R}
R可以用下式进行估计
R
^
=
1
N
∑
j
=
i
−
N
+
1
i
r
j
r
j
T
−
H
P
(
t
i
−
)
H
T
(3)
\hat{\mathbf{R}}=\frac{1}{N} \sum_{j=i-N+1}^{i} \mathbf{r}_{j} \mathbf{r}_{j}^{\mathrm{T}}-\mathbf{H} \mathbf{P}\left(t_{i}^{-}\right) \mathbf{H}^{\mathrm{T}}\tag{3}
R^=N1j=i−N+1∑irjrjT−HP(ti−)HT(3)
2. 应用
其中一种方法是把待估计的过程噪声方差矩阵建模为
Q
d
=
Q
d
0
+
a
[
Δ
Q
d
]
(4)
\mathbf{Q}_{\mathrm{d}}=\mathbf{Q}_{\mathrm{d} 0}+a\left[\mathbf{\Delta Q}_{\mathrm{d}}\right]\tag{4}
Qd=Qd0+a[ΔQd](4)
其中, Q d 0 \mathbf{Q}_{\mathrm{d} 0} Qd0是标称值, Δ Q d \Delta \mathbf{Q}_{\mathrm{d}} ΔQd是围绕标称值的小量, a a a是待调整参数。如果 a a a有有限个备选选项,则该问提会被转化为一个多模型滤波问题。不同滤波器的 a a a不同,结果选残差特性与理论值符合最好的 a a a所对应的滤波器结果。
心得,评注
- (1)式计算残差的实际协方差时,假定残差的期望是零。参考文献[2]则去实际计算残差的期望,因此更符合真实实际情况。其具体公式为
W ^ k = 1 k ∑ i = 1 k e ~ i e ~ i T , with e ~ i = e i − e ‾ i (5) \hat{\boldsymbol{W}}_{k}=\frac{1}{k} \sum_{i=1}^{k} \tilde{\boldsymbol{e}}_{i} \tilde{\boldsymbol{e}}_{i}^{T}, \quad \text { with } \quad \tilde{\boldsymbol{e}}_{i}=\boldsymbol{e}_{i}-\overline{\boldsymbol{e}}_{i}\tag{5} W^k=k1i=1∑ke~ie~iT, with e~i=ei−ei(5)
其中
e
ˉ
i
=
1
i
∑
j
=
1
i
e
j
(6)
\bar{e}_{i}=\frac{1}{i} \sum_{j=1}^{i} e_{j}\tag{6}
eˉi=i1j=1∑iej(6)
为 i i i时刻残差期望的估计值。公式(6)表明,对于 i i i时刻残差期望的估计,使用== i i i时刻及以前的所有残差的平均值==。
参考文献
- Stochastic models, estimation, and control VOLUME 2 10.10小节
- Motion and Parameter Estimation of Space Objects Using Laser-Vision Data