Automatica 46 (2010) 708–715
Contents lists available at ScienceDirect
Automatica journal homepage: www.elsevier.com/locate/automatica
Brief paper
On the global optimality of unbiased minimum-variance state estimation for systems with unknown inputsI Chien-Shu Hsieh ∗ Department of Electrical Engineering, Ta Hwa Institute of Technology, Hsinchu 30740, Taiwan, ROC
article
info
Article history: Received 27 April 2009 Received in revised form 29 September 2009 Accepted 10 January 2010 Available online 23 February 2010 Keywords: Global optimality Optimal estimation Unknown inputs Unbiased minimum-variance estimation
abstract In this paper, a globally optimal filtering framework is developed for unbiased minimum-variance state estimation for systems with unknown inputs that affect both the system state and the output. The resulting optimal filters are globally optimal within the unbiased minimum-variance filtering over all linear unbiased estimators. Globally optimal state estimators with or without output and/or input transformations are derived. Through the global optimality evaluation of this research, the performance degradation of the filter proposed by Darouach, Zasadzinski, and Boutayeb [Darouach, M., Zasadzinski, M., & Boutayeb, M. (2003). Extension of minimum variance estimation for systems with unknown inputs. Automatica, 39, 867–876] is clearly illustrated and the global optimality of the filter proposed by Gillijns and De Moor [Gillijns, S., & De Moor, B. (2007b). Unbiased minimum-variance input and state estimation for linear discrete-time systems with direct feedthrough. Automatica, 43, 934–937] is further verified. The relationship with the existing literature results is addressed. A unified approach to design a specific globally optimal state estimator that is based on the desired form of the distribution matrix from the unknown input to the output is also presented. A simulation example is given to illustrate the proposed results. © 2010 Elsevier Ltd. All rights reserved.
1. Introduction Unknown input filtering (UIF) has played a significant role in many applications, e.g., bias compensation (Friedland, 1969; Hsieh & Chen, 1999; Ignani, 1990), geophysical and environmental applications (Kitanidis, 1987), fault detection and isolation problems (Chen & Patton, 1996), and functional filtering (Hsieh, 2007). Except for the above-mentioned first application, the UIF problem is often solved by assuming that no prior information about the unknown input is available, and this is the main concern of this paper. A common approach for solving this more general UIF problem for stochastic systems is to produce the unknown input decoupled state estimation. Researchers devoted to this research area have proposed various different kinds of optimal state estimators (Borisov & Pankov, 1994; Chen & Patton, 1996; Cheng, Ye, Wang, & Zhou, 2009; Darouach & Zasadzinski, 1997; Darouach, Zasadzinski,
I This work was supported by the National Science Council, Taiwan under Grant NSC 97-2221-E-233-003. The material in this paper was partially presented at CACS 2009 International Automatic Control Conference, Nov. 27–29, Taipei, Taiwan. This paper was recommended for publication in revised form by Associate Editor Tongwen Chen under the direction of Editor Ian R. Petersen. The author would like to thank the Associate Editor and the anonymous referees for their insightful comments and suggestions. ∗ Corresponding address: Ta Hwa Institute of Technology, 1, Ta Hwa Road, Chiunglin, Hsinchu 30740, Taiwan, ROC. Tel.: +886 35927700; fax: +886 35927085. E-mail addresses:
[email protected],
[email protected].
0005-1098/$ – see front matter © 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.automatica.2010.01.029
& Boutayeb, 2003; Gillijns & De Moor, 2007a,b; Hou & Patton, 1998; Hsieh, 2000, 2009, in press; Kitanidis, 1987; Valcher, 1999). However, although some insights have been gained (Hsieh, 2009, in press), the relationship between the various optimal filters remains unclear, and the global optimality of existing optimal filters has not been fully addressed. All the above-mentioned issues have motivated this research. It is well known that the globally optimal state estimator for stochastic systems without direct feedthrough of the unknown input to the output can be obtained as given in Darouach and Zasadzinski (1997), Gillijns and De Moor (2007a), Hsieh (2000), Kerwin and Prince (2000), and Kitanidis (1987). On the other hand, the results concerning the global optimality of the optimal state estimator for systems with direct feedthrough of the unknown input to the output are limited. To the best of our knowledge, Darouach et al. (2003) first addressed this issue in their proposed optimal estimator filter (OEF) and optimal predictor filter (OPF). Specifically, a decorrelation constraint for the existence of the optimal filter design was proposed. However, a performance degradation problem of the OEF was observed in our previous study (Hsieh, 2006a). In trying to derive a globally optimal state estimator, we proposed the optimal unbiased minimum-variance filter (OUMVF) (Hsieh, 2006b) as a potential remedy for the problem encountered in the OEF. Unfortunately, only a restricted solution was obtained due to the decorrelation constraint proposed by Darouach et al. (2003). Recently, Gillijns and De Moor (2007b) proposed a recursive threestep filter (RTSF) to potentially solve the problem, having found
C.-S. Hsieh / Automatica 46 (2010) 708–715
that, through a joint input and state estimation, the optimal state estimator could be obtained in a less restricted way as compared to the OEF. However, the RTSF is only applicable when the distribution matrix from the unknown input to the output is of full column rank due to the rank condition for the existence of an unbiased input estimator; moreover, the global optimality of the RTSF was not addressed. An extension of the RTSF (Hsieh, 2009) which took into consideration the possibly biased input estimate, i.e., by trying to remove the restriction of the full rank of the distribution matrix, and the new modified unbiased minimum-variance filter (Hsieh, in press), which was intended to improve the filtering performance of the OUMVF, have also been proposed. However, the global optimality of these filters was not addressed. More recently, Cheng et al. (2009) addressed the global optimality issue in their proposed recursive optimal filter (ROF). They showed that the ROF is indeed a globally optimal solution to this particular problem. Nevertheless, it should be stressed that the ROF was obtained through output and unknown input transformations by applying the singular value decomposition (SVD) to the distribution matrix, which rendered the obtained gain matrices to satisfy the decorrelation constraint. Thus, a globally optimal solution without resorting to the decorrelation constraint, i.e., one requiring neither output and/or unknown input transformations, for systems with direct feedthrough of the unknown input to the output still remained open. The aims of this paper are as follows: (1) to present new globally optimal state estimators via the modifications of the OUMVF; (2) to illustrate the essence of the globally optimal state estimator within the unbiased minimum-variance filtering framework; (3) to address the global optimality of the new derived optimal filters in the sense of unbiased minimum variance over all linear unbiased estimators; and (4) to present a unified approach to design a specific globally optimal state estimator for systems with unknown inputs. Through the proposed results, the performance degradation problem of the OEF and the suboptimal filtering performance of the OUMVF have been clearly illustrated. Moreover, the proof of the global optimality of the RTSF has been provided. The relationship with the existing literature results has also been addressed. This paper is organized as follows. In Section 2, the statement of the problem is addressed. Section 3 recalls the design method of the OUMVF, through which a globally optimal state estimator, known as the global OUMVF (GOUMVF), can be constructed, and the relationship between the GOUMVF and the ROF is addressed. In Section 4, a new globally optimal state estimator without output or unknown input transformations, called the refined OUMVF (ROUMVF), is derived. The relationship between the ROUMVF and the existing literature results is also addressed and the global optimality of the ROUMVF is evaluated. Moreover, a compact form of the ROUMVF is presented. Section 5 summarizes a unified approach to design a specific globally optimal state estimator for a system of interest. An illustrative example is given to show the effectiveness of the proposed results. Section 6 has the conclusions. 2. Problem formulation Consider the following linear discrete-time stochastic timevarying system: xk+1 = Ak xk + Gk dk + wk ,
(1)
yk = Ck xk + Hk dk + vk ,
(2)
where xk ∈ R , dk ∈ R , yk ∈ R , wk ∈ R , and vk ∈ R are the state vector, deterministic unknown input, output, process noise, and measurement noise, respectively. Ak , Gk , Ck , Hk are deterministic known matrices with appropriate dimensions. wk and vk are uncorrelated white noises with covariance matrices Qk = E [wk wkT ] ≥ 0 and Rk = E [vk vkT ] > 0, respectively. The initial state x0 is with mean xˆ 0 and covariance P0 , and is independent of wk and vk . To n
m
p
n
p
709
simplify the discussions and without loss of generality, the known input vector uk is not included in the system model. The problem is to design a globally optimal state estimator of xk in the sense of unbiased minimum variance over all linear unbiased estimators for system (1)–(2) given measurements up to time instant k and without any information concerning the unknown input dk . The filter considered in this paper is given by the following (full-order) modified unbiased minimum-variance filter (MUMVF) (Hsieh, in press): xˆ k = Nk xˆ k−1 + Kk y˜ k−1 + Lk yk ,
(3)
where Nk , Kk , Lk , and y˜ k−1 are determined to achieve E [ek ] = 0 and tr (E [ek eTk ]) is minimized, in which ek = xk − xˆ k . Note that if vector y˜ k−1 in (3) is determined as yk−1 , then the MUMVF (3) reduces to the conventional full-order stochastic Luenberger observer-type filter considered in Cheng et al. (2009), Darouach et al. (2003) and Hsieh (2006b). The main purpose of this paper is to explore the globally optimal solution to the aforementioned state estimation problem in light of the filter structure (3). Two categories of optimal filter solutions are presented: one is designed for y˜ k−1 = yk−1 , and the other for y˜ k−1 6= yk−1 , which are addressed in Sections 3 and 4, respectively. 3. Optimal state estimator design for y˜ k −1 = yk −1 3.1. The optimal unbiased minimum-variance filter The previously proposed OUMVF (Hsieh, 2006b) for system (1)–(2) is obtained via consideration of the following filter form: xˆ k = Nk xˆ k−1 + Kk yk−1 + Lk yk ,
(4)
where Nk = (I − Lk Ck )Ak−1 − Kk Ck−1 ,
(5)
and the gain matrices Kk and Lk are determined to satisfy the following unbiasedness conditions:
(I − Lk Ck )Gk−1 − Kk Hk−1 = 0,
Lk Hk = 0.
(6)
As shown in Darouach et al. (2003), the optimal filter design in light of the filter structure (4)–(6) can be achieved by finding the gain matrices Lk and Kk+1 in satisfying the following decorrelation constraint: Nk+1 Lk Rk KkT+1 = 0.
(7)
Based on (7), two designs of the OUMVF were proposed in Darouach et al. (2003): one is the OEF design, and the other the OPF design, which are obtained by considering Kk+1 = 0 and Lk = 0, respectively. It should be stressed that these results are obtained regardless of the values of Nk+1 and Rk . Moreover, it is known that, for the estimator filter design problem, the choice with Kk+1 = 0 may not yield an optimal state estimator due to the greater restriction of the obtained unbiasedness constraint (Hsieh, 2009). In fact, a less restrictive solution for the OUMVF design can be obtained via finding sufficient conditions for (7) to be satisfied only for all Nk+1 , i.e., in trying to achieve Lk Rk KkT+1 = 0. This is the case that has been overlooked by Darouach et al. (2003) in deriving their proposed optimal filters. By exploring the relationship between the gain matrices Kk and Lk and treating the former as a design parameter, the following OUMVF was proposed (Hsieh, 2006b):
¯k = 0 G
U¯ k−1
Hk
T¯k = αk I − Hk
R˜ k = Ck P˜ k CkT + Rk ,
Ck U¯ k−1
Ck U¯ k−1
Ď
Hk
,
(8) Ck U¯ k−1
Ď
,
(9) (10)
710
C.-S. Hsieh / Automatica 46 (2010) 708–715
¯ k R˜ k )T¯kT (T¯k R˜ k T¯kT )−1 , K¯ k = (P˜ k CkT − G
(11)
¯ k + K¯ k T¯k , Lk = G
(12)
xˆ k = x¯ k + Lk (yk − Ck x¯ k ),
(13)
Pk = (I − Lk Ck )P˜ k (I − Lk Ck )T + Lk Rk LTk ,
(14)
x¯ k+1 = Ak xˆ k + Λk (yk − Ck xˆ k ), P˜ k+1 =
Φk Pk ΦkT Ď
+ Qk +
Λk Rk ΛTk
(15)
+
Φk Lk Rk ΛTk
+
Λk Rk LTk ΦkT
, (16)
where M is the Moore–Penrose pseudo-inverse of M, matrix parameter αk must be chosen so that T¯k is of full row rank, and U¯ k = Gk − Λk Hk ,
Λk =
Φk = Ak − Λk Ck ,
(17)
Ck U¯ k−1 = rank
0 Hk
U¯ k−1 Ck U¯ k−1
= rank[Hk ] + rank[U¯ k−1 ].
(19)
Note that if matrix Λk in (18) is given by Λk = 0, i.e., the design parameter Kk+1 is chosen as Kk+1 = 0, then the obtained OUMVF will be equivalent to the OEF. Ď For the case of Kk+1 6= 0, since Kk+1 = (I − Lk+1 Ck+1 )Gk Hk and Lk is in the left null space of Hk due to (6), a sufficient condition to achieve Lk Rk KkT+1 = 0 can be given as follows:
(I − Hk HkĎ )Rk (HkĎ )T = 0,
(20)
which may not generally hold, and hence the OUMVF with Ď Λk = Gk Hk may not be an optimal state estimator. However, in Section 3.2 we shall show that the OUMVF can be slightly tailored to derive a globally optimal state estimator. Remark 1. A necessary requirement for the rank condition (19) to hold is that the columns of Hk and U¯ k−1 must be independent. For the OEF, which can be viewed as a special case of the OUMVF, rank condition (19) becomes
rank Hk
Ck Gk−1 = rank[Hk ] + rank[Gk−1 ].
(21)
3.2. A global OUMVF design Ď Gk Hk
As stated in Section 3.1, the OUMVF with Λk = may become an optimal filter if condition (20) holds. Specifically, if matrix Rk can be chosen as Rk = ρk Ip×p , then (20) always holds. This dependence on the system covariance matrix may be too restrictive to be applied in practical applications. In the following, a more practical design to achieve (20) is presented, whereby one can derive a specific globally optimal state estimator. The main idea behind deriving a globally optimal state estimator in light of the filter structure (4) is to find the specific forms of Hk and Rk such that condition (20) is satisfied. In showing this, the following specific Hk and Rk are considered:
Hk =
¯k Σ 0
0 , 0
Rk =
R1,k 0
0 , R2,k
(22)
¯ k is of full row rank, R1,k > 0, and R2,k > 0. Using (22) where Σ
0 = 0
0 I
R1,k 0
0 R2,k
Ď T ¯k) (Σ 0
0 = 0. 0
the measurement and the unknown input vector, respectively, such that the original system (1)–(2) can be rewritten as the following equivalent system: xk+1 = Ak xk + (Gk Tkd )d¯ k + wk ,
(23)
(T˜k0 yk ) = (T˜k0 Ck )xk + (T˜k0 Hk Tkd )d¯ k + (T˜k0 vk ),
(24)
where
¯ ˜Tk0 Hk Tkd = Σk
¯ ,
Tkd dk
0 . 0
0
(25)
Normally, T˜k0 and Tkd can be determined by renumbering the output equation and the unknown input vector, respectively. Otherwise, it can always be done by applying the SVD to distribution matrix Hk , as shown in Cheng et al. (2009). Then, define the following notations: G2,k ,
Gk Tkd = G1,k
d¯ k =
,k ˜Tk0 Rk (T˜k0 )T = V11 T V12,k
1 d¯ k , d¯ 2k
V12,k . V22,k
(26)
Since Rk > 0 and T˜k0 6= 0, there always exists the following nonsingular transformation matrix T˜k1 : T˜k1 =
I 0
−1 −V12,k V22 ,k
I
,
(27)
such that the equivalent system (23)–(24) can be further represented as follows: xk+1 = Ak xk + G1,k d¯ 1k + G2,k d¯ 2k + wk ,
(28)
Σk d1k
z1,k = C1,k xk + ¯ ¯ + v1,k ,
(29)
z2,k = C2,k xk + v2,k ,
(30)
where z1,k = T˜k1 T˜k0 yk , z2,k
Comparing (21) with (19), it is clear that the existence of the OEF may be more stringent than that of the OUMVF since the columns of Hk and Gk−1 are more likely to be dependent than those of Hk and U¯ k−1 due to rank[Gk ] ≥ rank[U¯ k ]. This illustrates that the OUMVF may have better filtering performance than the OEF (Hsieh, in press). An illustrative example is given in Section 5.
Let T˜k0 and Tkd be two transformation matrices associated with
(18)
The existence of the above OUMVF is given by the following rank condition: rank Hk
(I −
Ď Ď Hk Hk )Rk (Hk )0
dk =
for Kk+1 = 0 for Kk+1 6= 0.
0 Ď Gk Hk
in (20) yields
C1,k C2,k
v1,k = T˜k1 T˜k0 vk , v2,k
= T˜k1 T˜k0 Ck ,
and v1,k and v2,k are uncorrelated and, respectively, have the following variances: −1 T R1,k = V11,k − V12,k V22 R2,k = V22,k . (31) ,k V12,k , Finally, applying the OUMVF (8)–(18) to the equivalent system (28)–(30), we obtain the following GOUMVF:
¯ k = U¯ k−1 (C2,k U¯ k−1 )Ď , G U¯ k = G1,k (I −
¯ k) ¯ kĎ Σ Σ
(32) G2,k ,
(33)
¯ k ), T¯k = αk (I − C2,k G ¯k + L 2 ,k = G
× T¯kT
(34)
˜
−¯ (
+ R2,k ) −1 T¯k , T¯k (C2,k P˜ k C2T,k + R2,k )T¯kT Pk C2T,k
Gk C2,k Pk C2T,k
˜
(35)
xˆ k = x¯ k + L2,k (z2,k − C2,k x¯ k ),
(36)
Pk = (I − L2,k C2,k )P˜ k (I − L2,k C2,k )T + L2,k R2,k LT2,k ,
(37)
x¯ k+1 = Ak xˆ k + P˜ k+1 =
¯ kĎ (z1,k G1,k Σ
− C1,k xˆ k ),
¯ kĎ C1,k )Pk (Ak − G1,k Σ ¯ kĎ C1,k )T (Ak − G1,k Σ ¯ kĎ R1,k (G1,k Σ ¯ kĎ )T , + Qk + G1,k Σ
(38)
(39)
C.-S. Hsieh / Automatica 46 (2010) 708–715
where (36) can also be represented in light of the filter structure (4) by choosing: Ď
¯ k−1 C1,k−1 ), Nk = (I − L2,k C2,k )(Ak−1 − G1,k−1 Σ ¯ kĎ−1 I Kk = (I − L2,k C2,k )G1,k−1 Σ 1 0 Lk = L2,k 0 I T˜k T˜k .
rank C2,k U¯ k−1 = rank
U¯ k−1 C2,k U¯ k−1
(49)
where Ď
= rank[U¯ k−1 ].
Remark 2. The above GOUMVF is a globally optimal state estimator, as highlighted in Remark 8, which has a strong connection to ¯ k is nonsingular, the ROF designed by Cheng et al. (2009). If matrix Σ then U¯ k in (33) will reduce to U¯ k = G2,k . Moreover, if the following transformations are used, T˜k0 = UkT ,
Then, we are in place to determine matrix Xk . Using the constraints in (47), we chose Xk as follows: Xk = X¯ k T¯k ,
0 T˜k1−1 T˜k0−1 ,
The existence of the GOUMVF is given by the following rank condition due to (19):
711
Tkd = Vk ,
T¯k = αk (I − Sk Sk ),
Sk = Hk
Ck U¯ k∗−1 .
Using (44)–(50), we have P¯ k∗+1 (≡ E [˜ek+1 e˜ Tk+1 ])
= Aˆ k Pk Aˆ Tk + Qk + Gk HkĎ Rk (Gk HkĎ )T + X¯ k T¯k R˜ ∗k T¯kT X¯ kT +Aˆ k (L∗k R˜ ∗k − P¯k∗ CkT )T¯kT X¯ kT + Gk HkĎ Rk T¯kT X¯ kT + X¯ k T¯k (L∗k R˜ ∗k − P¯k∗ CkT )T Aˆ Tk + X¯ k T¯k Rk (Gk HkĎ )T + Aˆ k L∗k Rk (Gk HkĎ )T + Gk HkĎ Rk (Aˆ k L∗k )T ,
(51)
where R˜ ∗k = Ck P¯ k∗ CkT + Rk .
(52)
where Uk and Vk are unitary matrices obtained by applying the SVD to the distribution matrix Hk , i.e., Hk = Uk diag{Σk , 0}VkT , where Σk is a diagonal matrix, then the above GOUMVF will become the ROF.
Solving (51) for X¯ k that minimizes the trace of P¯ k∗+1 , we obtain
4. Optimal state estimator design for y˜ k −1 6= yk −1
where we assume
Although the OUMVF serves as an effective filter structure to derive globally optimal state estimators, as given in the previous section, some transformations are needed to facilitate the derivation. In this section, we have shown how to refine the OUMVF to relax those transformations in deriving a new globally optimal state estimator.
(L∗k R˜ ∗k − P¯k∗ CkT )T¯kT X¯ kT = 0,
xˆ k = x¯ k + Lk (yk − Ck x¯ k ), x¯ ∗k+1
∗
= Ak xˆ k +
∗
Ď Gk Hk (yk
(40)
− Ck xˆ k ) + Xk y˜ k ,
(41)
where y˜ k = yk − Ck x¯ ∗k and L∗k and Xk are matrices to be determined. Defining x˜ k+1 = x¯ ∗k+1 + U¯ k∗ dk ,
e˜ k+1 = xk+1 − x˜ k+1 ,
(42)
y˜ k+1 = Ck+1 e˜ k+1 + Hk+1 dk+1 + Ck+1 U¯ k∗ dk + vk+1 . (43) Next, we derive the unbiasedness conditions of the state estimator (40)–(41). Using (1)–(2) and (40)–(43), we have ek = (I − L∗k Ck )(˜ek + U¯ k∗−1 dk−1 ) − L∗k vk − L∗k Hk dk ,
(44)
Ď Gk Hk vk
− Xk (Ck e˜ k + Hk dk + Ck U¯ k∗−1 dk−1 + vk ),
(45)
Ď
where Aˆ k = Ak − Gk Hk Ck . From (44) and (45), we obtain the following unbiasedness conditions: Lk Hk = 0,
(I −
Xk Hk = 0,
Xk Ck U¯ k∗−1 = 0.
∗
Lk Ck )U¯ k∗−1 ∗
= 0,
∗ Ck U¯ k∗−1 = rank[Hk ] + rank[U¯ k−1 ].
(55)
) ¯ (T¯k Rk Tk ) Tk ,
(56)
Pk = (I − L∗k Ck )P¯ k∗ (I − L∗k Ck )T + L∗k Rk (L∗k )T ,
(57)
x¯ ∗k+1
¯∗
¯∗ T
= Ak xˆ k +
P¯ k∗+1 =
¯ ∗ ˜∗
˜ ∗ ¯ T −1 ¯
Gk Rk TkT
Ď Gk Hk (yk
− Ck xˆ k ) + w ˜ k,
(58)
Ď Ď Aˆ k Pk Aˆ Tk + Qk + Q˜ k + Gk Hk Rk (Gk Hk )T + Aˆ k L∗k Rk (HkĎ )T GTk + Gk HkĎ Rk (L∗k )T Aˆ Tk ,
(59)
where Ď
¯ ∗k = Γk Sk , G w ˜k =
U¯ k∗−1 ,
Γk = 0
−Gk HkĎ Rk T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k (yk
(60)
− Ck x¯ k ), ∗
(61)
Q˜ k = −Gk Hk Rk T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k Rk (Gk Hk )T .
(62)
Ď
Ď
Using (53) and (56), we have
× T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k (Gk HkĎ Rk )T = 0, which verified (54). Moreover, using the following relationship: yk − Ck xˆ k = (I − Ck L∗k )˜yk , (55) can be represented in light of the filter structure (3) by choosing: Nk = (I − L∗k Ck )Ak−1 , Ď
(47)
˜ k = Ξk − Ck L∗k , Ξ
(48)
Lk = L∗k ,
˜ k−1 , Kk = (I − L∗k Ck )Gk−1 Hk−1 Ξ where
Using (46), the existence of the gain matrix Lk is given by the following rank condition:
xˆ k = x¯ ∗k + L∗k (yk − Ck x¯ ∗k ),
(46) ∗
rank Hk
(54)
(L∗k R˜ ∗k − P¯k∗ CkT )T¯kT X¯ kT = (G¯ ∗k R˜ ∗k − P¯k∗ CkT ) I − T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k R˜ ∗k
Ď
where U¯ k∗ = Gk (I − Hk Hk ), one has
e˜ k+1 = Aˆ k ek + wk −
(53)
which will be shown later. Finally, using (40)–(41) and (46)–(54), the proposed optimal state estimator, i.e., the ROUMVF, is given by
Lk = Gk + (Pk Ck −
The main idea in deriving a new optimal state estimator which does not need any transformations is to further minimize the trace of the matrix P˜ k+1 of the OUMVF, i.e., (16), and hence the resulting trace of the estimation error covariance (14) would be globally minimized. In deriving the new filter, we first modify (13) and (15), respectively, to give: ∗
Ď
X¯ k = −Gk Hk Rk T¯kT (T¯k R˜ ∗k T¯kT )−1 ,
∗
4.1. A refined OUMVF design
(50)
Ξk = I − Rk T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k .
(63)
Remark 3. As shown in the simulation example in Hsieh (2009), if U¯ k∗ = Gk then rank condition (48) will be equivalent to (21).
712
C.-S. Hsieh / Automatica 46 (2010) 708–715
Then, one has (I − L∗k Ck )Gk−1 = 0, and hence (58) and (59) can be simplified as follows: x¯ ∗k+1
= Ak xˆ k ,
P¯ k∗+1
=
Ak Pk ATk
+ Qk .
Thus, the ROUMVF will reduce to the OEF.
Remark 5. Rewrite (58) and (59), respectively, as follows: x¯ ∗k+1 = Ak xˆ k + Gk Mk∗ (yk − Ck x¯ ∗k ),
(64)
Pk (L∗k R˜ ∗k − P¯k∗ CkT )(Mk∗ )T = Ak Gk • Mk∗ R˜ ∗k (Mk∗ )T T × Ak Gk + Qk ,
(65)
where
¯ ∗k ) I − R˜ ∗k T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k . Mk∗ = Hk (I − Ck G
(66)
¯ ∗k and T¯k are replaced, respectively, by (1) If the above G ¯ ∗k = U¯ k∗−1 (T˘k Ck U¯ k∗−1 )Ď T˘k , G Ď
Note, that if the following relationship holds,
(I − L∗k Ck ) L∗k −P k −Mk∗ Ck Mk∗ −Pkdx P¯ k∗ 0 0 Rk In CkT × 0 HkT −(U¯ k∗−1 )T 0
Sk =
(R˜ ∗k )−1 Sk
Ď
P¯ k∗+1 = Ak
∗ x¯ k
Ď
Mk∗ = Hk Hk
Kk = P¯ k∗ CkT (R˜ ∗k )−1 ,
∗ 0 Sk ,
and, hence, the obtained ROUMVF will be equivalent to the recently developed ERTSF in Hsieh (2009). Furthermore, if rank[Hk ] = m, i.e., Hk is of full column rank, one has
Γk = 0 , L∗k = Kk (I − Hk Mk∗ ), −1 Mk∗ = HkT (R˜ ∗k )−1 Hk HkT (R˜ ∗k )−1 . Thus, the obtained ROUMVF will be equivalent to the RTSF in Gillijns and De Moor (2007b). Note that in this case, the unbiased input estimate dˆ k and its corresponding error variance Pkd can be obtained as follows: dˆ k = Mk∗ (yk − Ck x¯ ∗k ),
Pkd = Mk∗ R˜ ∗k (Mk∗ )T .
Remark 6. Using (55)–(57) and (64)–(66), the state estimate xˆ k can be expressed alternatively as follows:
Gk
=
I Ck
0 Hk 0 0 0
−U¯ k∗−1
Ď
0 0 0 0
In 0
0 0 0, Im 0
Pkxd Pkd
Pk
Pkdx
0 Hk
+
L∗k = Kk + (Γk − Kk Sk )Sk∗ ,
In Ck 0 0 0
0 0
0 Im
0 0
(70)
ATk GTk
+ Qk ,
(71)
then the proposed ROUMVF will be equivalent to the descriptor Kalman filter (DKF) (Nikoukhah, Willsky, & Levy, 1992) that is obtained by using the maximum likelihood estimation to the following system:
(R˜ ∗k )−1 ,
then one has
0 0 In 0 0
−Pkxd 0 = 0 −Pkd
where
yk SkT
0 Ip 0 0 0
0
(67)
R˜ ∗k T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k = I − Sk Sk∗ ,
(69)
= (L∗k R˜ ∗k − P¯k∗ CkT )(Mk∗ )T , Pkd ≡ E [(Φk dk − dˆ k )(Φk dk − dˆ k )T ] = Mk∗ R˜ ∗k (Mk∗ )T .
In
where SkT
(68)
0
Pkxd ≡ E [(xk − xˆ k )(Φk dk − dˆ k )T ]
0 × 0 0
in which the matrix parameter βk must be chosen so that T˘k is of full row rank, then the obtained ROUMVF will be equivalent to the optimal linear minimum-variance estimator (OLMVE) (Hsieh, 2006b, 2009), which is obtained by applying the innovations filtering technique (Hou & Patton, 1998) to the considered system (1)–(2). (2) If the following relationships hold:
− −
Ď
¯ ∗k ), T¯k = T˘k (I − Ck G
T˘k = βk (I − Hk Hk ),
∗
−P k −Pkdx
Lk Mk∗
∗ x¯ k yk 0,
combination of the unknown input Φk dk , in which Φk = Hk Hk (Hsieh, 2009), (Pkdx )T = Pkxd , and
where
Ď Sk = Sk∗ ,
Pkxd Pkd
where dˆ k (≡Mk∗ (yk − Ck x¯ ∗k )) is an unbiased estimate of the linear
Ď
∗
x¯ ∗k+1 = Ak xˆ k + Gk dˆ k ,
Remark 4. Comparing (59) with (16), it is clear that the ROUMVF may have better filtering performance than the OUMVF with Λk = Ď Gk Hk since Q˜ k ≤ 0. In other words, the predicted state estimation error covariance of the former is further minimized as compared to that of the latter. This is the rationale that the proposed ROUMVF is an optimal filter without resorting to the constraint (7) proposed by Darouach et al. (2003). As shown in Remark 3, the ROUMVF will be equivalent to the OUMVF if rank condition (21) holds.
P¯ k∗+1
xˆ k (I − L∗k Ck ) = dˆ k −Mk∗ Ck
Ak−1 0
−U¯ k∗−1
" xk #
0
dk dk−1
+
−wk−1 vk
Gk−1 νk−1 , 0
T
where νk ≡ − (xk − xˆ k )T (Φk dk − dˆ k )T is zero-mean Gaussian and independent of wk and vk+1 . Specifically, if matrix Hk is of full column rank, then one has U¯ k∗−1 = 0 and Φk = I, and, hence,
the obtained dˆ k is the unbiased estimate of the unknown input dk . Moreover, it can be checked that the matrices L∗k , Mk∗ , Pk , Pkxd , and Pkd in (70) are obtained as those given in the RTSF. 4.2. Global optimality evaluation In this subsection, we show that the proposed ROUMVF is globally optimal among all possibly linear filters. To show this, we first replace (55) by xˆ k = x¯ ∗k + L∗k y˜ k +
k−1 X
Mi y˜ i .
(72)
i=0
For the last term in (72) to be unbiased regardless of all past unknown input vectors, i.e., E [Mi y˜ i ] = 0, the gain matrix Mi must be constrained as follows due to (43):
C.-S. Hsieh / Automatica 46 (2010) 708–715
Mi Hi
Ci U¯ i∗−1 = 0,
i = 0, . . . , k − 1.
(73)
713
Proof. Assume that for some mi satisfying (74), condition (75) does not hold, i.e.,
Remark 8. It can be checked that the GOUMVF is, in fact, a subset of the ROUMVF with the specific forms of Hk and Rk in (22), which satisfy relationship (20). Specifically, it is noted that the Ď relationship will always guarantee the condition Hk Rk T¯kT = 0 since the matrix T¯k is in the left null space of Hk . Hence, one has Xk = 0 Ď and L∗k Rk (Hk )T = 0, which renders P¯ k∗+1 in (59) to the simplified one in (39). This shows that the predicted state estimation error covariance of the GOUMVF, i.e., the P˜ k+1 in (39), indeed has a minimized trace. In light of this, constraint (20) for guaranteeing the existence of the optimal filter in the form of (4) has the same role as the matrix Xk to minimize the predicted state estimation error covariance. Thus, by using the same approach in Lemma 1 and Theorem 1, the global optimality of the GOUMVF can be easily verified.
E [(xk − x˜ k )(mi y˜ i )] = b 6= 0,
4.3. A compact form of the ROUMVF
Then, following a similar approach given in Lemma 1 of Kerwin and Prince (2000), we have the following lemma. Lemma 1. Let y˜ 1 , . . . , y˜ k−1 have finite variances and let x˜ k given by (42) be the linear unbiased minimum-variance predicted estimate of xk given by y1 , . . . , yk−1 . If mi is any 1 × p vector satisfying
mi Hi
Ci U¯ i∗−1 = 0,
(74)
where i ≤ k − 1, then E [(xk − x˜ k )(mi y˜ i )] = 0.
(75)
where b is an n × 1 vector. Now consider another unbiased predicted estimate of xk : xˇ k = x˜ k −
mi y˜ i
γ2
b, 2
E [kxk − xˇ k k2 ] = E [kxk − x˜ k k2 ] −
G2,k d¯ k + wk ,
xk+1 = Ak xk + G1,k
where γ = E [(mi y˜ i ) ], then one has 2
As shown in Section 3.2, the ROUMVF can also be implemented more efficiently as below. Let system (1)–(2) be rearranged by renumbering the unknown input vector dk into d¯ k as follows:
bT b
γ2
< E [kxk − x˜ k k ]. 2
This contradicts the fact that x˜ k is obtained by solving the minimization problem: finding Xk that minimizes the trace of P¯ k∗+1 (51). Thus, the lemma is proved. Finally, the desired result is given in the following theorem. Theorem 1. For xˆ k given by (72) with constraints (46) and (73), the mean square error achieves a minimum when M0 = M1 = · · · = Mk−1 = 0. Proof. By application of Lemma 1 to each of the rows of Mi , we obtain E [˜ek (Mi y˜ i )T ] = 0 (i ≤ k − 1).
(76)
xˆ k = x¯ ∗k + L∗k (yk − Ck x¯ ∗k ),
(79)
¯ ∗k + (P¯k∗ CkT − G¯ ∗k R˜ ∗k )T¯kT (T¯k R˜ ∗k T¯kT )−1 T¯k , L∗k = G
(80)
Pk = (I − Lk Ck )P¯ k∗ (I − L∗k Ck )T + L∗k Rk (L∗k )T ,
(81)
˜ kĎ Ξ˜ k (yk G1,k Σ
(82)
∗
x¯ ∗k+1
= Ak xˆ k +
− Ck x¯ k ), ∗
˜ kĎ Ξk Rk (G1,k Σ ˜ kĎ )T P¯ k∗+1 = Aˆ k Pk Aˆ Tk + Qk + G1,k Σ
˜ kĎ )T + G1,k Σ ˜ kĎ Rk (Aˆ k L∗k )T , + Aˆ k L∗k Rk (G1,k Σ
(83)
˜ k and Ξk are given in (63), and where R˜ ∗k is given by (52), Ξ Ď
ek = fL + gM ,
(78)
˜ k is of full column rank. where Σ Applying the algorithm of ROUMVF (55)–(62) to (77)–(78) and using (63), we obtain the following compact ROUMVF (CROUMVF):
˜ k Ck , Aˆ k = Ak − G1,k Σ
Next, using (44), (46) and (72), we have
(77)
0 d¯ k + vk ,
˜k yk = Ck xk + Σ
T¯k = αk (I −
Ď Sk Sk ),
Ď
G2,k−1 Sk ,
¯ ∗k = 0 G
(84)
Ck G2,k−1 .
˜k Sk = Σ
(85)
where fL = (I − L∗k Ck )˜ek − L∗k vk ,
gM = −
k−1 X
Mi y˜ i .
i=0 T Using (76), one can easily verify that E [fL gM ] = 0. Thus, the mean square error achieves a minimum when gM = 0, which occurs for M0 = M1 = · · · = Mk−1 = 0. This completes the proof of the theorem.
Remark 7. As addressed in Kerwin and Prince (2000), the above Lemma 1 does not indicate that the (predicted) estimation error and the data sequence are orthogonal. However, the lemma does lead to the fact that (76) holds for any Mi satisfying (73), which is the key point in showing the global optimality of the ROUMVF. Specifically, using Hk = 0 in (45) yields e˜ k = Ak−1 ek−1 + wk−1 , thereby Lemma 1 will imply the following relationships: E [(xk−1 − xˆ k−1 )(mi y˜ i )] = 0,
E [˜yk (Mi y˜ i )T ] = 0,
which are the main results given in Lemmas 1 and 2 of Kerwin and Prince (2000).
Remark 9. The above CROUMVF has a strong connection to the ˜ k can be further ROF designed by Cheng et al. (2009). If matrix Σ ˜ kT = ΣkT 0 , where Σk is a nonsingular rearranged as Σ diagonal matrix, and the measurement noise covariance Rk can be pre-processed as a diagonal matrix through some output transformations, then one has
˜ kĎ Ξk = Σk−1 Σ
0 ,
Ď
˜ k )T = 0, L∗k Rk (Σ
L∗k = 0
L∗2,k .
Then, the CROUMVF will reduce to the ROF. 5. An illustrative example To illustrate the proposed results, the numerical example given in Cheng et al. (2009) was considered and slightly modified, where the system parameters of (1) and (2) are given as follows: 0.5 0 Ak = 0 0 0
2 0.2 0 0 0
0 1 0.3 0 0
0 0 0 0.7 0
0 1 1 , 1 0.1
1 1 Gk = 0 0 0
0 0 0 0 0
−0.3 0 0 , 0 0
714
C.-S. Hsieh / Automatica 46 (2010) 708–715
Table 1 y
y
Matrices (Tk Hk Tkd ), Tk , and Tkd of the ROUMVF, CROUMVF, GOUMVF, ROF, OUMVF, OEF, and DKF. Filter
(Tky Hk Tkd )
ROUMVF
Hk
CROUMVF
Tkd
I5
I3
1 0 0 0 0
0 0 2 0 0
0 0 0 0 0
I5
0 2 0 0 0
0 0 0 0 0
0 1 0 0 0
1
0 0 0
GOUMVF
y
Tk
0
2
ROF
0 0 0
OUMVF OEF DKF
Hk Hk Hk
0
0 0 0 0 0
1
0 0 0 0 0
1 0 0 0
0 0 1 0 0 0 0 −1 0 0
−0.5
0 1 0 0 0 0 0 0 0
0 1 0
1 0 0
0 1 0
1 0 0
0 0 0 0 1
0 0 1 0
−1
0 0 1
0
−0.5 0 1 0
0
0 1
0 0 0 0 1
0 −1 0
I5 I5 I5
0 0 1
1 0 0
I3 I3 I3
Table 2 Performance of the ROUMVF, CROUMVF, GOUMVF, ROF, OUMVF, OEF, and DKF. Filter
rmse(e1k )
rmse(e2k )
rmse(e3k )
rmse(e4k )
rmse(e5k )
trace(Pk )
ROUMVF CROUMVF GOUMVF ROF OUMVF OEF DKF
1.2783 1.2783 1.2783 1.2783 1.2868 3.4766 1.2783
0.3023 0.3023 0.3023 0.3023 0.3023 0.3023 0.3023
0.0521 0.0521 0.0521 0.0521 0.0521 0.0521 0.0521
0.0604 0.0604 0.0604 0.0604 0.0604 0.0604 0.0604
0.0317 0.0317 0.0317 0.0317 0.0317 0.0317 0.0317
0.1941 0.1941 0.1941 0.1941 0.1947 0.0677 0.1941
0 0 Hk = 0 0 0
0 0 2 0 0
1 0 0 , 0 0
0.5
1
0.3
1
Rk = 10−2 × 0.5
1
0.3
1
1
Qk = 10−4 ×
,
1
1 0.5
0.5 1 1
,
C k = I5 .
1 The unknown input and other simulation settings are the same as those given in the original paper. In the simulation, the ROUMVF (55)–(62), the CROUMVF (79)–(85), the GOUMVF (32)–(39), the ROF in Cheng et al. (2009), the OUMVF (8)–(18), the OEF in Darouach et al. (2003), and the DKF given by (68)–(71) were considered. In outlining a unified procedure to derive a globally optimal filter within the unbiased y minimum-variance filtering framework, we define Tk and Tkd as two transformation matrices associated with the output and the unknown input, respectively. Thus, the original system (1)–(2) can be transformed into the following equivalent system:
(
)=(
y Tk Ck
y Tk Hk Tkd
)xk + (
y Tk
y Tk k
)d¯ k + ( v ),
¯ and is nonsingular. Then, based on the desired where dk = y form of the matrix (Tk Hk Tkd ) of the above equivalent system, one can design a specific globally optimal state estimator. The Tkd dk
rank Hk
Ck Gk−1 = 3
6= rank[Hk ] + rank[Gk−1 ] = 4. The simulation clearly shows the suboptimal essence of the OUMVF and the smaller estimation error covariance of the ROUMVF since trace (Pk (OUMVF)) = 0.1947
> trace (Pk (ROUMVF)) = 0.1941.
xk+1 = Ak xk + (Gk Tkd )d¯ k + wk , y Tk yk
transformation matrices used in deriving the considered filters have been summarized in Table 1, from which it is clear that the designs of the filters: ROUMVF, OUMVF, OEF, and DKF do not need any transformations, while those of the others do. However, it should be stressed that these extra transformation matrices may reduce the computational complexity of the derived filter since the dimensions of some of the parameter matrices have been somehow reduced. We have listed the root-mean-square-errors (rmse) in the state estimates of the aforementioned filters, as well as the traces of their steady-state estimation error covariance, in Table 2, from which we obtained the following results: (1) all the first four filters yield the same filtering performance; (2) the DKF can yield the same filtering performance as that obtained by the first four filters; (3) the filtering performance of the OUMVF is slightly worse than that of the first four filters; this is due to the fact that condition (20) is not satisfied in this simulation case; and (4) the OEF has the worst filtering performance as rank condition (21) is not satisfied because
The simulation result also shows that the OEF is, in fact, more restrictive than the others, due to the fact that rank condition (21) is more stringent than that in (48), which is given as
rank Hk
∗ Ck U¯ k∗−1 = rank[Hk ] + rank[U¯ k−1 ] = 3.
C.-S. Hsieh / Automatica 46 (2010) 708–715
6. Conclusion
715
References
Gillijns, S., & De Moor, B. (2007a). Unbiased minimum-variance input and state estimation for linear discrete-time systems. Automatica, 43, 111–116. Gillijns, S., & De Moor, B. (2007b). Unbiased minimum-variance input and state estimation for linear discrete-time systems with direct feedthrough. Automatica, 43, 934–937. Hou, M., & Patton, R. J. (1998). Optimal filtering for systems with unknown inputs. IEEE Transactions on Automatic Control, 43, 445–449. Hsieh, C.-S. (2000). Robust two-stage Kalman filters for systems with unknown inputs. IEEE Transactions on Automatic Control, 45, 2374–2378. Hsieh, C.-S. (2006a). Optimal minimum-variance filtering for systems with unknown inputs. In Proceedings of the 6th world congress on intelligent control and automation (pp. 1870–1874). Hsieh, C.-S. (2006b). Optimal filtering for systems with unknown inputs via unbiased minimum-variance estimation. In Proceedings of IEEE Tencon 2006. Hsieh, C.-S. (2007). The unified structure of unbiased minimum-variance reducedorder filters. IEEE Transactions on Automatic Control, 52, 1716–1721. Hsieh, C.-S. (2009). Extension of unbiased minimum-variance input and state estimation for systems with unknown inputs. Automatica, 45, 2149–2153. Hsieh, C.-S. (2010). On the optimality of two-stage Kalman filtering for systems with unknown inputs. Asian Journal of Control (in press). Hsieh, C.-S., & Chen, F.-C. (1999). Optimal solution of the two-stage Kalman estimator. IEEE Transactions on Automatic Control, 44, 194–199. Ignani, M. B. (1990). Separate bias Kalman estimator with bias state noise. IEEE Transactions on Automatic Control, 35, 338–341. Kerwin, W. S., & Prince, J. L. (2000). On the optimality of recursive unbiased state estimation with unknown inputs. Automatica, 36, 1381–1383. Kitanidis, P. K. (1987). Unbiased minimum-variance linear state estimation. Automatica, 23, 775–778. Nikoukhah, R., Willsky, A. S., & Levy, B. C. (1992). Kalman filtering and riccati equations for descriptor systems. IEEE Transactions on Automatic Control, 37, 1325–1342. Valcher, M. E. (1999). State observers for discrete-time linear systems with unknown inputs. IEEE Transactions on Automatic Control, 44, 397–401.
Borisov, A. V., & Pankov, A. R. (1994). Optimal filtering in stochastic discretetime systems with unknown inputs. IEEE Transactions on Automatic Control, 39, 2461–2464. Chen, J., & Patton, R. J. (1996). Optimal filtering and robust fault diagnosis of stochastic systems with unknown disturbances. IEE Proceedings Control Theory and Applications, 143, 31–36. Cheng, Y., Ye, H., Wang, Y., & Zhou, D. (2009). Unbiased minimum-variance state estimation for linear systems with unknown input. Automatica, 45, 485–491. Darouach, M., & Zasadzinski, M. (1997). Unbiased minimum variance estimation for systems with unknown exogenous inputs. Automatica, 33, 717–719. Darouach, M., Zasadzinski, M., & Boutayeb, M. (2003). Extension of minimum variance estimation for systems with unknown inputs. Automatica, 39, 867–876. Friedland, B. (1969). Treatment of bias in recursive filtering. IEEE Transactions on Automatic Control, 14, 359–367.
Chien-Shu Hsieh received the B.S., M.S., and Ph.D. degrees in Electrical and Control Engineering from National Chiao-Tung University, Taiwan, in 1982, 1987, and 1999, respectively. From 1987 to 1992, he was an Assistant Researcher at the Chung-Shan Institute of Science and Technology. In 1998, he joined the Mechanical Industry Research Laboratories of Industry Technology Research Institute, Hsinchu, Taiwan, where he was a Researcher. Since 2000, he has been with the Department of Electrical Engineering, Ta Hwa Institute of Technology, Hsinchu, Taiwan, where he is currently an Associate Professor. His research interests include Kalman filtering, state estimation, reliable control, robust control, and intelligent control.
It has been shown that the originally proposed OUMVF can be refined to obtain globally optimal state estimators through further minimizing the predicted state estimation error covariance. To achieve this, two methods are proposed in this paper: one is based on the transformation-based approach to satisfy the decorrelation constraint (7), and the other is the modification of the conventional observer-type filter (4) by using (3). The resulting globally optimal state estimators have been given in three equivalent forms: the GOUMVF, the ROUMVF, and the CROUMVF. Specifically, it has been shown that the existing globally optimal filter ROF is, in fact, a special case of the GOUMVF and the CROUMVF. Through the proposed results, the performance degradation problem of the OEF and the suboptimal essence of the originally developed OUMVF have been clearly illustrated. The relationship with the existing literature results has also been addressed, and indicated that the optimal filters developed by Gillijns and De Moor (2007b), Hou and Patton (1998), Hsieh (2009), and Nikoukhah et al. (1992) are all globally optimal in addition to that proposed by Cheng et al. (2009). An equivalent system description obtained through output and input transformations serves as a unified approach to design a specific globally optimal state estimator for systems with unknown inputs. Simulation results verify the proposed results.