Multi-rate stochastic H∞ filtering for networked multi-sensor fusion

Automatica 46 (2010) 437–444

Contents lists available at ScienceDirect

Automatica journal homepage: www.elsevier.com/locate/automatica

Brief paper

Multi-rate stochastic H∞ filtering for networked multi-sensor fusionI Yan Liang a , Tongwen Chen b,∗ , Quan Pan a a

School of Automation, Northwestern Polytechnical University, Xi’an, 710072, China

b

Department of Electrical and Computer Engineering, University of Alberta, Edmonton, Alberta, Canada T6G2V4

article

info

Article history: Received 9 January 2009 Received in revised form 28 October 2009 Accepted 12 November 2009 Available online 5 December 2009 Keywords: Multi-rate systems Multi-sensor fusion Networked control systems H∞ filtering Unknown input observers

abstract This paper presents a multi-rate filtering problem for a class of networked multi-sensor fusion systems with packet dropouts (PDs): the state evolves according to a linear discrete-time model with normbounded unknown inputs (UIs), and its underlying period is h; the p sensors are distributively deployed with different sampling periods n1 h, . . . , np h; multi-rate sensor measurements, corrupted by UIs, are subject to stochastic PDs in the transmission to a fusion center for state estimation; the estimation is updated at the period mh. Different from the single-rate estimator design with PDs which are treated as stochastic parameters, a UI observer is proposed where PDs are represented as zero-mean white input noises of the linear time-variant estimation error system. The results on the existence of a stable observer are proposed. Due to insufficient design freedom for absolute error decoupling, we turn to designing an observer-based stochastic H∞ filter. A numerical example of distributive multi-sensor target tracking is given to illustrate the proposed filter. © 2009 Elsevier Ltd. All rights reserved.

1. Introduction Networked control systems (NCSs) have gained attention during the last few years (e.g., see Basin & Martinez-Zuniga, 2004; Hespanha, Naghshtabrizi, & Xu, 2007; Hu & Zhu, 2003; Luck & Ray, 1990; Matveev & Savkin, 2003; Nilsson, Bernhardsson, & Wittenmark, 1998; Wang, Yang, Ho, & Liu, 2006; Zhang, Basin, & Skliar, 2007, and references therein). When compared with using the conventional point-to-point system connection, using an NCS has advantages like easy installation and reduced setup, wiring, and maintenance costs. In an NCS, data travel through the communication channels from sensors to the controller and/or from the controller to the actuators. For the estimation problem, an immediate concern arising from the introduction of these network media is that the network properties, such as random transmission delay and packet dropout (PD), give rise to parameter uncertainties. Most research conducted on NCSs considers the issues of random transmission delay and uncertain observations where the associated systems can be described by stochastic systems (see e.g., Basin & Martinez-Zuniga, 2004; Nilsson et al., 1998; Wang

I A part of the work in this paper will be presented at the joint 48th IEEE Conference on Decision and Control and 28th Chinese Control Conference, December 16–18, 2009, Shanghai, China. This paper was recommended for publication in revised form by Associate Editor Hideaki Ishii under the direction of Editor Ian R. Petersen. ∗ Corresponding author. Tel.: +1 780 492 3940; fax: +1 780 492 1811. E-mail addresses: [email protected] (Y. Liang), [email protected] (T. Chen), [email protected] (Q. Pan).

0005-1098/$ – see front matter © 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.automatica.2009.11.019

et al., 2006; Zhang et al., 2007, and references therein). The closely related random PD has been the focus of some research studies in the last few years. The stability of a time-varying Kalman filter in relation to a PD rate was investigated (Sinopoli, Schenato, Franceschetti, Jordan, & Sastry, 2004). Through transforming the PD rate into a stochastic parameter in the system’s representation, stochastic H2 and H∞ norms of an estimation error system were defined, and thus H2 and H∞ filters were proposed to deal with the possible delay of one sampling period, uncertain observations and multiple PDs under a unified framework (Sahebsara, Chen, & Shah, 2007, 2008). Considering a linear system with stochastic parameters due to PDs, Sun, Xie, Xiao, and Soh (2008a) proposed the optimal linear minimum-variance (LMV) estimators, including filters, predictors and smoothers via innovation analysis approaches. Furthermore, a reduced-order filter was designed (Sun, Xie, Xiao, & Soh, 2008b) and the maximum successive PD number was considered (Sun, Xie, Xiao, & Soh, 2008c). In Gao and Chen (2007), a robust H∞ estimation scheme for a class of linear systems with polytypic uncertain parameters was studied subject to limited communication capacity; here, as the main properties of limited communication channels, measurement quantization, data transmission delay, and PDs were considered simultaneously in the estimator design. One common step of the above methods is to augment the state using the transmitted measurement and/or received input to capture the dynamic properties introduced by PDs, and then design filters for linear systems with stochastic parameters. In many complex systems, it is often unrealistic or sometimes impossible to guarantee all physical signals operating at one single rate (Chen & Francis, 1995). State estimation with dual-rate

438

Y. Liang et al. / Automatica 46 (2010) 437–444

