Distributed Covariance Intersection Fusion Estimation ... - IEEE Xplore

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/TAC.2016.2539221, IEEE Transactions on Automatic Control 1

Distributed Covariance Intersection Fusion Estimation for Cyber-Physical Systems with Communication Constraints Bo Chen, Guoqiang Hu, Member, IEEE, Daniel W. C. Ho, Senior Member, IEEE, Li Yu, Member, IEEE

Abstract—This paper addresses the distributed covariance intersection fusion estimation problem for a class of cyber-physical systems (CPSs), where the raw measurements are preprocessed in each sink node to obtain the local state estimates of a CPS. When each local estimate is transmitted to a remote fusion center (FC), the communication network between the sink node and the FC is subject to finite bandwidth, varying delays and packet dropouts. Based on a new communication constraint model covering all these network phenomena, a recursively distributed fusion estimator (DFE) is designed for the addressed CPSs. Since the estimation performance directly impacts the stability of control operation in CPSs, stability condition is derived such that the mean square error of the designed DFE is bounded, and it is proved that the DFE is independent of the initial values under certain conditions. The design procedure is illustrated using a numerical example. Index Terms—Cyber-Physical Systems, Distributed Fusion Estimation, Kalman Filtering, Communication Constraints, Stability Analysis

I. I NTRODUCTION Advances in embedded computing, communication, and related hardware technologies have recently brought the paradigm of cyberphysical system (CPSs) as a new research frontier [1]. As one of important issues in CPSs, real-time state estimation plays a vital role in provisioning real-time monitoring and control [2]–[4]. For example, estimating the real voltage from sensor information must be completed before taking certain actions to regulate the voltage with in some desired range in a power grid [5]. Moreover, the accuracy of state estimation has an important impact on computing control commands for safe and efficient operation of a CPS [3]. Since multi-sensor information fusion estimation can potentially improve estimation accuracy and enhance reliability and robustness against sensor failures [6]–[12], it is of theoretical significance and practical relevance to investigate the problem of information fusion estimation for the CPSs. In this case, the distributed information fusion estimation problem is investigated in this paper for a class of CPS architecture (see Fig.1 in [13]), where a typical example of such CPS is the smart grid communication systems [14]. When the local estimates are transmitted to the fusion center (FC) over communication networks, however, three key issues must be taken into account [15]: i) Bandwidth constraint; ii) Packet dropout; iii) Communication delay. Notice that information loss is inevitable under the constraints (i)–(iii), and such a fusion estimation with incomplete information may degrade the estimation performance. In this sense, how to develop fusion estimation algorithms in the presence of the above communication constraints is essential for real-time monitoring and control of the CPSs. This work was supported in part by Singapore MOE AcRF Tier 1 Grant RG60/12(2012-T1-001-158), in part by the National Natural Science Funds of China (61403345, 61273117), in part by the Hong Kong Scholars Program (XJ2014013), in part by GRF grant from HKSAR (CityU 11300415), and in part by the Zhejiang Province Postdoctoral Research Project (BSH1502151). B. Chen was with the School of Electrical and Electronic Engineering, Nanyang Technological University, 639798 Singapore, is currently with the Department of Mathematics, City University of Hong Kong, 999077 Hong Kong, and is also with the College of Information Engineering, Zhejiang University of Technology, HangZhou 310023, China (email:[email protected]). G. Hu is with the School of Electrical and Electronic Engineering, Nanyang Technological University, 639798 Singapore (email: [email protected]). D. W. C. Ho is with the Department of Mathematics, City University of Hong Kong, Hong Kong, 999077, China (email:[email protected]). L. Yu is with the College of Information Engineering, Zhejiang University of Technology, HangZhou 310023, China (email:[email protected]).

Under the centralized fusion framework, the Kalman filter over a lossy network was designed in [16] for a class of spatially distributed CPSs with packet losses. Additionally, the centralized fusion estimation methods in [17], [18] may solve the estimation problem for CPSs with communication delays and packet dropouts. Compared with the centralized fusion structure, the distributed fusion structure is generally more robust, reliable, and fault-tolerant [6]–[8], and thus can provide an attractive alternative to fusion estimation problems for CPSs. Under this case, the distributed fusion estimators in [19]–[21] may be applicable to the CPSs with delays and packet dropouts. Notice that the number of state variables for a CPS may be large, and thus it is unrealistic in the distributed fusion structure that each local estimate for the state “x(t)” is completely transmitted to the FC over a bandwidth constrained communication channel. In this sense, bandwidth constraint in CPSs is the primary consideration for designing a distributed state fusion estimator. Due to this constraint, the distributed fusion estimation algorithms in [19]–[21] cannot be directly applied to the CPSs. Most of the existing results on fusion estimation with bandwidth constraint resort to the dimensionality reduction method (see [22]–[24] and the references therein) and the quantization method (see [25], [26] and the references therein). Moreover, for a high-dimensional signal, the dimensionality reduction method may be efficient in traffic reduction as compared with the quantization method [27], and is more suitable to solve the problem of bandwidth constraint in CPSs. Therefore, the distributed dimensionality reduction fusion estimation algorithm was developed in [13] for the CPSs with constant communication delays. Without the assumptions on the statistical property of delays and packet dropouts in [18] and [21], it is considered in this paper that the communication delays and the number of successive packet losses (NSPL), which can be identified by each packet’s timestamp, are time-varying from the sink nodes to the FC. The main contributions are summarized as follows: i) A novel communication constraint model with compensation strategy is proposed to describe the limited bandwidth, varying communication delays and packet dropouts, and then the recursive distributed fusion estimator (DFE) is designed based on covariance intersection (CI) fusion criterion; ii) A delay-dependent and probability-dependent condition is derived such that the mean square error (MSE) of the designed DFE is bounded, and it is proved that the estimation performance of the designed DFE is independent of the initial values under certain conditions. Moreover, when the upper bounds of the delays and NSPL are known in advance, the probability selection criterion is presented for the stochastic dimensionality reduction strategy. Finally, the advantage of the DFE is demonstrated by comparing with the fusion estimators in [16] and [22] by a numerical example. Notations: The superscript ′ T′ represents the transpose, and E{·} is the mathematical expectation. Im represents the identity matrix of size m × m, while diag{·} stands for a block diagonal matrix. Prob{A} means the occurrence probability of the event A, while Tr(B) denotes the trace of the matrix B. ||A||1 , ||A||2 and ||A||∞ represent the 1-norm, 2-norm and infinity-norm of the matrix A, respectively. x⊥y denotes that x and y are orthogonal vectors. Moreover, τ1 > τ2 , it ∏ ∑τif 2 will be specified that ττ2=τ1 F (τ ) = Im and τ =τ1 G(τ ) = 0, where F (τ ) ∈ Rm×m and G(τ ) ∈ Rn×n represent different matrix functions with respect to the variable τ .

