Technical report from Automatic Control at Linköpings universitet
On parametric smoothing Cramér-Rao bounds
Carsten Fritsche, Umut Orguner
Division of Automatic Control
E-mail: carsten@isy.liu.se, umut@metu.edu.tr
23rd March 2017
Report no.: LiTH-ISY-R-3097
Address:
Department of Electrical Engineering Linköpings universitet
SE-581 83 Linköping, Sweden
WWW: http://www.control.isy.liu.se
AUTOMATIC CONTROL REGLERTEKNIK LINKÖPINGS UNIVERSITET
Technical reports from the Automatic Control group in Linköping are available from http://www.control.isy.liu.se/publications.
1
On parametric smoothing Cramér-Rao bounds
Carsten Fritsche and Umut Orguner
††
Department of Electrical & Electronics Engineering, Middle East Technical University, 06800, Ankara, Turkey
1
Parametric Filtering Cramér-Rao bound
We are interested in deriving a parametric CRLB for unbiased estimators of the current state xk, for the
conditional MSE matrix
M(ˆx0:k(y0:k)|x0:k) ≥ [J0:k(x0:k)]−1, (1)
where ˆxT
0:k(y0:k) = [ˆxT0(y0:k), . . . , ˆxTk(y0:k)] is any unbiased estimator of the state sequence x0:k, and J0:k(x0:k)
is the auxiliary Fisher information matrix of the state sequence x0:k, defined as
J0:k(x0:k) = Ey0:k,z1:k{−∆
x0:k
x0:klog(p(y0:k, z1:k|x0:k))|x0:k}. (2)
With the introduction of the auxiliary Fisher information matrix it is possible to incorporate information from the deterministic state dynamics into the bound calculations. In particular, since the state sequence x0:k is
deterministic, it can be rewritten as a set of equality constraints: 0 = xi+1− fi(xi) − Giwi, i = 0, . . . , k − 1.
These again can be interpreted as “perfect” measurements zi+1 = xi+1− fi(xi) − Giwi with zi+1 = 0 ∀i. In
order to stay in a probabilistic framework, we add to these equations zero-mean Gaussian noise with covariance
M, which is later on set to zero to recover the equality constraint. Thus, it is possible to establish the following recursion
p(y0:k, z1:k|x0:k) = p(yk|xk)p(zk|xk, xk−1)
× p(y0:k−1, z1:k−1|x0:k−1). (3)
The auxiliary Fisher information matrix J0:k−1(x0:k−1) can be then partitioned as follows:
J0:k−1(x0:k−1) = −Ey0:k−1, z1:k−1 ( ∆x0:k−2 x0:k−2 ∆ xk−1 x0:k−2 ∆x0:k−2 xk−1 ∆ xk −1 xk−1 log(p(y0:k−1, z1:k−1|x0:k−1)) x0:k−1 ) = A11 k−1 A12k−1 (A12 k−1)T A22k−1 . (4)
The auxiliary Fisher information submatrix Jk−1(x0:k−1) is computed as the inverse of the n × n lower-right
partition of [J0:k−1(x0:k−1)]−1, given by Jk−1(x0:k−1) = A11k−1− (A 12 k−1)T(A 22 k−1)−1A 12 k−1. (5)
Similarly, by making use of the recursion (3) it can be easily verified that the matrix J0:k(x0:k) can be simplified
as follows J0:k(x0:k) = −Ey0:k, z1:k ∆x0:k−2 x0:k−2 ∆ xk−1 x0:k−2 ∆ xk x0:k−2 ∆x0:k−2 xk−1 ∆ xk−1 xk−1 ∆ xk xk−1 ∆x0:k−2 xk ∆ xk−1 xk ∆ xk xk log(p(y0:k, z1:k|x0:k)) x0:k = A11 k−1 A12k−1 0 (A12 k−1)T A22k−1+ L11k L 12 k 0 (L12 k ) T L22 k (6) with elements L11k = Ey0:k,z1:k{−∆ xk−1 xk −1log(p(zk|xk, xk−1))|x0:k}, (7a) L12k = Ey0:k,z1:k{−∆ xk xk−1log(p(zk|xk, xk−1))|x0:k}, (7b) L22 k = Ey0:k,z1:k{−∆ xk xklog(p(zk|xk, xk−1))|x0:k} + Ey0:k,z1:k{−∆ xk xklog(p(yk|xk))|x0:k}. (7c)
Since J0:k(x0:k) has a block tri-diagonal structure, a recursive computation of the (auxiliary) Fisher information
submatrix Jk(x0:k) is possible. By noting that this submatrix is computed as the inverse of the n × n lower-right
partition of [J0:k(x0:k)]−1, block matrix inversion yields
Jk(x0:k) = L22k −0 (L 12 k ) T A11 k−1 A12k−1 (A12 k−1)T A22k + L 11 k −1 0 L12 k = L22 k − L 21 k Ak−1+ L11k − (A 12 k ) T(A11 k )−1A 12 k−1 −1 L12k . (8)
Inserting the definition of Jk(x0:k) given in (5) into (11), yields the desired recursive formula for computing
Jk(x0:k) = L22k − L 21 k L 11 k + Jk(x0:k) −1 L12k . (9)
For nonlinear additive Gaussian systems, the pdfs required to evaluate the L-terms in (7) are given by p(zk|xk, xk−1) =
N (zk; xk− fk−1(xk−1) − Gk−1wk−1, M) and p(yk|xk) = N (yk; hk(xk), Rk), respectively. The corresponding
L−terms then simplify to
L11k = F T k−1(xk−1)M−1Fk−1(xk−1), (10a) L12k = F T k−1(xk−1)M−1, (10b) L22k = H T k(xk)R−1k Hk(xk) + M−1, (10c)
where we have introduced the Jacobian matrices Fk−1(xk−1) = [∇xk−1f
T
k−1(xk−1)]T and Hk(xk) = [∇xkh
T k(xk)]T.
The auxiliary Fisher information matrix recursion then simplifies to
Jk(x0:k) = HkT(xk)R−1k Hk(xk) + M−1− M−1Fk−1(xk−1)FkT−1(xk−1)M−1Fk−1(xk−1) + Jk(x0:k) −1
× FT
k−1(xk−1)M−1. (11)
Applying the matrix inversion lemma to the above expression yields
Jk(x0:k) = (M + Fk−1(xk−1)Jk−1−1(x0:k−1)F
T
k−1(xk−1))−1+ Hk(xk)TR−1k Hk(xk). (12)
Finally setting M = 0 gives the desired recursion
Jk(x0:k) = (Fk−1(xk−1)Jk−1−1(x0:k−1)FkT−1(xk−1))−1+ Hk(xk)TR−1k Hk(xk). (13)
2
Parametric Smoothing Cramér-Rao bound
The presentation of the derivation of the parametric smoothing Cramér-Rao bound in this technical report follows closely the derivation of the posterior smoothing Cramér-Rao bound presented in [1, 2].
The parametric smoothing Cramér-Rao bound for a smoothing state estimate ˆxℓ(y1:k), 0 ≤ ℓ ≤ k − 1, is based
on the auxiliary FIM J0:k(x0:k). Note, that [J0:k(x0:k)]−1 is the lower bound of the conditional MSE matrix
M(ˆx0:k(y0:k)|x0:k) of the estimate trajectory consisting of the smoothing estimates ˆx0T(y0:k), . . . , ˆxTk−1(y0:k),
and the filtering estimate ˆxTk(y0:k). Thus, smoothing Cramér-Rao bounds should be implicitly contained in
J0:k(x0:k). We decompose J0:k(x0:k) into blocks which correspond to time instances 0, . . . , ℓ and ℓ+1, ℓ+2, . . . , k,
yielding J0:k(x0:k) ∆ = Γℓ|ℓ Kℓ KT ℓ Lℓ|k = A11 ℓ A12ℓ A12ℓ T A22ℓ + L 11 ℓ+1 L 12 ℓ+1 L12 ℓ+1 T L22 ℓ+1+ L11ℓ+2 L12ℓ+2 L12 ℓ+2 T . .. . .. . .. L22 k−1+ L11k L12k L12 k T L22 k (14)
3
The inverse [J0:k(x0:k)]−1 contains the smoothing bounds Pℓ|k, ℓ = 0, 1, . . . , k − 1 and the filtering bound Pk|k
on the main diagonal
[J0:k(x0:k)] −1 = P0|k . .. Pℓ|k Pℓ+1|k . .. Pk|k ∆ = " [J0:k(x0:k)]−11,1 [J0:k(x0:k)]−11,2 [J0:k(x0:k)]−12,1 [J0:k(x0:k)]−12,2 # . (15)
In order to obtain the smoothing bound Pℓ|k, we have to extract the lower-right block of[J0:k(x0:k)]−1]
1,1 by Pℓ|k= [0, Inx][J0:k(x0:k)] −1 1,1[0, Inx] T, (16)
where Inx is an nx× nx identity matrix and 0 demotes a zero matrix of appropriate size. In order to obtain [J0:k(x0:k)]−11,1, block-matrix inversion will be applied to [J0:k(x0:k)]−1, yielding
[J0:k(x0:k)]−1 1,1= [Γℓ|ℓ] −1+ [Γ ℓ|ℓ]−1KℓLℓ|k− KℓT[Γℓ|ℓ]−1Kℓ −1 KℓT[Γℓ|ℓ]−1 = [Γℓ|ℓ]−1+ [Γℓ|ℓ]−1Kℓ[J0:k(x0:k)]−12,2KℓT[Γℓ|ℓ]−1, (17) with [Γℓ|ℓ]−1= " ∗ ∗ ∗ hA22 ℓ + L 11 ℓ+1− A 12 ℓ T A11 ℓ −1 A12 ℓ i−1 # (18)
Using the intermediate result [0, Inx][Γℓ|ℓ] −1[0, I nx] T=hA22 ℓ + L 11 ℓ+1− A 12 ℓ T A11 ℓ −1 A12ℓ i−1 =hPℓ−1|ℓ + L11 ℓ+1 i−1 , (19)
we can finally evaluate (16)
Pℓ|k= h Pℓ−1|ℓ + L11 ℓ+1 i−1 +hPℓ−1|ℓ + L11 ℓ+1 i−1 L12ℓ+1Pℓ+1|k L12ℓ+1 Th Pℓ−1|ℓ + L11 ℓ+1 i−1 (20)
Matrix inversion (A + BCD)−1 with A =hP−1
ℓ|ℓ + L 11 ℓ+1 i−1 , B = DT =h Pℓ−1|ℓ + L11 ℓ+1 i−1 L12ℓ+1 and C = Pℓ+1|k
gives the smoothing FIM Jℓ|k
Jℓ|k= Jℓ|ℓ+ L11ℓ+1− L 12 ℓ+1 h L12ℓ+1 T Jℓ|ℓ+ L11ℓ+1 −1 L12ℓ+1+ Jℓ+1|k i−1 L12ℓ+1 T (21) We replace L12 ℓ+1= F T k−1(xk−1)M−1 and L11ℓ+1= F T k−1(xk−1)M−1Fk−1(xk−1), yielding h L12ℓ+1 T Jℓ|ℓ+ L11ℓ+1 −1 L12ℓ+1+ Jℓ+1|k i−1 =hM−1Fk−1(xk−1)Jℓ|ℓ+ FkT−1(xk−1)M−1Fk−1(xk−1) −1 FkT−1(xk−1)M−1+ Jℓ+1|k i−1 =hJℓ+1|k+ M−1− h M−1− M−1Fk−1(xk−1)Jℓ|ℓ+ FkT−1(xk−1)M−1Fk−1(xk−1) −1 FkT−1(xk−1)M−1 ii−1 = Jℓ+1|k− h M+ Fk−1(xk−1)Jℓ|ℓ −1 FkT−1(xk−1) i−1 + M−1 −1 , (22)
where we again used the matrix inversion lemma. Similarly, we can rewrite (21)
Jℓ|k = Jℓ|ℓ+ FkT−1(xk−1) M−1− M−1h L12 ℓ+1 T Jℓ|ℓ+ L11ℓ+1 −1 L12ℓ+1+ Jℓ+1|k i−1 M−1 Fk−1(xk−1) (23)
Inserting (22) into (23) yields Jℓ|k= Jℓ|ℓ+ FkT−1(xk−1) " M−1− M−1 Jℓ+1|k− h M+ Fk−1(xk−1)Jℓ|ℓ −1 FkT−1(xk−1) i−1 + M−1 −1 M−1 # Fk−1(xk−1) = Jℓ|ℓ+ FkT−1(xk−1) " M + Jℓ+1|k− h M + Fk−1(xk−1)Jℓ|ℓ −1 FkT−1(xk−1) i−1−1 #−1 Fk−1(xk−1) (24)
Setting M = 0 gives the recursion for the parametric smoothing FIM
Jℓ|k= Jℓ|ℓ+ FkT−1(xk−1) Jℓ+1|k− h Fk−1(xk−1)Jℓ|ℓ −1 FkT−1(xk−1) i−1 Fk−1(xk−1) (25)
and its inverse Pℓ|k=Jℓ|k
−1
gives the parametric smoothing Cramér-Rao lower bound.
References
[1] N. Bergman, “Recursive Bayesian estimation: Navigation and tracking applications,” Ph.D. dissertation, Linköping University, Linköping, Sweden, 1999.
[2] M. Simandl, J. Královec, and P. Tichavský, “Filtering, predictive, and smoothing Cramér-Rao bounds for discrete-time nonlinear dynamic systems,” Automatica, vol. 37, no. 11, pp. 1703–1716, Nov. 2001.
Avdelning, Institution
Division, Department
Division of Automatic Control Department of Electrical Engineering
Datum Date 2017-03-23 Språk Language ✷ Svenska/Swedish ✷ Engelska/English ✷ ⊠ Rapporttyp Report category ✷ Licentiatavhandling ✷ Examensarbete ✷ C-uppsats ✷ D-uppsats ✷ Övrig rapport ✷ ⊠
URL för elektronisk version
http://www.control.isy.liu.se
ISBN
—
ISRN
—
Serietitel och serienummer
Title of series, numbering
ISSN
1400-3902
LiTH-ISY-R-3097
Titel
Title
On parametric smoothing Cramér-Rao bounds
Författare
Author
Carsten Fritsche, Umut Orguner
Sammanfattning
Abstract
In this report, the parametric Cramér-Rao lower bound for the smoothing problem is derived.
Nyckelord