sensors was proposed in Andiusani and Gau (1987). Its dynamic system was decomposed into dual subsystems corresponding to the dual-rate sensors; the filtering residuals of one Kalman filter for the fast-rate subsystem were fused with slow-rate sensor measurements by another Kalman filter based on the slow-rate subsystem. Since different dynamic models always have different frequency properties in multiple model systems, the fast-rate sensor measurement was thus compressed to be a slow-rate substitution with little or no accuracy degradation in low-frequency models, and cost-effective multi-rate interacting multiple model estimators were proposed (Hong, 1999). In the case that the updating rate of state estimates is a multiple of the measurement sampling rate, fast-rate estimation was proposed by extending traditional single-rate H2 or H∞ filters (Sheng, Chen, & Shah, 2005). In fast-rate fault detection for sampled-data systems with slow-rate measurements, fast-rate residuals were generated (Izadi, Zhao, & Chen, 2005; Zhong, Ye, Ding, & Wang, 2007). Recently, we proposed the LMV estimator for a four-rate estimation problem, which includes the state updating rate in the model, the measurement sampling rate, the estimate updating rate, and the estimate output rate (Liang, Chen, & Pan, 2009). In general, all the above multi-rate filters cannot deal with stochastic time-varying parameters. Thus, designing a multi-rate estimator under stochastic PDs is not trivial by the fact that state augmentation, as the common step in the existing single-rate estimators with PDs, inevitably leads to stochastic parameters. Here, we present the multi-rate filtering problem for multisensor fusion with PDs: sensors are distributively deployed and their sampling rates are not identical; unknown inputs (UIs) drive the evolvement of the interested state and corrupt sensor measurements; multi-rate sensor measurements are subject to PDs in the transmission to the fusion center for state estimation. To the best of the authors’ knowledge, this is the first consideration of multirate filtering with PDs. Instead along the idea of filter design with stochastic parameters, the problem of multi-rate sensor fusion is transformed into the design of an unknown input observer (UIO). It is interesting that the resultant estimation error system of the observer is a linear time-invariant system with deterministic parameters, and stochastic PDs are presented as the zero-mean white input noises. The conditions for the existence of a stable UIO are proposed. It is founded that such conditions are dependent on the system matrix, measurement matrix, sampling rate of each sensor, and estimation rate, but independent of the non-zero probability of successful packet transmission. Based on the resultant UIO, we design a multi-rate stochastic H∞ estimator. A numerical example of distributive multi-sensor target tracking is illustrated. It is worth mentioning that the multirate nature gives rise to higher computation complexity in two aspects: first, the lifting process, a common step in treating multirate systems, leads to higher dimensional parameter matrices to be optimized; second, the resultant estimator design is non-convex, and thus requires a computation-intensive iterative parameter optimization. Fortunately, such computational issues would not affect the estimator implementation because all estimator parameters are time-invariant and can be optimized off-line. The rest of this paper is organized as follows. The problem under investigation is formulated in Section 2. The UIO is designed and analyzed in Section 3. Based on the resultant UIO, the multi-rate stochastic H∞ filter is developed in Section 4. A simulation example is given in Section 5 and some conclusions are drawn in Section 6. Throughout this paper, the superscripts ‘‘−1’’ and ‘‘T ’’ represent the inverse and transpose operations, respectively; the symbols ‘‘I’’ and ‘‘0’’, respectively, represent identity and zero matrices with proper dimensions; diag{·} denotes a block diagonal matrix; prob{·} denotes a probability measure; E {·} denotes the operator of mathematical expectation; trace{·} denotes the operator of matrix trace; ↓ m represents down-sampler by m; q−1 is the discrete-time unit delay operator; rank{·} represents the matrix rank.

Fig. 1. Multi-rate filtering problem.

2. Problem formulation As shown in Fig. 1, we present the following linear multi-rate multi-sensor system: xk+1 = Axk + Bwk , yi,k = Ci xni k + Di wni k ,

(1)

(i = 1, 2, . . . , p),

zmk = Lxmk ,

(2) (3)

where x is the state vector evolving with a period h; h is a positive real number; yi is the ith channel sensor measurement vector with a sampling period ni h; z is the vector to be estimated with an updating period mh; both the sampling rate ni (i = 1, 2, . . . , p) and the estimation rate m are positive integers; w is a normbounded UI vector; the matrices A, B, Ci , Di and L are known with proper dimensions. In (1)–(3), there are p + 2 periods: the state updating period h, the p sensor sampling periods n1 h, . . . , np h, and the estimation output period mh. Moreover, the sensors and a fusion unit are considered to be spatially separated and connected via network media pulsed buffers. Each buffer records only the latest received measurement in the corresponding channel. Due to PDs, the buffered measurement could be different from the original measurement as follows: y∗i,k = i,k yi,k + (1 − i,k )y∗i,k−1

(4)

where the stochastic parameter i,k takes the value of 0 or 1, representing PD or not. Similar to Sahebsara et al. (2007) and Sun et al. (2008c), i,k is unknown in the estimator design and satisfies the following assumption. Assumption 1. The stochastic parameters i,k ’s are Bernoulli distributed white sequences satisfying prob i,k = 1 = αi ,





0 < αi ≤ 1,

(5)

where αi is the probability of successful packet transmission for the ith channel sensor measurement. The main requirements of multi-rate estimator design for (1)– (5) under Assumption 1 in the fusion unit are:

• the estimator parameters can be calculated off-line, which is desirable in real-time applications;

Y. Liang et al. / Automatica 46 (2010) 437–444

• the estimator should be robust to two kinds of uncertainties: one is w from the original system and sensors; the other is  from data transmission; • the estimator is causal, i.e., no received measurements later than kmh are utilized in estimating zmk . Remark 1. Now, a multi-rate filtering problem for networked multi-sensor fusion systems with PDs is proposed. To the best of our knowledge, this is the first attempt at multi-rate estimation with PDs. As mentioned in the introduction, the proposed problem is non-trivial since the common idea in single-rate estimation is to transform PDs into stochastic parameters, while the existing multirate estimators cannot deal with stochastic parameters. As shown later, we turn to designing one UIO with causality constraints so that PDs are represented as zero-mean white input noises of the linear time-invariant estimation error system; then, we define a stochastic H∞ norm and further design a UIO-based robust filter. 3. Unknown input observer Define the least common multiple of m and ni , i = 1, 2, . . . , p by N = LCM{m, ni , i = 1, 2, . . . , p}.

(6)

Using the lifting technique (Chen & Francis, 1995) on (1), (2) and (4), we have the single-rate lifted system with a slow period Nh: xk+1 = Axk + B w k , y

i,k

= C i xk + Di wk ,

y∗ = αi y i,k

i,k

i = 1, 2, . . . , p