0018-9286 (c) 2015 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.


where σ~i i (t)(~i = 1, 2, · · · , ∆i ) are required to satisfy ∑∆i σ~i i (t)σ~i 0 (t) = 0(~i ̸= ~0i ) and σ~i i (t) = 1

II. P ROBLEM F ORMULATION A. System Model

~i =1

i

Consider a class of physical processes described by the following discrete state-space model: x(t + 1) = Ax(t) + Γw(t)

(1)

where x(t) ∈ Rn is the current state of the process, w(t) ∈ Rm is the system noise, while A and Γ are constant matrices. When the measurements from each sensor are sent to the sink nodes, the ith sink node’s measurement yi (t) ∈ Rqi is modeled by: yi (t) = Ci x(t) + vi (t)(i = 1, 2, · · · , L)

(2)

where Ci ∈ Rqi ×n is the measurement matrix, and vi (t) ∈ Rqi is the measurement noise. Moreover, w(t) and vi (t) are uncorrelated zero-mean Gaussian white noises satisfying E{[wT (t) viT (t)]T [wT (t1 ) vjT (t1 )]} = δt,t1 diag{Qw , δi,j Qvi }, where δt,t1 = 0 if t ̸= t1 and δt,t1 = 1 otherwise. Each sink node in Fig.1 of [13] is a gateway node, which is responsible for receiving measurements, computing the local optimal (linear minimum variance sense) estimate (LOE) of the system (1) and sending the LOE to the FC. Then, based on the measurements {yi (1), yi (2), · · · , yi (t)}, the LOE x î (t) is given by the well-known Kalman filter: x î (t) = GKi (t)Aˆ xi (t − 1) + Ki (t)yi (t) ∆

(3)

∆

where GKi (t) = In − Ki (t)Ci . Define x ˜i (t) = x(t) − x î (t). Then the optimal gain matrix Ki (t) and the local estimation error covari∆ ance matrix Pii (t) = E{˜ xi (t)˜ xT i (t)} are calculated by  ∗ T ∗ T −1  Ki (t) = Pii (t)Ci [Ci Pii (t)Ci + Qvi ] Pii (t) = GKi (t)Pii∗ (t) (4)  ∗ Pii (t) = APii (t − 1)AT + ΓQw ΓT B. The Model of Communication Constraints Under the distributed fusion framework, each LOE x î (t) must be sent to the FC to design an optimal fusion estimator. Notice that the number of state variables in (1) is large in many large-scale CPSs, however, any communication network can only carry a finite amount information per unit time. Thus, it is unrealistic to send the complete information of x î (t) to the FC over communication networks. To reduce communication traffic, only ri (1 ≤ i < n) components of the ith LOE x î (t) are allowed to be transmitted to the FC, and the allowed sending components (ASC) of each x î (t) are selected in a random way. The above strategy can reduce the size of the data packet containing messages from the sink node, and thus satisfy the finite communication bandwidth. According to this communication strategy, the/ ASC of x î (t) has ∆i possible cases, ri∏ −1 ri ∏ n r where ∆i = (n − ℓi ) ℓi . For the ith LOE x î (t), let ℓn i =0

ℓr i =1

us define the selected ASC (SASC) by x ˆsi (t) ∈ Rri that will be transmitted to the FC, while one group of ASC is denoted by x ˆ~sii (t) ∈ Rri . Then the SASC x ˆsi (t) only takes one ASC at time t from the finite set Ssi (t) = {ˆ x~sii (t)|~i = 1, 2, · · · , ∆i }. Suppose that each SASC x ˆsi (t) can be received by the FC at time t, it is specified that the unselected components of x î (t)(i ∈ {1, · · · , L}) with respect to x ˆsi (t) are replaced by “0” in the FC. In this case, according to the set Ssi (t), the reorganized state estimate (RSE) of the system state x(t) at the FC side, which is denoted by x âi (t), only takes one signal from the following set: Si (t) = {H~i i x î (t)|~i = 1, 2, · · · , ∆i } H~i i

(5)