B = [AN −1 B

Di

 Si,ni −1 .. .

D∗i =  

Si,N −ni −1

Si,ni −2

··· ···

Si,N −ni −2

··· ···

0

.. .

, o

p,k

y∗d = col y∗d , y∗d , . . . , y∗d 1,k

k

2 ,k

p,k

, o

n

y∗k = col y∗1,m1,k , y∗2,m2,k , . . . , y∗p,mp,k ,

n

o

y∗k d = col y∗1,m1,k −1 , y∗2,m2,k −1 , . . . , y∗p,mp,k −1 , where E1 , . . . , Eq , F1 , . . . , Fq+1 , F1∗ , . . . , Fq∗+1 , G1 , . . . , Gq , H1 , . . . , Hq+1 , H1∗ , . . . , Hq∗+1 are the periodically time-varying parameter matrices to be determined. Denoting

e xk = xk −b xk ,

(14)

∆j,k = (i,kN +N − αi )(yj,mj,k −

y∗j,mj,k −1 ),

(15)

e zk = zk −b zk ,

(16)

we transform the UIO in (10)–(13) into the following expression through decoupling the estimation error with the state and its measurements (see the Appendix):

b xkN +mi = Amib xkN + Fi 1y∗ , k

i = 1, . . . , q − 1,

(17)

(19)

AB

B],

∗

∗

i = 1, . . . , q

0 ,



···

1y∗j,k = y∗j,k − (1 − αj )y∗j,dk − αj C jb xkN , 1y∗j,mj,k = y∗j,mj,k − (1 − αj )y∗j,mj,k −1 − αj Cj Ab xkN .

Si,0

..

. ···

Di

..

. ···

··· 0

..

Hi,j = diag

Si,0



..  . .



(10)

b xkN +N = Eqb xkN + Fq y∗ + Fq∗ y∗d + Fq+1 y∗k + Fq∗+1 y∗k d , k k

(11)

b zkN +mi = Gib xkN + Hi y∗ + Hi∗ y∗d , k k

(12)

i ∈ {1, 2, . . . , q − 1},

b zkN +N = Gqb xkN + Hq y∗ + Hq∗ y∗d + Hq+1 y∗k + Hq∗+1 y∗k d , k k

li,j blocks