where denotes a diagonal matrix that contains ri diagonal elements “1” and n − ri diagonal elements “0”. To describe the RSE x âi (t) in a simple way, it is assumed that the ∆i elements of the set Ssi (t) are indexed from 1 to ∆i , then we introduce the following indication functions: { 1 if x ˆsi (t) = x ˆ~sii (t) i σ~i (t) = (6) 0 if x ˆsi (t) ̸= x ˆ~sii (t)

(7)

such that the SASC x ˆsi (t) ∈ Rri only takes one ASC from the set Ssi (t), i.e., ∑∆i x ˆsi (t) = σ~i i (t)ˆ x~sii (t) (8) ~i =1

From (5) and (8), the RSE x âi (t) is modeled by: x âi (t) = Hi (t)ˆ xi (t) where ∑∆i Hi (t) =

~i =1

σ~i i (t)H~i i = diag{γ1i (t), · · · , γni (t)}

(9) (10)

∑n i which yields that γℓi (t) ∈ {0, 1} and ℓ=1 γℓ (t) = ri (i = 1, · · · , L). It is noted that when each SASC x ˆsi (t) is transmitted to the FC through communication networks, however, some data packets not only suffer time-varying communication delay but also, even worse, may be lost during transmission in the unreliable communication networks. Under this condition, the RSE model (9) will not be true, and the FC may receive a data packet containing the information x ˆsi (t − di (t)) at time t, where di (t) is an integer representing the communication delay between the ith sink node and the FC. On one hand, to identify the time-varying communication delay at the FC side, a time stamp is added to each SASC x ˆsi (t) before it is transmitted into the network medium. On the other hand, there may be no data packet from the ith sink node arriving at the FC at a particular time, and then the NSPL is denoted by τi (t) at time t. To establish the mathematical relationship between the communication delay di (t) and the NSPL τi (t), the following recursive function is introduced: { di (t) if the FC receives a data packet (11) dpi (t) = dpi (t − 1) + 1 otherwise Then the NSPL τi (t) is determined by τi (t) = min{(t − t∗ )|dpi (t∗ ) − di (t∗ ) = 0(t∗ ≤ t)}, and thus one has dpi (t) = di (t − τi (t)) + τi (t). In practical applications, it is reasonable to consider that the time-varying delay di (t) and the NSPL τi (t) are bounded, i.e., dli ≤ di (t) ≤ dui , 0 ≤ τi (t) ≤ τi

(12)

Then, let the ith compensating state estimate (CSE) of the state x(t) under the bandwidth constraints, varying delays and packet dropouts be x ˆci (t), and then each CSE is given by: p

x ˆci (t) = Adi (t) Hi (t − dpi (t))ˆ xi (t − dpi (t)) p p di (t) +A [In − Hi (t − di (t))]Aˆ xci (t − dpi (t) − 1)

(13)

where Hi (t−dpi (t)) is determined by (10–11), and the compensation strategy in (13) is designed based on the dpi (t) − step prediction information of x ˆci (t − dpi (t)). Moreover, the SASC x ˆsi (t) can be only observed at most once in the FC for a practical system, and thus one has x î (t − di (t)) ̸= x î (t1 − di (t1 )) and Hi (t − di (t)) ̸= Hi (t1 − di (t1 )) (t ̸= t1 ). Remark 1. According to the time stamps, if the buffer in the FC is updated only when the arrived data packet is more recent than the existing one, then there must be “t1 − dpi (t1 ) ≤ t2 − dpi (t2 )(t1 < t2 )”. By such a mechanism, the time sequence of x î (t− dpi (t)) in (13) is correct, and thus it is easy to find a fixed linear ∆ span of the vectors {˜ xci (1), · · · , x ˜ci (t1 )} (˜ xci (t) = x(t)−ˆ xci (t)) such that the random variables w(t) and vi (t) are orthogonal to this span, which means that the recursive structure of the fusion estimator can be easily determined. However, the RSE model (9) implies that the variable dpi (t) in (13) only represents the communication delay and NSPL of the transmitted components, while the delay and NSPL for the untransmitted components are unknown. In this sense, the above storing mechanism is conservative for dealing with the case of bandwidth constraints. By taking this factor into account, it is specified in this paper that the arrived data packet is directly stored



into the buffer to design the DFE at a particular time, and thus the time sequence of x î (t − dpi (t)) is out of order, which adds new difficulties to the design of the recursive DFE. C. Problem of Interest When a group of ASC x ˆ~sii (t)(~i ∈ {1, 2, · · · , ∆i }) in the set Ssi (t) is randomly selected as the SASC x ˆsi (t) at time t, let each stochastic process {σ~i i (t)}(~i ∈ {1, 2, · · · , ∆i }) be independent and identically-distributed (i.i.d) when designing the stochastic dimensionality reduction strategy. Then, the occurrence probabilities of the cases σ~i i (t) = 1 and σ~i i (t) = 0 are given by Prob{σ~i i (t) = 1} = π~i i and Prob{σ~i i (t) = 0} = 1 − π~i i , where the selection probability π~i i ≥ 0 is a given scalar satisfying: ∑∆i π~i i = 1 (i ∈ {1, 2, · · · , L}) (14) ~i =1

Therefore, it is concluded from (10) that the binary variables γℓi (t)(ℓ = 1, · · · , n) are independent Bernoulli dis∆ tributed white noise sequences with Prob{γℓi (t) = 1} = ∆ γℓi and Prob{γℓi (t) = 0} = 1 − γℓi . Moreover, it is specified that each random variable σ~i i (t) is uncorrelated with w(t) and vi (t)(i = 1, · · · , L). On the other hand, according to the CSEs x ˆci (t)(i = 1, · · · , L) in the FC, the DFE for the addressed CPSs is given by: x ˆ(t) =

∑L i=1

Ωi (t)ˆ xci (t)

(15)

∑L where ˆ(t) is unbiased i=1 Ωi (t) = In , and thus the DFE x when E{ˆ xci (−dni )} = E{x(−dni )} = E{ˆ xi (−dni )}(dni = 0, 1, · · · , dui + τi ). Consequently, the problems to be solved are described as follows: 1) For a given arbitrary group of random binary variables ςi (t) = i (σ1i (t), · · · , σ∆ (t))(i = 1, · · · , L) satisfying (7) and (14), the aim i is to design the optimal weighting matrices Ω1 (t), · · · , ΩL (t) such that the MSE of DFE x ˆ(t) is minimal, i.e., x ˆ(t) = arg min E{[x(t) − x ˆ∗ (t)]T [x(t) − x ˆ∗ (t)]} x ˆ∗ (t)

(16)

where x ˆ∗ (t) denotes an arbitrary group of convex linear combination with respective to the CSEs x ˆci (t)(i = 1, 2, · · · , L). 2) Find stability conditions, which are dependent on the selection probability π~i i in (14), the communication delay di (t) and the NSPL τi (t), such that the MSE of the designed DFE x ˆ(t) is bounded and independent of the initial values, i.e., J = lim E{[x(t) − x ˆ(t)]T [x(t) − x ˆ(t)]} < p t→∞

(17)

and J is independent of the initial values. Remark 2. When the sensor messages are transmitted to the FC over communication networks, the communication delays and NSPL are always time-varying in a practical system. Under this case, there are mainly two different strategies designing the fusion estimators: S.i) By introducing certain extent of conservatism, each time-varying delay and each NSPL are prolonged to the upper bound dui and τi at the FC side, and then the fusion estimator is designed based on the constant delays and bandwidth constraints; S.ii) At each time step, the fusion estimator gains are designed according to the time-varying communication delay and NSPL that are identified at the FC side, and this method does not introduce the conservatism. In this paper, we will adopt the second strategy to design the DFE with communication constraints. It is noted that though the distributed fusion estimator was designed in [13] for CPSs with constant communication delays and bandwidth constraints, however, the method in [13] cannot solve the addressed problem in this paper by using (S.i) because the dimensionality reduction strategy in [13] at the sink node side was required to know communication delay at each time step, but only the FC knows the time-varying communication delays in practical applications.

III. D ISTRIBUTED F USION E STIMATOR F OR CPS S Denote the local estimation error covariance matrix and the local estimation error cross-covariance matrix of x ˆci (t) and x ˆcj (t) as ∆ ∆ c c T c c T Σii (t) = E{˜ xi (t)[˜ xi (t)] } and Σij (t) = E{˜ xi (t)[˜ xj (t)] }(i ̸= j), ∆ where x ˜ci (t) = x(t) − x ˆci (t). If each optimal weighting matrix Ωi (t) in (15) is designed by using the similar equation (21) in [13], one must establish the recursive relationship between Σij (t) and the variable E{˜ xci (t − dpi (t) − 1)[˜ xcj (t − dpj (t) − 1)]T } before obtaining the recursive form of Σij (t). However, the variables dpi (t) and dpj (t) are time-varying, and thus one cannot guarantee that there must exist positive integers ~1 > 0 and ~2 > 0 at time t such that ∆ ~1 ~2 fio (t) = fjo (t)(i ̸= j), where fi (t) = t − dpi (t) − 1, ~ 0 fio (t) = fi (fi (· · · (fi (t)) · · · )) and fio (t) = t. This implies that the | {z } ∆

∆

~ times

DFE (15) cannot be designed by using the optimal fusion criterion weighted by matrices. In this case, the covariance intersection (CI) fusion criterion is adopted to design the DFE x ˆ(t) for the communication constraint model (13), where the CI fusion criterion is not required to know the cross-covariance matrix Σij (t) [9]. In this fusion framework, the weighting matrices Ωi (t)(i = 1, 2, · · · L) is given by: Ωi (t) = P (t)ηi (t)Σ−1 (18) ii (t)(i = 1, 2, · · · , L) ∑L where 1 (ηi (t) ≥ 0) and P (t) = i=1 ηi (t) ) = (∑ −1 L −1 . Then the optimal coefficients ηi (t)(i = i=1 ηi (t)Σii (t) 1, 2, · · · , L) are obtained by solving the following problem: { min Tr(P (t)) ηi (t)(i=1,··· ,L) (19) ∑L s.t. : i=1 ηi (t) = 1, 0 ≤ ηi (t) ≤ 1 Though (19) is a nonlinear optimization problem, it can be solved by the function “fmincon” in Matlab Optimization Toolbox when the local estimation error covariance matrices Σii (t)(i = 1, 2, · · · , L) are given. Therefore, the recursive form of Σii (t) will be presented in the following theorem. Theorem 1: Define  ∆ ∆  Φ (t) = GKi (t)A, Hi = E{Hi (t)}   Ki ∏  ∆   ΨΦi (t, ξ) = ξφi =0 ΦKi (t − φi )  ∆ (20) Λi = E{Hi (t) ⊙ Hi (t)}   ∆ V =  E{Hi (t) ⊙ [In − Hi (t)]} i    ∆ Wi = E{[In − Hi (t)] ⊙ [In − Hi (t)]} where GKi (t) is defined in (3), while the operators “⊙” and “⊗” are given by Lemma 1 [22]. Then, each local estimation error covariance matrix Σii (t) is calculated by:  p p  Σii (t) = Adi (t) [Wi ⊗ (AΣii (fi (t))AT )][Adi (t) ]T    ∆   +∆Σii (t) (fi (t) = t − dpi (t) − 1)   p  d (t)  ∆Σii (t) = A i [Λi ⊗ Pii (fi (t) + 1) + Wi    ⊗(ΓQw ΓT ) + Vi ⊗ [GKi (fi (t) + 1)ΓQw ΓT (21)  +ViT ⊗ (ΓQw ΓT GT Ki (fi (t) + 1))   p  12 T di (t) T  +Σ12 ] + Σ22  ii (t) + (Σii (t)) ][A ii (t)   p ∑dpi (t) dp (t)−1  22 T di (t)−1 T  Σii (t) = i A ΓQ Γ (A )  w θ=1   12 Σii (t) = Vi ⊗ [ΨΦi (fi (t) + 1, 0)Φxi (fi (t))AT ] ∆

where Pii (t) and GKi (t) are determined by (4); Φxi (t) = E{˜ xi (t)[˜ xci (t)]T } is computed by:  p p d (t) T   Φxi (t) = ΨΦi (t, di (t))Φxi (fi (t))[In − Hi ][A i ]    +∆Φxi (t)  ∆Φxi (t) = {ΨΦi (t, dpi (t) − 1)GKi (fi (t) + 1)ΓQw ΓT [In (22) p  p   −Hi ] + ΨΦi (t, di (t) − 1)Pii (fi (t) + 1)Hi }[Adi (t) ]T  p  ∑ d (t)  T + αii =1 ΨΦi (t, αi − 2)GKi (t − αi + 1)ΓQw ΓT (Aαi −1 ) Moreover, the estimation performances of the local CSEs (13) and the DFE (15) have the following relationship: Tr(Po (t)) ≤ Tr(P (t)) ≤ Tr(Σii (t))(i ∈ {1, · · · , L})

(23)



∆

where Po (t) = E{˜ x(t)˜ xT (t)} denotes the fusion estimation error covariance matrix of the DFE x ˆ(t). Proof: It follows from (3) and (13) that x ˜i (t) = ΨΦi (t, t − t2 − 1)˜ xi (t2 ) ∑ 2 + t−t Ψ Φi (t, αi − 2)GKi (t − αi + 1)Γw(t − αi ) (24) α =1 ∑ i 2 −1 − t−t ΨΦi (t, αi − 1)Ki (t − αi )vi (t − αi ) αi =0 p

x ˜ci (t) = Adi (t) Hi (t − dpi (t))˜ xi (t − dpi (t)) p p di (t) +A {[In − Hi (t − di (t))]A˜ xci (fi (t)) +[In − Hi (t − dpi (t))]Γw(fi (t))} + Fw (dpi (t), t)

(25)

p ∆ ∑di (t) θ−1 where Fw (dpi (t), t) = θ=1 A Γw(t − θ), while ΨΦi (t, ξ) is defined in (20). Meanwhile, according to the geometric meaning of x ˜i (t) and x ˜ci (t), one has by the statistical properties of w(t), vi (t) and Hi (t) that  ˜i (t)⊥w(t1 )(t1 ≥ t), x ˜i (t)⊥vi (t1 )(t1 > t)  x  x ˜ci (t)⊥w(t1 )(t1 ≥ t), x ˜ci (t)⊥vi (t1 )(t1 > t) (26) x ˜i (t − dpi (t))⊥Fw (dpi (t), t)    c p p x ˜i (t − di (t) − 1)⊥Fw (di (t), t)

Then it is derived from (24–26) that  p p xi (t)˜ xT  E{˜ i (t − di (t))} = ΨΦi (t, di (t) − 1)Pii (fi (t) + 1) p c T (27) E{˜ xi (t)[˜ xi (fi (t))] } = ΨΦi (t, di (t))Φxx (fi (t))  E{˜ xi (t)wT (fi (t))} = ΨΦi (t, dpi (t) − 1)GKi (fi (t) + 1)ΓQw Therefore, (22) is obtained from (24–25) and (27). On the other hand, the recursive equation (21) can be deduced from (22), (24– 26) and Lemma 1 in [22], and the detailed derivation is omitted due to the page limitation. Moreover, the result (23) can be obtained from Theorem 2 in [12], which implies that though the estimation error cross-covariance matrix Σij (t) cannot obtained for the communication constraint model (13), the estimation performance of the designed DFE x ˆ(t) is still better than that of each CSE x ˆci (t). Based on Theorem 1, the computation procedures for the DFE x ˆ(t) can be summarized by Algorithm 1.

controllable and observable, i.e., rank([Γ, AΓ, · · · , An−1 Γ]) = n and rank([CiT , (Ci A)T , · · · , (Ci An−1 )T ]T ) = n (i = 1, · · · , L). Theorem 2: For the addressed CPSs (1–2) with communication constraint model (13), if the system parameters, selection probabilities, communication delays and NSPL satisfy the following condition:  ∆ χ χ l l   ζ∞,i = max{(||A ||1 ||A ||∞ )|χ = di , di + 1, · · · ,   u 1/2 di + τi }||(In − Hi ) A||1 ||(In − Hi )1/2 A||∞ < 1 (28) ∆ χ+1   ζ2,i = max{(||ΦKi ||2 ||Aχ (In − Hi )A||2 )   χ = dli , dli + 1, · · · , dui + τi } < 1 (i = 1, · · · , L) where dli , dui , τi and Hi are respectively defined in (12) and (20), and ∆ ΦKi = lim ΦKi (t), then each estimation error covariance matrix t→∞ Σii (t) calculated by (21) can guarantee: lim ||Σii (t)||1 ≤ Σci , lim Σii,1 (t) = lim Σii,2 (t)

t→∞

t→∞

t→∞

(29)

where Σii,1 (t) and Σii,2 (t) are any estimation error covariance matrices of the CSE x ˆci (t) with different initial conditions. In this case, the MSE of the designed DFE x ˆ(t) will be bounded and independent of the initial conditions, i.e., ∆

J = lim Tr{Po (t)} ≤ p

(30)

t→∞

where J is independent of the initial conditions. Proof: According to the well-known stability condition of Kalman filter, the controllable and observable systems (1–2) imply: lim Pii (t) = Pii , lim GKi (t) = GKi

t→∞

t→∞

lim ΦKi (t) = ΦKi , lim ΨΦi (t, ξ) = Φξ+1 Ki

t→∞

(31)

t→∞

where the limits in (31) are independent of the initial conditions, while ΦKi is a stable matrix. Then it is deduced from (31) that there must exist an integer NΨi > 0 such that (22) is equivalent to: p

d (t)+1

Φxi (t) = ΦKii

p

ˆ xi (t)(32) Φxi (fi (t))[In − Hi ][Adi (t) ]T + ∆Φ p

Algorithm 1 For a given arbitrary group of random binary variables ςi (t)(i = 1, 2, · · · , L) satisfying i.i.d and (14) 1: Calculate each LOE x î (t) by (3), while each local estimation error variance matrix Pii (t) is computed by (4); 2: Calculate Φxi (t) by (22); 3: Calculate Σii (t) by (21); 4: Determine η1 (t), · · · , ηL (t) by solving the problem (19); 5: Calculate Ω1 (t), · · · , ΩL (t) by (18); 6: Calculate the local CSEs x ˆc1 (t), · · · , x ˆcL (t) by (13); 7: Calculate the DFE x ˆ(t) by (15).

{ΦKii GKi ΓQw ΓT [In − Hi ] + p ∑dpi (t) T αi −1 di (t) T Φ Pii Hi }[A ] + GKi ΓQw ΓT (Aαi −1 ) . αi =1 ΦKi ˆ xx (t) in (32) is bounded, i.e, there Combining (12) yields that ∆Φ ˆ xi (t)||2 ≤ ∆Ψx , must exist a scalar ∆Ψxi such that lim ||∆Φ i t→∞ ˆ and ∆Φxi (t) is independent of the initial conditions. Therefore, it follows from (32) that

Remark 3. It is concluded from Theorem 1 that the computation procedures for the local estimation error covariance matrices Σii (t)(i = 1, 2, · · · , L) are independent of the sensor measurements, and thus the optimization problem (19) can be separately solved in the FC when the selection probabilities π~i i (i = 1, · · · , L; ~i = 1, · · · , ∆i ) and the initial conditions are given. In this case, if only each SASC x ˆsi (t) ∈ Rri (ri < n) is sent to the FC, Algorithm 1 will be easily implemented in practical applications, where x ˆsi (t) is determined by the binary variables σ~i i (t) (~i = 1, · · · , ∆i ). Notice that when each sink node knows the selection probabilities, these binary variables obeying the i.i.d will be randomly generated. On the other hand, when the estimation error cross-covariance Σij (t) is difficult to obtain, the standard CI fusion method is adopted in this paper to design the DFE, however, the distributed fusion estimator for the communication constraint model (13) can also been obtained by using the modified CI fusion methods in [10]–[12] based on Theorem 1.

~0 = min{~∗0 |fio0 (t) − NΨi ≥ 0}. Notice that t → ∞ ⇔ ~ → ∞, and thus it is concluded that when the second condition in (28) holds, Φxi (t) will be bounded. At the same time, let Φ1xi (t) and Φ2xi (t) denote any covariance matrices that are computed by (22) ~0 ~0 with different initial conditions Φ1xi (fio (t)) and Φ2xi (fio (t)). Then ˆ xi (t) that one has by (32) and the property of ∆Φ p ∏~ p κ d (f κ (t))+1 κ ˜ xi (t) = ˜ xi (fio Φ ΦKii io Φ (t))[In − Hi ]Adi (fio (t))(34)

IV. S TABILITY A NALYSIS FOR THE DFE

and ∆Σii (t) is independent of the initial conditions. Then, it follows from (35) that there must exist an integer NΣi > NΨi such that

In this section, we will discuss the performance of the designed DFE. Without loss of generality, it is assumed that (1) and (2) are

where

ˆ xi (t) ∆Φ

=

d (t)

p di (t) Ki

||Φxi (t)||2 ≤ ζ2,i ||Φxi (fi (t))||2 + ∆Ψxi (t > NΨi )

(33)

where ζ2,i is defined in (28). Then, ∑ one0 +1hasκ by (33) that ~−~0 ~0 ||Φxi (t)||2 ≤ ζ2,i ||Φxi (fio (t))||2 + ~−~ ζ2,i ∆Ψxi , where κ=0 ~∗

κ=0

∆ ˜ xi (t) = where Φ Φ1xi (t) − Φ2xi (t). Then, taking the second condition ˜ xi (t) = 0, and thus the (28) and (34) into account yields lim Φ t→∞ matrix variable Φxi (t) (22) is independent of the initial conditions. Since the matrix variable ∆Σii (t) in (21) is dependent on the variables Pii (t), GKi (t), ΨΦi (t, ξ) and Φxi (t), it can be concluded that there must exist a scalar ∆Σci > 0 satisfying:

lim ||∆Σii (t)||1 ≤ ∆Σci

(35)

t→∞

p

p

ˆ ii (t)(36) Σii (t) = Adi (t) [Wi ⊗ (AΣii (fi (t))AT )][Adi (t) ]T + ∆Σ


This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/TAC.2016.2539221, IEEE Transactions on Automatic Control 5 p

ˆ ii (t) = Adi (t) [Λi ⊗ Pii + Wi ⊗ (ΓQw ΓT ) + where ∆Σ ˆ 12 Vi ⊗ (GKi ΓQw ΓT ) + ViT ⊗ (ΓQw ΓT GT Ki ) + Σii (t) + p 12 T d (t) T 22 12 ˆ ii (t)) ][A i ] + Σii (t) and Σ ˆ ii (t) (Σ = Vi ⊗ [ΦKi Φxi (fi (t))AT ]. According to the statistical property “i.i.d” on the binary variable σ~i i (t), it is known from (7) and (10) that the diagonal elements of Hi (t) satisfy: √ (37) E{γℓi1 (t)γℓi2 (t)} ≤ γℓi1 γℓi2 ≤ γℓi1 γℓi2 (ℓ1 ̸= ℓ2 ) For a given matrix X ∈ Rn×n , one has by (37) and the definitions of “⊙” and “⊗” in Lemma 1 of [22] that ||Wi ⊗ X||1 ≤ ||(In − Hi )1/2 X(In − Hi )1/2 ||1

(38)

It follows from (35) and (36) that ||Σii (t)||1 ≤ ||Wi ⊗ (AΣii (fi (t))AT )||1 p p ×||Adi (t) ||1 ||[Adi (t) ]T ||1 + ∆Σci

(39) p

Then, it is derived from (38–39) and the fact ||[Adi (t) ]T ||1 = p di (t) ||A ||∞ that ||Σii (t)||1 ≤ ζ∞,i ||Σii (fi (t))||1 + ∆Σci ∑~−~1 +1 κ (40) ~−~1 ~1 ≤ ζ∞,i ||Σii (fio (t))||1 + κ=0 ζ∞,i ∆Σci ~∗

where ζ∞,i is defined in (28), and ~1 = min{~∗1 |fio1 (t) − NΣi ≥ 0}. Thus, when the first condition in (28) holds, one ∆Σci has lim ||Σii (t)||1 ≤ 1−ζ∞,i , and thus the local estimation error t→∞ variance matrix Σii (t) is bounded. Moreover, it is deduced from ˜ xi (t) = 0” that lim ||Σii,1 (t) − the similar derivation of “ lim Φ t→∞ t→∞ Σii,2 (t)||1 = 0, and thus the result (29) is obtained. On the other hand, it is known from (18) and (19) that each weighting matrix Ωi (t) is determined by the local estimation error covariance matrices Σii (t)(i = 1, 2, · · · , L), then it is concluded from the second equation in (29) that the weighting matrices are independent of the initial conditions. Based on this result, (30) is obtained from (23) and (29). This completes the proof. Remark 4. It is concluded from (13) and (21) that each CSE and the corresponding estimation error covariance matrix are designed by recursive forms, and thus the Algorithm 1 is required to know the initial conditions of the CSEs and the covariance matrices before being implemented. However, the initial values may be poorly known and are sometimes rather arbitrary in a practical system (e.g., the FC may difficultly know the initial conditions of the spatially distributed sensors). In this sense, the result in Theorem 2 has important engineering significance because the effects of the initial values decay as more and more data is processed if the condition (28) is satisfied for the communication constraint model (13). On the other hand, according to Theorem 2, when ∆Σii (t) in (21) and ∆Φxi (t) in (22) are respectively replaced by ˆ ii (t) of (36) and ∆Φ ˆ xi (t) of (32), a asymptotically optimal local ∆Σ covariance matrix, which is denoted by Σaii (t), can be obtained (i.e., lim Σaii (t) = lim Σii (t)), and the asymptotically optimal DFE t→∞ t→∞ (AODFE) can be designed based on Σaii (t). Moreover, the estimation performance of the AODFE will be equivalent to that of the DFE as t goes to ∞, while the computational complexity of the AODFE is lower than that of the DFE. Remark 5. The stability condition in Theorem 2 is dependent on the selection probabilities, the upper bounds of delays and NSPLs, which implies that the bandwidth constraints, delays and packet dropouts can degrade the estimation performance of the DFE. Furthermore, if the upper bounds of delays and NSPLs are known a priori, the selection probabilities defined in (14) for the dimensionality reduction strategy can be adjusted according to the stability condition in Theorem 2 such that the MSE of the designed DFE is bounded. V. N UMERICAL E XAMPLES Suppose that there are two sink nodes sending their information to the FC for estimating the state of a CPS, and the system parameters

in (1–2) are taken as follows:      1.0 0.1 0 0 0.2     0.2 0.1 0.2 0   0.5           A =  0 0.2 0.1 0.1  , Γ =  0.3     0.4  0 0 0.1  0.9  1010 1001     0 1 0 0 1 0 1 0       C1 =    0 1 1 0  , C2 =  0 0 0 1     1010 1010

(41)

where λmax (A) = 1.0228 > 1 means that the system (1) is unstable. The noises w(t) and vi (t)(i = 1, 2) are gaussian white noises, and their covariance matrices are taken as Qw = 1, Qv1 = diag{0.9, 0.6, 0.9, 0.4} and Qv2 = diag{0.3, 0.4, 0.5, 0.2}, respectively. It is known from (41) that the measurement information on the fourth component of “x(t)” cannot be obtained by the first sink node, while the measurement information on the second component of “x(t)” cannot be obtained by the second sink node. Moreover, one has by (41) that rank([Γ, AΓ, A2 Γ, A3 Γ]) = 4(i = 1, 2) and rank([CiT , (Ci A)T , (Ci A2 )T , (Ci A3 )T ]T ) = 4, which means that the limits GKi and ΦKi in (31) exist. For this example, it is considered that each sink node has sufficient processing capabilities to compute the LOE x î (t), however, only two elements of x î (t) are allowed to be transmitted to the FC for satisfying the finite bandwidth, i.e., r1 = r2 = 2. In this case, one has ∆1 = ∆2 = 6, and it is concluded from (5) that  i i  H1 = diag{1, 1, 0, 0}, H2 = diag{1, 0, 1, 0} i H = diag{1, 0, 0, 1}, H4i = diag{0, 1, 1, 0} (42)  3i H5 = diag{0, 1, 0, 1}, H6i = diag{0, 0, 1, 1} with i = 1, 2. Then it follows from (9) and (42) that each RSE x âi (t) is determined by x âi (t) = Hi (t)ˆ xi (t), where Hi (t) = diag{σ1i (t) + σ2i (t) + σ3i (t), σ1i (t) + σ4i (t) + σ5i (t), σ2i (t) + σ4i (t) + σ6i (t), σ3i (t)+σ5i (t)+σ6i (t)}. Then the binary variables σ~i i (t)(~i = 1, · · · , 6) satisfy (7), and each stochastic process {σ~i i (t)} is i.i.d, ∑6 i i and Prob{σ~i (t) = 1} =π~i and ~i =1 π~i i = 1. When the bandwidth constraints, time-varying communication delays and packet dropouts in a CPS are modeled by (13), the delays and NSPLs can be identified by each packet’s time stamp, where the bounds in (12) are taken as dli = 1, dui = 2, τi = 2(i = 1, 2). When the LOE x ˆ1 (t) is transmitted to the FC, the dimensionality reduction strategy is determined by (42), and the selection probabilities are given by π21 = 0.61, π31 = 0.39, πj1 = 0 (j = 1, 4, 5, 6) and π12 = 0.15, π22 = 0.25, π32 = 0.17, π42 = 0.12, π52 = 0.18, π62 = 0.13. Then it can be calculated that ς∞,1 = 0.9941 < 1, ς2,1 = 0.5475 < 1 and ς∞,2 = 0.9645 < 1, ς2,2 = 0.4867 < 1, which means that the condition (28) holds. Therefore, according to Theorem 2, the MSE of the DFE for this example is bounded, and the corresponding fusion estimation performance is independent of the initial values. As pointed out in Section III, the cross-covariance matrix Σ12 (t) cannot be derived when the communication delays are time-varying, and thus the theory MSE of the DFE is difficult to be obtained. In this case, the Monte Carlo method (ch.1, [7]) is adopted to approach the theory value, then, by using Algorithm 1, the practical MSEs (PMSEs) over an average of 500 runs of Monte Carlo method for the DFE x ˆ(t) and the CSEs x ˆci (t)(i = 1, 2) are plotted in (a) of Fig.1, which shows that the estimation performance of the DFE is better than that of each local CSE. This is in line with the result “Tr(Po (t)) ≤ Tr(Σii (t))” in (23). Meanwhile, (a) of Fig.1 also shows that the MSE of the DFE is bounded, which accords with the result (30). Moreover, the relationships between the local covariance matrices Σii (t)(i = 1, 2) of the DFE and the local covariance matrices Σaii (t)(i = 1, 2) of the AODFE are shown in (b) of Fig.1. It is seen from this figure that the error between the Σii (t) and Σaii (t) converges to zero for t ≥ 50, which is in line with the result in Remark 4. Furthermore, it implies that the fusion estimation performance of the DFE for this example is independent of the initial values, and thus the weighting matrices Ωi (t)(i = 1, 2) can



a.

R EFERENCES

1 PMSE

PMSE

1

0

50

0.1 0.05 0

PMSE for the DFE

2

100

b. ||

(t)—

11

a (t)|| 11 2

Fusion Estimator in [16]

PMSE PMSE

x 10 4 2 0 100

||

(t)—

22

50

c.

0.5

150 a (t)|| 22 2

DFE

50

4

100

d.

DFE

Fusion Estimator in [22]

100

150

200

250

300

t/step

Fig. 1. (a): The estimation performance of the local CSEs and the DFE, where RMSEi represents the estimation performance of the ith CSE; (b): The trajectories of ||Σii (t) − Σa ii (t)||2 (i = 1, 2); (c): Comparison of the estimation performance for the DFE and the fusion estimation method in [16], where the packet dropping rate is 0.2; (d): Comparison of the estimation performance for the DFE and the fusion estimation method in [22].

be separately calculated without known the initial conditions of the sink nodes and the FC. To demonstrate the advantage of the proposed method, when we only consider the case of packet dropouts, the fusion estimation precision of the DFE, which is assessed by its PMSE, is compared with that of the centralized fusion estimator in [16]. Then, the PMSEs over an average of 500 runs of Monte Carlo method are plotted in (c) of Fig.1, which shows that the performance of the DFE outperforms the fusion estimator in [16]. This is because the estimation performance of the centralized fusion estimator in [16] is dependent on the previous measuring information, while the estimation performance of the DFE is only dependent on the current local estimation information at each time step. On the other hand, when considering the finite bandwidth, varying transmission delays and packet dropouts, the comparison between the DFE and the distributed fusion estimator in [22] is shown in (d) of Fig.1, where the fusion estimator [22] can only deal with the problem of the bandwidth constraints It is seen from this figure that the estimation performance of the DFE is better than that of the fusion estimator in [22], and the PMSE of the fusion estimator in [22] becomes larger as time goes. This implies that the communication delays and packet dropouts can degrade the estimation performance seriously if the above two factors are not considered for the fusion estimation of CPSs. VI. CONCLUSIONS A unified framework was presented in this paper for studying the distributed fusion estimation problem for a class of CPSs with communication constraints. A novel communication model with compensation strategy was proposed to establish the relationship among bandwidth constraints, delays and packet dropouts. Since the local estimation error cross-covariance matrices are difficult to be derived, a recursive DFE was designed by resorting to the covariance intersection fusion criterion. Stability condition, which is dependent on the selection probabilities and the upper bounds of delays and NSPLs, has been derived such that the MSE of the designed DFE is bounded, and it has been proved that the estimation performance of the DFE is independent of the initial values under certain conditions. Finally, a numerical example was given to demonstrate the feasibility and advantage of the proposed methods. Compared with the single-sensor state estimator, it is shown from Theorem 1 and the simulations that the multi-sensor fusion estimation approach can effectively improve the accuracy of state estimation that directly impacts the stability of control operation in the CPSs.

[1] K.H. Johansson, G.J. Pappas, P. Tabuada, C.J. Tomlin, Special issue on control of cyber-physical systems, IEEE Transactions on Automatic Control, vol. 59, no. 12, 2014, pp. 3120-3121. [2] Y.F. Huang, S. Werner, J. Huang, N. Kashyap, V. Gupta, State estimation in electric power girds: meeting new challenges presented by requirements of the future grid, IEEE Signal Processing Magazine, vol. 29, no. 5, 2012, pp. 33-43. [3] S. Deshmukh, Estimation and control in spatially distributed cyber physical systems, Ph.D dissertation, Kansas State University, 2013. [4] X. Cao, P. Cheng, J. Chen, S.S. Ge, Y. Cheng, Y. Sun, Cognitive radio based state estimation in cyber-physical systems, IEEE Journal on Selected Areas in Communicaiton, vol. 32, no. 3, 2014, pp. 489-502. [5] H. Li, L. Lai, V. Poor, Multicast routing for decentralized control of cyber physical systems with an application in smart grid, IEEE Journal on Selected Areas in Communication, vol. 30, 2012, pp. 1097-1107. [6] J.A. Roecker, C.D. McGillen, Comparison of two-sensor tracking methods based on state vector fusion and measurement fusion, IEEE Transactions on Aerospace and Electronic Systems, vol. 24, no. 4, 1988, pp. 447-449. [7] Y. Bar-Shalom, X.R. Li, T. Kirubarajan, Estimation with applications to tracking and navigation, John Wilely and Sons, Inc., 2001. [8] S.L. Sun, Z.L. Deng, Muti-sensor optimal information fusion Kalman filter, Automatica, vol. 40, 2004, pp. 1017-1023. [9] S.J. Julier, J.K. Uhlmann, General decentralized data fusion with covariance intersection (CI), Handbook of Multisensor Data Fusion, CRC Press, 2001, ch. 12, pp. 1-25. [10] Y. Wang, X.R. Li, Distributed estimation fusion with unavaliable crosscorrelation, IEEE Transactoins on Aerospace and Electronic Systems, vol. 48, no.1, 2012, pp. 259-278. [11] J. Sijs, M. Lazar, State fusion with unknown correlation: ellipsoidal intersection, Automatica, vol. 48, 2012, pp. 1874-1878. [12] Z.L. Deng, P. Zhang, W.J. Qi, J.F. Liu, Y. Gao, Sequential covariance intersection fusion Kalman filter, Informaiton Science, vol. 189, 2012, pp. 293-309. [13] B. Chen, G. Hu, W.A. Zhang, L. Yu, Information fusion estimation for spatially distributed cyber-physical systems with communication delay and bandwidth constraints, Proceedings of the American Control Conference, Chicago, USA, 2015, pp. 5152-5157. [14] P.Y. Chen, S.M. Cheng, K.C. Chen, Smart attack in smart grid communication networks, IEEE Communication Magazine, vol. 50, no. 8, 2012, pp. 24-29. [15] J.P. Hespana, P. Naghshtabrizi, Y. Xu, A survey of recent results in networked control systems, Proceedings of the IEEE, vol. 95, no. 1, 2007, pp. 138-162. [16] S. Deshmukh, B. Natarajan, A. Pahwa, State estimation over a lossy network in spatially distributed cyber-physical systems, IEEE Transactions on Signal Processing, vol. 62, no. 15, 2014, pp. 3911-3923. [17] A. Chiuso, L. Schenato, Information fusion strategies and performance bounds in packet-drop networks, Automatica, vol. 47, no. 7, 2011, pp. 1304-1316. [18] B. Chen, W.A. Zhang, G. Hu, L. Yu, Networked fusion kalman filtering with multiple uncertainties, IEEE Transactions on Aerospace and Electronic Systems, vol. 51, no. 3, 2015, pp. 2332-2349. [19] J. Ma, S. Sun, Information fusion estimation for systems with multiple sensors of different packet dropout rates, Information Fusion, vol. 12, 2011, pp. 213-222. [20] Y. Xia, J. Shang, J. Chen, G.P. Liu, Networked data fusion with packet losses and variable delays, IEEE Transactions on Systems, Man and Cybernetics–Part B: Cybernetics, vol. 39, no. 5, 2009, pp. 1107-1120. [21] B. Chen, W.A. Zhang, L. Yu, Distributed fusion estimation with missing measurements, random transmission delays and packet dropouts, IEEE Transactions on Automatic Control, vol. 59, no. 7, 2014, pp. 1961-1967. [22] B. Chen, W.A. Zhang, L. Yu, Distributed finite-horizon fusion Kalman filtering for bandwidth and energy constrained wireless sensor networks, IEEE Transactions on Signal Processing, vol.62, no.4, 2014, pp.797-812. [23] H. Ma, Y.H. Yang, Y. Chen, K.J. Ray Liu, Q. Wang, Distributed state estimation with dimension reduction preprocessing, IEEE Transactions on Signal Processing, vol. 62, no. 12, 2014, pp. 3098-3110. [24] B. Chen, W.A. Zhang, L. Yu, G. Hu, H. Song, Distributed fusion estimation with communication bandwidth constraints, IEEE Transactions on Automatic Control, vol.60, no. 5, 2015, pp. 1398-1403. [25] A. Riberio, I.D. Schizas, S.I. Roumeliotis, G.B. Glannkis, Kalman filtering in wireless sensor networks, IEEE Control System Magazine, vol. 30, 2010, pp. 66-96. [26] X. Shen, P. K.Varshney, Y. Zhu, Robust distributed maximum likelihood estimation with dependent quantized data, Automatica, vol. 50, no. 1, 2014, pp. 169-174. [27] B. Chen, G. Hu, W.A. Zhang, L. Yu, Distributed mixed H2 /H∞ fusion estimation with limited communication capacity, IEEE Transactions on Automatic Control, vol. 61, no. 3, 2016, pp. 805-810.