0 · · 0} | ·{z

 (20)

N /nj −li,j blocks

Fi = Ki,1 Hi,1 · · · Ki,p Hi,p ,

Di

i ∈ {1, 2, . . . , q − 1},

I| ·{z · · }I

with li,j = bmi/nj c + 1, where bmi/nj c represents the largest integer not larger than mi/nj . As shown in the following remark, Fi can be represented by

0 0

.

As shown in (8), the measurements at different time instants are lifted to be a new measurement vector. It gives rise to the causality constraint that the future measurements should not be utilized in the current estimate in a physically realizable observer. In other words, some elements (weights of future measurements) of Fi in (17) should be zero. Denote the causality constraint matrix



Si,j = Ci Aj B

···

n o 1y∗k = col 1y∗1,k , 1y∗2,k , . . . , 1y∗p,k , n o 1y∗k = col 1y∗1,m1,k , 1y∗2,m2,k , . . . , 1y∗p,mp,k , with

The lifted system in (7) and (8) evolves at the slow period Nh, while the estimate should be updated at the fast period mh. At every time (kN + mi)h, 0 < i ≤ q, the observer updates the state estimate by using the available measurements in the kth frame period (kNh, (k + 1)Nh]. As shown in (4), the received measurement y∗i,k has the dynamic property due to stochastic PDs. Hence, y∗i,k and y∗i,k−1 are expected to appear simultaneously in the following UIO to compensate such dynamics:

b xkN +mi = Eib xkN + Fi y∗ + Fi∗ y∗d , k k

2 ,k

n

b zkN +mi = Lb xkN +mi ,

C i A N −n i



1,k

k

o

(8)

i



y∗ = col y∗ , y∗ , . . . , y∗

(18)

e i,k = diag{(i,kN − αi )I , (i,kN +ni − αi )I , . . . , (i,kN +N −ni − αi )I }  C  ∗ Di = Di

mj,k = (k + 1)N /nj ,

n

b xkN +N = Ab xkN + Fq 1y + Fq+1 1yk , k

(9)

··· wk = col {wkN , wkN +1 , . . . , wkN +N −1 } ,  y = col yi,kN , yi,kN +ni , . . . , yi,kN +N −ni , i,k  y∗ = col y∗i,kN , y∗i,kN +ni , . . . , y∗i,kN +N −ni , i,k  y∗d = col y∗i,kN −ni , y∗i,kN , . . . , y∗i,kN +N −2ni , i,k   ∆i,k = e i,k yi,k − y∗i,dk ,

n  Ci A i   Ci =  , ..  . 

q = N /m,

(7)

where A = AN ,

with

where

+ (1 − αi )y∗i,dk + ∆i,k

xk = xkN ,

439

(13)



(21)

where Ki,j , i ∈ {1, 2, . . . , q − 1}, j ∈ {1, 2, . . . , p} is a free matrix to be designed. Remark 2. Using (21), we have Fi 1y∗ =

p X

k

Ki,j (Hi,j 1y∗ )

j =1

j ,k

with Hi,j 1y∗ = col{y∗j,kN − (1 − αj )y∗j,kN −nj , . . . , y∗j,kN + j ,k

− (1 −

bmi/nj cnj

αj )y∗j,kN +bmi/n −1cn , 0, . . . , 0} j j

− αj Hi,j C jb xkN .

(22)

As shown in (22), no future measurements are utilized in estimating the state at every time instant (kN + mi)h.

440

Y. Liang et al. / Automatica 46 (2010) 437–444

Constructing a matrix F with its ith row sub-block being Fi in (21), and putting (A.3) and (A.5) into (A.1), (A.4)–(A.6) into (A.2), (A.9) into (A.7) and (A.10) into (A.8), we obtain the estimation error system of the observer in (17)–(19)

 T :



H = H1

ξk+1 = Aξk + B1 wk + B2 ∆k ek = C ξk + D1 wk∗ + D2 ∆∗k , ∗

∗

(23)

B1 = B − Fq D − Fq+1 C B



T C = α1 C T1 α2 C T2 · · · αp C Tp ,  T C 2 = α1 C1T α2 C2T · · · αp CpT ,  T D1 = α1 DT1 α2 DT2 · · · αp DTp ,  T D2 = α1 DT1 α2 DT2 · · · αp DTp , h
T
T iT e C 0 = LT , LAm · · · LAN −m  0 ··· · · · LB 0 · · ·  LAm−1 B  .. .. .. .. ..  e D10 =  . . . . .  (q−1)m−1 LA B ··· LB | {z } N −m blocks

e D20 =  

··· ..

0 ,



C∗



0 0

.. .. ..   . . . ,  0 · · · 0 | {z } m blocks

  . 

L

|

{z

q−1 blocks

,



(24)

Thus, there exist such gain matrices Fq and Fq+1 that all eigenvalues of the matrix A are asymptotically stable if and only if pair (C ∗ , AN ) is detectable. Re-scaling (C ∗ )T , and further rearranging its column blocks, we obtain

0

.

T

Fq+1 C ∗ .

A = AN − Fq



2



0 L

(C 2 AN )T



with



 (Ci AN )T .

we rewrite the matrix A in (23) as

−Fq+1 D ,   B2 = −Fq −Fq+1 , C =e C0 − e D20 F C 1 ,    D1 = e D2 = −e D10 − e D20 F D1 0 , D20 F 2

···

C ∗ = (C 1 )T

A = A − Fq C 1 − Fq+1 C 2 A, 1

(Ci Ani )T

,

Proof. Denoting

∆k = col{∆1,k , ∆2,k , . . . , ∆p,k },



Hp

with



ξk = e xkN , ek = col{e zkN ,e zkN +m , . . . ,e zkN +N −m }, ∗ ∗ wk = col{w k , wkN +N }, ∆k = col{∆k , ∆k }, ∆k = col{∆1,k , ∆2,k , . . . , ∆p,k },

T

···

H2

H i = CiT

where

1

Theorem 1. There exist such gain matrices Fq and Fq+1 that the estimation error system (23) of the observer (17)–(19) for (1)–(5) is stable, if and only if pair (H , AN ) is detectable, where

}

It is clear that ∆∗k is a zero-mean white noise from (5) and Assumption 1. We further have that ∆∗k is independent of ξk and wk∗ . Remark 3. Different from the existing single-rate estimation schemes where PDs are directly treated as stochastic parameters, the proposed UIO represents PDs as the additional zero-mean white noise inputs ∆∗k to the linear time-invariant error system in (23). It will result in different filter analysis and design. Remark 4. It is easily verified that it is impossible to design B1 , B2 , D1 , and D2 as zero matrices with proper dimensions through choosing Fq,j , Fq+1,j , and Ki,j with i = 1, . . . , q − 1 and j = 1, . . . , p. In other words, there is no sufficient design freedom to decouple the estimation error e in (23) with wk∗ and ∆∗k so that the estimation error converges to zero as k → ∞. Thus, the following design cannot follow the idea of classical observer design that first absolutely decouples the estimation error with UIs and then implements the pole placement (Valcher, 1999). In the following part, we will first explore the condition that the estimation error system (23) is asymptotically stable, and further optimize the estimator parameters to minimize the gain from the UIs to the estimation error.

h
T = α1 C 1


T


T · · · αp C p
T
T i · · · α1 C1 AN · · · αp Cp AN h

T
T T ∼ C1 · · · Cp · · · C1 AN ··· h



T i T T T ∼ C1 C1 A N · · · Cp Cp AN   ∼ H1 · · · Hp .

Cp A N


T i

It is known that rearranging or re-scaling the measurement vector does not change the detectability, namely, the detectability of pair (C ∗ , AN ) and that of pair (H , AN ) are equivalent. Thus, we complete the proof.  Proposition 1. There does not exist the stable estimation error system (23) of the observer (17)–(19) for (1)–(5), if and only if there exists at least one common unobservable and unstable mode among the pairs (C1 , AN ), (C1 An1 , AN ), . . . , (C1 AN , AN ), . . . , (Cp , AN ), (Cp Anp , AN ), . . . , (Cp AN , AN ). Based on a proof by contradiction, this proposition is easily obtained from Theorem 1.  Theorem 2. There exist such gain matrices Fq and Fq+1 that all eigenvalues of the estimation error system (23) of the observer in (17)–(19) for (1)–(5) can be arbitrarily assigned if and only if rank{M } = dx , where dx is the dimension of x and



M = M1

M2

···

Mp



with Mi = CiT



(Ci Ani )T

···

 (Ci (Ani )dx −1 )T .

Proof. First, we examine the rank of C ∗ as follows rank C ∗





= rank

n

= rank



= rank



···

(Cp AN )T  (Cp AN )T

···

C T1

(C1 AN )T

···

C Tp

C1T

(C1 An1 )T

···

(C1 AN )T  (Cp AN )T ,

C Tp

C Tp

(C1 AN )T

C T1

(Cp Anp )T

···

···

T o

(25)

where the first equality uses the fact that αj > 0, j = 1, . . . , p; the second equality is derived through rearranging the columns; the last equality is obtained through substituting C i in (8). Then, we further check the rank of the observability Gramian matrix

Y. Liang et al. / Automatica 46 (2010) 437–444

rank



(C ∗ )T 

(C1 A )



···

(Cp A )

np T

C Tp

···

···

n1 T

C1T

= rank

= rank

(C ∗ A)T

M1

(C1 A

) 

Ndx T

(Cp A  Mp ,

)

Ndx T

···

···

M2

 (C ∗ AN (dx −1) )T

(26)

where the first equality uses (25) and the second equality uses the Cayley–Hamilton theorem so that the dx × dx matrix (Anj )i , for arbitrary i ≥ dx , is a linear combination of I , Anj , . . . and (Anj )dx −1 , and thus the redundant column blocks are deleted without changing the rank. From (24), it is concluded that there exist such gain matrices Fq and Fq+1 that all eigenvalues of the estimation error system in (23) can be arbitrarily assigned, if and only if pair (C ∗ , AN ) is observable, i.e. rank



∗ N (dx −1) T

(C ) (C A ) · · · (C A ∗ T

∗ N T

)



= dx .

Combining it with (26), we complete the proof.



Remark 5. Theorem 1 and Proposition 1 are established based on detectability, which is a weaker condition than observability. Based on the observability, Theorem 2 supplies a sufficient condition for the existence of the stable estimation error system in (23). As shown in Theorems 1 and 2 and Proposition 1, it is concluded:

• the existence of the stable estimation error system (23) depends on system matrix A in (1), measurement matrix Ci in (2), and sampling rate ni in (2); • the existence of the stable estimation error system (23) is independent of the causality constraint and PDs (if the probability of successful packet transmission in each channel αi 6= 0); • the least common multiple N in (6) appears in the sufficient and necessary conditions proposed by Theorem 1, but does not appear in the sufficient condition proposed by Theorem 2. • if N in (6) equals LCM{ni , i = 1, 2, . . . , p}, then m in (3) will have nothing to do with the stability conditions proposed by Theorem 1 or 2. Remark 6. In Theorems 1 and 2 and Proposition 1, the existence of a stable observer is independent of PDs. In the following part, we will show that in some cases, severe PDs can render the UIO design impractical in the case that A is not stable. If A is not stable, then AN will not be stable, and pair (0, AN ) will not be detectable. Denoting



K = K1

···

K2

Kp ,

Ki = αi Fq,i





Fq+1,i ,



we represent A in (23) as A − K H, where H is defined in Theorem 1. If the pair (H , AN ) is detectable, then we can find the suitable gain K for the stable UIO. Denote the set k = {1, 2, . . . , p}. Its subset k1 = {l1 , l2 , . . . , lf } is defined as one critical channel set if the pair ∗

where



H = H l1 H

∗∗

absolutely depends on the critical measurement channel sets, and thus each component of K related to the critical channel set should not be zero. Otherwise, it means that the feasible K can be designed based on a proper subset of k1 , which conflicts with the definition of k1 . Thus, for any measurement channel i ∈ k1 , we obtain its filter gain:



Fq,i

Fq+1,i = Ki /αi 6= 0.



(27)

As shown in (27), the filter gain related to channel i ∈ k1 is inversely proportional to αi in the case that absolute PDs happen in the trivial channel set. In such a case, as αi → 0 (the probability of successful packet transmission in any channel i in the critical set k1 approaches zero), its non-zero filter gain in (27) will become impractical in the implementation. 4. Multi-rate stochastic H∞ filter design Before proceeding to the main theorem, we consider the following definitions adapted from Sahebsara et al. (2008): Definition 1. The estimation error system T as defined in (23) is exponentially stable in the mean-square sense or internally stable if there exists β > 0, and 0 < τ < 1 such that with zero inputs (∆∗k = 0 and wk∗ = 0),

E {kξk k22 } < βτ k E {kξ0 k22 },

∀k > 0.

(28)

In the deterministic parameter case, the H∞ norm of a linear discrete time-invariant system is the maximum bound on the 2-norm of the output over all inputs of the unit 2-norm. For a compromise in combining the deterministic and stochastic inputs in (23), we introduce a known weighting factor ρ > 0, and define a weighted stochastic H∞ norm of T in (23) by kT k∞ as follows: ∞ P

kT k2∞ = sup

E kek k22





k=0

n

2 o . 2 E ∆∗ + ρ w ∗

∞ P

k

k=0

2

k

(29)

2

Remark 7. The weighting factor ρ is utilized to trade off the effect of deterministic and stochastic inputs on the estimation error. If ρ is a very large number, then the gain from wk∗ to ek will be attenuated with priority. Otherwise, minimizing the gain from ∆∗k to ek will be prioritized. In the case of ρ = 1, minimizing the proposed H∞ norm is optimal in the minimum mean-square error sense. Definition 2. For ρ > 0, the estimation error system T in (23) achieves an H∞ performance level γ > 0, i.e., ∞ ∞ n X X

2 o  2 E kek k22 ≤ γ 2 E ∆∗k 2 + ρ wk∗ 2 . k=0

(30)

k=0

∗∗

(H , AN ) is detectable while any pair (H , AN ) is not detectable, ∗

441



= H t1

H l2 H t2

··· ···

H lf

T

H ts

, T

,

and {t1 , t2 , . . . , ts } is a proper subset of k1 . If the pair (H , AN ) is detectable and A is not stable, there must exist a critical channel set from the fact that (H , AN ) is detectable and (0, AN ) is not detectable. It is worth mentioning that the critical channel set may not be unique. Denote the trivial channel set k2 = k −k1 . Consider the case that all probabilities of successful packet transmission related to the trivial channel set approach zero, i.e., their packets are absolutely lost in such channels. In such a case, the observer design

Theorem 3. For γ > 0, the estimation error system T in (23) is exponentially stable in the mean-square sense and achieves an H∞ performance level γ if there exists a positive-definite matrix P such that

 Φ11 Φ1 = ∗

 Φ12 < 0, Φ22

T

T

Φ2 = B2 PB2 + D2 D2 − γ 2 I < 0 with T

T

Φ11 = A PA + C C − P , T

T

Φ12 = A PB1 + C D1 , Φ22 =

T B1 PB1

+

T D1 D1

− ργ 2 I .

(31) (32)

442

Y. Liang et al. / Automatica 46 (2010) 437–444

Proof. First, we verify that the estimation error system T in (23) is T

exponentially stable in the mean-square sense. As C C ≥ 0, from (31), we have T

A PA − P < 0,

(33)

and thus there exists a small enough positive number 1 > ε0 > 0 satisfying T

A PA − P < −ε0 P ,

(34)

by (33) and P > 0. In the case of zero inputs, for any non-zero vector ξk , we have

n



T



E ξkT P ξk = E ξkT−1 P ξk−1 + E ξkT−1 A PA − P ξk−1









(35)

where the first inequality uses (34), and the second inequality is obtained through iteratively using the first inequality. Since P > 0, there exist ε1 > 0 and ε2 > 0 such that ε1 I ≤ P ≤ ε2 I, namely,

  ε1 E kξi k22 ≤ E ξiT P ξi ≤ ε2 E {kξi k22 },

∀i = 0, 1, . . . .

(36)

Combining (35) and (36), we have (28) with β = ε2 /ε1 > 0 and 0 < τ = 1 − ε0 < 1. Now, we verify that the estimation error system T as defined in (23) achieves kT k∞ < γ . Using the positive-definite matrix P, we define a Lyapunov function as Vk = E ξkT P ξk ,





(37)

and thus have

1Vk = Vk+1 − Vk n 

2 o 2 = E − kek k22 + γ 2 ∆∗k 2 + ρ wk∗ 2 h   
T
T i ξ + E ξkT Φ1 k∗ + ∆∗k Φ2 ∆∗k , wk∗ wk



n



T

T



o

(38)

E ξkT

+

T D1 D2



o

∆∗k = 0.

Putting (31) and (32) into (38), we have

n



o

2 2 E kek k22 < E γ 2 ∆∗k 2 + ρ wk∗ 2







− 1Vk .

(39)

Summing up both sides of (39) for k = 0, 1, . . . , ∞, we get ∞ ∞ n X X

2 o  2 E kek k22 ≤ γ 2 E ∆∗k 2 + ρ wk∗ 2 + V0 − V∞ , (40) k =0

k=0

and thus obtain (30) under zero initial conditions.



with

∀i = 1, 2,



B Ψ2 = 2 , D2

Γ1 = diag(−P , −ργ 2 I ),

Π = diag(−Q , −I ),

Γ2 = −γ 2 I ,

PQ = I .

Using the Schur complement, we have the following equivalent expression of (41):



Π ∗

 Ψi < 0, Γi

i = 1, 2.

(42)

Here, the above H∞ filter design problem is equivalent to the following optimization problem: min

γ subject to (42) and PQ = I .

(43)

It is worth mentioning that the optimization in (43) includes a non-convex constraint, that is, PQ = I. Even in the case of no PD, namely, αi = 1, this non-convex constraint would still exist. Thus, it is concluded that such a non-convex constraint comes from the multi-rate nature. The optimization in (43) may be tackled by the product reduction algorithm in Oliveira and Geromel (1997) as follows: given a large enough γ , choose initial values P (0) and Q (0) satisfying (42), and iteratively implement the following convex optimization: min

Fq ,Fq+1 ,Ki,j ,P =P (k+1) ,Q =Q (k+1)

 subject to (42) and

P I

trace(P (k) Q + Q (k) P )

I Q

 ≥0

(44)

until convergence; then, try a smaller γ and repeat the above iterative optimization until the smallest γ is obtained.

Remark 9. In designing the multi-rate estimator with PDs, the multi-rate nature gives rise to design complexity in the following ways. causality constraint matrices, the causality is guaranteed in the estimator design, not affecting the LMI-based optimization. • The second is the non-convex optimization, which has to be implemented via a computation-intensive iterative optimization in the estimator design. • The third is high dimensions in parameter matrices. As shown in (7)–(9) and (23), lifting results in larger-sized parameter matrices, such as Bi , C , and Di . These drawbacks impact only the parameter optimization process (which is off-line), not the estimator implementation, since the estimator parameters are time-invariant. 5. Numerical example

To minimize the effect of deterministic and stochastic inputs on the estimation error, we state the following filtering problem: Multi-rate stochastic H∞ filtering problem: Design a filter as in (17)–(19) satisfying (31) and (32) with the minimum value γ such that the filtering error system in (23) is exponentially stable in the mean-square sense and the H∞ condition in (30) is satisfied. Notice in (23) that A, B1 , B2 , C , D1 and D2 , are all affine in the design parameter matrices (Fq , Fq+1 , and Ki,j ), and that the inequalities in (31) and (32) are not convex. By using matrix manipulations, we represent (31) and (32) by

− ΨiT Π −1 Ψi + Γi < 0,



• The first is the causality constraint. Through introducing the

E ξkT A PB1 + C D1 ∆∗k = 0, T B1 PB2



B1 , D1

Remark 8. The proposed estimator is implemented as follows: check the stability conditions proposed by Theorem 1; if the conditions are satisfied, iteratively optimize the filtering parameters by using (44); obtain the recursive estimate from the UIO in (11)–(13).

by using the fact that

n

A Ψ1 = C

Fq ,Fq+1 ,Ki,j ,P ,Q

o

 < (1 − ε0 )E ξkT−1 P ξk−1  < (1 − ε0 )k E ξ0T P ξ0 ,



(41)

Consider an example of distributive target tracking with two sensors. According to Newton’s force principle, the movement of a maneuvering target can be modeled as a constant-velocity nominal system with a bounded but unknown acceleration input. One sensor (called S1) samples the target position (range) with the period 2h, while the other (called S2) samples the target velocity (Doppler) with the period 3h. The measurements of both sensors are corrupted by bounded errors, and transmitted to the centralized processing unit (CPU) for the state estimation in both position and velocity at a fast rate with the period h. Stochastic PDs are independently Bernoulli distributed with the parameters being α1 and α2 in the S1-to-CPU and S2-to-CPU channels, respectively.

Position and its estimates

Y. Liang et al. / Automatica 46 (2010) 437–444 150 100 12 actual value method 1 method 2

50 0

0

20

40

60

80

10

100

8

120

6 γ

Time instant (0.5s) Velocity and its estimates

443

15

actual value method 1 method 2

10

4 2

5 0

0 0.1

-5

0.2

-10 0

20

40

60

80

100

0.3

0.4

0.5

120

0.6

Time instant (0.5s)

0.7

0.8

α2

Fig. 2. Method comparison with α1 = 0.7, α2 = 0.7, and ρ = 0.75.

This situation can be modeled by (1)–(5) with the following parameters: h , 1

h2 /2 B= h

C2 = 0

1 ,

D1 = 0

L = I,

p = 2,



1 A= 0













0 0 1

n1 = 2,

0 , 0 0 ,



0 ,





C1 = 1



D2 = 0

n2 = 3,

0

1 ,



m = 1.

Here, the estimate updating period h is set to be 0.5 s. It is easily verified that a stable UIO for this distributive target tracking exists, and its optimal parameters can be obtained via solving the above iterative optimization problem. To the best of our knowledge, no results exist in the literature which solve this type of filter design problem. In this simulation, we compare the proposed multi-rate H∞ filter design with given α1 and α2 (Method 1) with a multi-rate H∞ filter design ignoring PD (Method 2) which is equivalent to the proposed method with α1 = α2 = 1. The UI vector wk in (1) and (2) includes three components: the first models the acceleration input; the second and third represent measurement errors of S1 and S2, respectively. In this simulation, the actual acceleration input is a sine signal with amplitude 0.125 and frequency 0.05 rad/s; independent measurement errors of the two sensors obey a truncated zero-mean Gaussian distribution with variance 1 and the truncating interval [−3, 3]. As shown in Fig. 2, Method 1 obtains satisfactory estimates while Method 2 leads to a significant estimate deviation from the actual position. Fig. 3 shows the changes in γ when both α1 and α2 change from 0.1 to 0.9 with a step size 0.1 for ρ = 0.75. It is concluded that a larger γ will result if PDs are more severe.

0.9

0.5

0.6

0.7

0.8

0.3

0.2

0.1

α1

Fig. 3. γ vs. α1 and α2 with ρ = 0.75.

Acknowledgements This research was supported in part by NSERC, the National Natural Science Foundation of China (60634030), the 111 Project (B08015), the 973 Project, the Scientific and Technological Innovation Foundation of NPU, and the Program for New Century Excellent Talents of University, China. Appendix. Derivation of the proposed UIO Putting (7), (10) and (11) into (14), we have mi

e xkN +mi = Eie xkN + A −

p X

! Fi,j αj C j − Ei

xk

j =1

−

p p X X (Fi,j (1 − αj ) + Fi∗,j )y∗j,dk − Fi,j ∆j,k j =1

+



j =1

mi−1

A

···

B

AB

B



0 −

!

p X

Fi,j αj Dj

wk ,

(A.1)

j=1

e xkN +N =

A−

p X

! Fq,j αj C j + Fq+1,j αj Cj A − Eq




xk

j =1

−

p X

Fq,j ∆j,k + Fq+1,j ∆j,k −




p X

j =1

+

Fq+1,j αj Dj wkN +N

j=1

B−

p X

! wk

Fq,j αj Dj + Fq+1,j αj Cj B




j =1

+ Eqe xkN −

p X

Fq,j (1 − αj ) + Fq∗,j y∗d




j =1

6. Conclusions This paper presented the multi-rate state estimation problem with norm-bounded UIs and stochastic PDs. To the best of our knowledge, this is the first consideration of state estimation using multi-rate measurements with PDs. Through constructing a UIO with a causality constraint, state estimation was transformed into the design of a stable UIO. Conditions were proposed for the existence of a stable UIO. Based on the resultant UIO, a multi-rate stochastic H∞ filter was designed. A numerical example to track a maneuvering target was given for illustration. In this paper, the probabilities of successful packet transmission for sensor measurements in (5) are assumed to be time-invariant. If they are time-varying, as in some networks, our proposed estimator design (for time-invariant parameters) may not be applied, because of computational issues; how to design an estimator in this case is an interesting problem for future work.

0.9

0.4

j ,k

p

−

X

Fq+1,j (1 − αj ) + Fq∗+1,j y∗j,mj,k −1 .




(A.2)

j =1

To decouple e xkN +mi and e xkN +N with xk , y∗d , and y∗j,mj,k −1 , we must j ,k

have Ei = Ami −

p X

Fi,j αj C j ,

(A.3)

j=1

Eq = A −

p X

Fq,j αj C j + Fq+1,j αj C j A ,




(A.4)

j =1

Fi∗,j = −(1 − αj )Fi,j ,

Fq∗,j = −(1 − αj )Fq,j ,

Fq∗+1,j = −(1 − αj )Fq+1,j . Putting (3) and (19) into (16), we have

(A.5) (A.6)

444

Y. Liang et al. / Automatica 46 (2010) 437–444

e zkN +mi = Gie xkN +mi + (L − Gi )xkN +mi −

p X

Hi,j αj C j xk

j =1

−

p X

Hi,j (1 − αj ) + Hi,j y∗d


∗

j =1 p

−

X

Hi,j αj Dj w k −

j =1

p X

j ,k

Hi,j ∆j,k ,

(A.7)

j =1

e zkN +N = Gqe xkN +N −

p X

Hq,j αj C j xk

j =1

+

L − Gq −

p X

! Hq+1,j αj Cj

xk+1

j =1

− −

p X j =1 X

Hq,j (1 − αj ) + Hq∗,j y∗d




j,k

p Hq+1,j (1 − αj ) + Hq∗+1,j y∗j,mj,k −1




j =1

−

p X

Hq,j (αj Dj w k + ∆j,k )

j =1 p

−

X

Hq+1,j (αj Dj wkN +N + ∆j,k ).

(A.8)

j =1

To decouple e zkN +mi with xkN +mi , xk , y∗d , w k , and ∆j,k , we must j ,k choose Gi = L,

Hi,j = Hi∗,j = 0.

(A.9)

To decouple e zkN +N with xk+1 , xk , y∗d , y∗j,mj,k −1 , w k , w k+1 , ∆j,k , and

∆j,k+1 , we must choose Gq = L,

Sahebsara, M., Chen, T., & Shah, S. L. (2008). Optimal H∞ filtering in networked control systems with multiple packet dropouts. System and Control Letters, 57(9), 696–702. Sheng, J., Chen, T., & Shah, S. L. (2005). Optimal filtering for multi-rate systems. IEEE Transactions on Circuits Systems II: Express Briefs, 52(4), 228–232. Sinopoli, B., Schenato, L., Franceschetti, M., Jordan, K., & Sastry, S. (2004). Kalman filtering with intermittent observations. IEEE Transactions on Automatic Control, 49(9), 1453–1464. Sun, S., Xie, L., Xiao, W., & Soh, Y. C. (2008a). Optimal linear estimation for systems with multiple packet dropouts. Automatica, 44(5), 1333–1342. Sun, S., Xie, L., Xiao, W., & Soh, Y. C. (2008b). Optimal full-order and reducedorder estimators for discrete-time systems with multiple packet dropouts. IEEE Transactions on Signal Processing, 56(8), 4031–4038. Sun, S., Xie, L., Xiao, W., & Soh, Y. C. (2008c). Optimal filtering for systems with multiple packet dropouts. IEEE Transactions on Circuits and Systems II: Express Briefs, 55(7), 695–699. Valcher, M. E. (1999). State observer for discrete-time linear systems with unknown inputs. IEEE Transactions on Automatic Control, 44(2), 397–401. Wang, Z., Yang, F., Ho, D., & Liu, X. (2006). Robust H∞ filtering for stochastic time-delayed systems with missing measurements. IEEE Transactions on Signal Processing, 54(7), 2579–2587. Zhang, H., Basin, M., & Skliar, M. (2007). Ito–Volterra optimal state estimation with continuous, multi-rate, randomly sampled, and delayed measurements. IEEE Transactions on Automatic Control, 52(3), 401–416. Zhong, M., Ye, H., Ding, S., & Wang, G. (2007). Observer-based fast rate fault detection for a class of multi-rate sampled-data systems. IEEE Transactions on Automatic Control, 52(3), 520–525.

j ,k

Hq,j = Hq∗,j = Hq+1,j = Hq∗+1,j = 0.

(A.10)

Putting (A.6) into (10) and (11), and (A.9) and (A.10) into (12) and (13), we obtain (17)–(19). References Andiusani, D., & Gau, C. (1987). Estimation using a multi-rate filter. IEEE Transactions on Automatic Control, 32(7), 653–656. Basin, M., & Martinez-Zuniga, R. (2004). Optimal linear filtering over observations with multiple delays. International Journal of Robust and Nonlinear Control, 14(8), 685–696. Chen, T., & Francis, B. (1995). Optimal sampled-data control systems. New York: Springer-Verlag. Gao, H., & Chen, T. (2007). H∞ estimation for uncertain systems with limited communication capacity. IEEE Transactions on Automatic Control, 52(11), 2070–2084. Hespanha, J. P., Naghshtabrizi, P., & Xu, Y. (2007). A survey of recent results in networked control systems. Proceedings of the IEEE, 95(1), 138–162. Hong, L. (1999). Multi-rate interacting multiple model filtering for target tracking using multi-rate models. IEEE Transactions on Automatic Control, 44(7), 1326–1340. Hu, S., & Zhu, Q. (2003). Stochastic optimal control and analysis of stability of networked control systems with long delay. Automatica, 39(11), 1877–1884. Izadi, I., Zhao, Q., & Chen, T. (2005). An optimal scheme for fast rate fault detection based on multi-rate sampled data. Journal of Process Control, 15(3), 307–319. Liang, Y., Chen, T., & Pan, Q. (2009). Multi-rate optimal state estimation. International Journal of Control, 82(11), 2059–2076. Luck, R., & Ray, A. (1990). An observer-based compensator for distributed delays. Automatica, 26(5), 903–908. Matveev, A. S., & Savkin, A. V. (2003). The problem of state estimation via asynchronous communication channels with irregular transmission times. IEEE Transactions on Automatic Control, 48(4), 670–676. Nilsson, J., Bernhardsson, B., & Wittenmark, B. (1998). Stochastic analysis and control of real-time systems with random time delays. Automatica, 34(1), 57–64. Oliveira, M.C., & Geromel, J.C. (1997). Numerical comparison of output feedback design methods. In Proceedings of American control conference, Albuquerque, USA, (pp. 72–76). Sahebsara, M., Chen, T., & Shah, S. L. (2007). Optimal H2 filtering with random sensor delay, multiple packet dropout and uncertain observations. International Journal of Control, 80(2), 292–301.

Yan Liang received the B.E. degree from the School of Marine Engineering at Northwestern Polytechnical University (NPU), China, in 1993, and the M.E. and Ph.D. degrees from the Department of Automation at NPU in 1998 and 2001, respectively. He was a Post-Doctoral Fellow at Tsinghua University, China, from April 2001 to July 2003, a Research Fellow at the Hong Kong Polytechnic University from May 2006 to August 2006, and a Visiting Scholar at the University of Alberta from December 2007 to December 2008. Since 2009, he has been a Professor in the School of Automation at NPU. His research interests are in the areas of state estimation, information fusion and target tracking. His e-mail is [email protected]. Tongwen Chen is presently a Professor of Electrical and Computer Engineering at the University of Alberta, Edmonton, Canada. He received the B.Eng. degree in Automation and Instrumentation from Tsinghua University (Beijing) in 1984, and the M.A.Sc. and Ph.D. degrees in Electrical Engineering from the University of Toronto in 1988 and 1991, respectively. His research interests include computer and networkbased control systems, and their applications to the process and power industries. He is an IEEE Fellow and a Fellow of the Engineering Institute of Canada. He has served as an Associate Editor for several international journals, including IEEE Transactions on Automatic Control, Automatica, and Systems and Control Letters. He is a registered Professional Engineer in Alberta, Canada. Quan Pan was born in Shanghai in 1961. He received his bachelor’s degree from Huazhong Institute of Technology in Wuhan in 1991, and his master’s and Ph.D. degrees from the School of Automation at Northwestern Polytechnical University (NPU) in Xi’an in 1991 and 1997, respectively. His research interests are target tracking, information fusion, hybrid system estimation theory, multi-scale estimation theory, and image processing. He is a Member of IEEE, a Member of the International Society of Information Fusion (ISIF), a Board Member of the Chinese Association of Automation, and a Member of Chinese Association of Aeronautics and Astronautics. He was the Duty Dean of Graduate School of NPU from 1998 to 2002, and the Duty Dean of Management School of NPU from 2002 to 2004. He is the Director of Office of Development and Planning of NPU from 2004, the Director of Research Institute of Control and Information of NPU from 1993, and a Professor of the School of Automation of NPU from 1998 onwards. As a Visiting Scholar, he visited Wright State University and Colorado School of Mines in 1992, University of New Orleans in 1993, California State University, Los Angeles in 1997, University of Oxford, University of Cambridge and University of Twente in 1999, University of California Berkeley, Stanford University, University of California Los Angeles, Columbia University, Harvard University and Massachusetts Institute of Technology in 2002, University of Auckland, University of Canterbury, University of Melbourne and University of New South Wales in 2005, and Waterloo University and Ottawa University in 2006. He has published 2 books, almost 20 international journal papers and more then 40 international conference papers. He obtained the 6th National Youth Award for Outstanding Contribution to Science and Technology in 1998 and the Chinese National New Century Excellent Professional Talent in 2000. His e-mail is [email protected].