- Email: [email protected]
Optimal Estimation Algorithm Design under Event-based Sensor Data Scheduling*
6th IFAC Workshop on Distributed Estimation and Control in 6th Networked Systemson 6th IFAC IFAC Workshop Workshop on Distributed Distributed Estimation Estimation and and Control Control in in 6th IFAC Workshop on Distributed Estimation and Control in Networked Systems Available online at www.sciencedirect.com September 8-9, 2016. Tokyo, Japan Networked Systems Networked Systems September 8-9, 2016. Tokyo, Japan September September 8-9, 8-9, 2016. 2016. Tokyo, Tokyo, Japan Japan
ScienceDirect
IFAC-PapersOnLine 49-22 (2016) 157–162
Optimal Estimation Algorithm Design Optimal Estimation Algorithm Design Optimal Algorithm Design underEstimation Event-based Sensor Data under Event-based Sensor Data ⋆ under Event-based Sensor Data Scheduling ⋆⋆ Scheduling Scheduling Lidong He ∗∗ Yifei Qi ∗∗ Jiming Chen ∗∗ Lidong He ∗∗ Yifei Qi ∗∗ Jiming Chen ∗∗ Lidong Lidong He He Yifei Yifei Qi Qi Jiming Jiming Chen Chen ∗ State Key Laboratory of Industrial Control Technology, College of ∗ ∗ Key Laboratory of Industrial Control Technology, College of ∗ State State Key of Control Technology, College Science and Engineering, Zhejiang University State Control Key Laboratory Laboratory of Industrial Industrial Control Technology, College of of Control Science and Engineering, Zhejiang University Control Science and Engineering, Zhejiang University (E-mail: [email protected], [email protected], Control Science and Engineering, Zhejiang University (E-mail: [email protected], [email protected], (E-mail: [email protected], [email protected]). (E-mail: [email protected], [email protected], [email protected], [email protected]). [email protected]). [email protected]). Abstract: This paper investigates a remote state estimation problem over a communication Abstract: This investigates a state over Abstract: Thisis paper paper investigates a remote remotesignificance state estimation estimation problem over aaa communication communication channel, which of theoretical and practical in the problem study of cyber-physical systems Abstract: This paper investigates a remote state estimation problem over communication channel, which is of theoretical and practical significance in the study of cyber-physical systems channel, which is of theoretical and practical significance in the study of cyber-physical systems (CPS). consideration of bothand bandwidth estimation an event-based channel,Inwhich is of theoretical practical limitation significanceand in the study ofquality, cyber-physical systems (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based schedule is proposed to decide whether to transmit the current measurement to the remote (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based schedule is proposed to decide whether to transmit the current measurement to the remote schedule is proposed to decide whether to transmit the current measurement to the estimator or not. To overcome the difficulties induced by the bilaterally hidden-truncation schedule is proposed to decide whether to transmit the current measurement to the remote remote estimator or not. To overcome the difficulties induced by the bilaterally hidden-truncation estimator or event-based not. To To overcome overcome the difficulties induced by the the bilaterally hidden-truncation process from scheme,the we difficulties introduce the generalized closed skew normal distribution estimator or not. induced by bilaterally hidden-truncation process from event-based scheme, we introduce the generalized closed skew normal distribution process from scheme, we the closed skew normal to accurately portray the state distribution. Furthermore, an exact MMSE algorithm process from event-based event-based scheme, we introduce introduce the generalized generalized closed skewestimation normal distribution distribution to accurately portray the state distribution. Furthermore, an exact MMSE estimation algorithm to accurately portray the state distribution. Furthermore, an exact MMSE estimation is designed, which does not rely on the Gaussian approximation assumption. Numerical examples to accurately portray the state distribution. Furthermore, an exact MMSE estimation algorithm algorithm is designed, which does not rely on the Gaussian approximation assumption. Numerical examples is designed, which does not rely on the Gaussian approximation assumption. Numerical are provided to show the effectiveness of the proposed algorithm. is designed, which does not rely on the Gaussian approximation assumption. Numerical examples examples are provided to show the effectiveness of the proposed algorithm. are provided provided to to show show the the effectiveness effectiveness of of the the proposed proposed algorithm. algorithm. are © 2016, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: CPS, remote state estimation, event-based scheduling, MMSE estimation Keywords: CPS, CPS, remote state state estimation, event-based event-based scheduling, MMSE MMSE estimation Keywords: Keywords: CPS, remote remote state estimation, estimation, event-based scheduling, scheduling, MMSE estimation estimation 1. INTRODUCTION assumed to be Gaussian, which benefits of obtaining the 1. INTRODUCTION INTRODUCTION assumed to be Gaussian, benefits of obtaining the 1. assumed Gaussian, which benefits obtaining the closed-form expression of which the estimation error covariance 1. INTRODUCTION assumed to to be be Gaussian, which benefits of of obtaining the closed-form expression of the estimation error covariance closed-form expression of the estimation error covariance the packet is kept saving due to its small innovation. With the integration of communication, computation and when closed-form expression of the estimation error covariance With the the integration of communication, communication, computation and when the the packet packet is is kept kept saving saving due due to to its its small small innovation. innovation. With integration of and control elements, cyber-physical systemscomputation (CPS) enable a when when packet is kept saving due to its small innovation. With the integration of communication, computation and In a the slightly different problem setting, Sukhavasi and control elements, cyber-physical systems (CPS) enable a control elements, cyber-physical systems (CPS) wealth new applications, including environment moni-aa In a slightly different problem setting, Sukhavasi and control of elements, cyber-physical systems (CPS) enable enable In a slightly different problem setting, Sukhavasi Hassibi (2013) have indicated that the hidden-truncation a slightly different problem setting, Sukhavasi and and wealth and of new new applications, including environment moni- In wealth applications, including environment monitoring control, smart home, unmanned aerial vehicle (2013) have indicated that the hidden-truncation wealth of of new applications, including environment moni- Hassibi Hassibi (2013) have indicated that the hidden-truncation process caused quantization is generally non-Gaussian. Hassibi (2013) by have indicated that the hidden-truncation toring and control, control, smart(Chen home, et unmanned aerial vehicle toring and smart home, unmanned aerial vehicle and modern battlefields al., 2011). Despite of process caused by quantization is generally non-Gaussian. toring and control, smart home, unmanned aerial vehicle Such caused by is non-Gaussian. a finding us the remark of ˚ Astr¨ om and Bernprocess caused recalls by quantization quantization is generally generally non-Gaussian. and exciting modern application battlefields (Chen (Chen et al., al., 2011). Despite of process and et Despite of ˚ the prospect, some2011). technical issues aa finding recalls us the remark of A str¨ omisand Bern˚ and modern modern battlefields battlefields (Chen et al., 2011). Despite of Such Such finding recalls us the remark of A str¨ Bernhardsson (2002): “An interesting conclusion that the ˚ a finding recalls us the remark of Astr¨oom m and and Bernthe exciting exciting application prospect, some technical issues Such the prospect, technical issues must be fullyapplication analyzed before CPS some become commonplace (2002): “An interesting conclusion is that the the exciting application prospect, some technical issues hardsson hardsson (2002): “An interesting conclusion is that steady state probability distribution of the control error (2002): “An interesting conclusion is that the the must be be fullyetanalyzed analyzed before CPSthem, become commonplace must before CPS commonplace (Johansson al., 2014). Among much effort has hardsson steady state probability distribution of the control error must be fully fully analyzed before CPS become become commonplace steady state probability distribution of the control error (for event-based sample) distribution is non-Gaussian even if theerror disstate probability of the control (Johansson et to al.,examining 2014). Among Among them, much estimation effort has has steady (Johansson et al., 2014). them, much effort been devoted how to improve (for sample) isInstead non-Gaussian even if the dis(Johansson et al., 2014). Among them, much effort has turbances (for event-based event-based sample) even if are Gaussian.” of Gaussian, event-based sample) is is non-Gaussian non-Gaussian even Sukhavasi if the the disdisbeen devoted devoted to to examining examining how to to improve estimation been how estimation performance, the consideration limited sensing and (for turbances are Gaussian.” Instead of Gaussian, Sukhavasi been devotedwith to examining how toofimprove improve estimation turbances are Gaussian.” Instead of Gaussian, Sukhavasi and Hassibi (2013) showed that such a process can be propturbances are Gaussian.” Instead of Gaussian, Sukhavasi performance, with the consideration of limited sensing and performance, with the consideration of limited sensing and communication resources, which motivates research (2013) that such aaClosed processSkew can be propperformance, with the consideration of limitedthe sensing and and and Hassibi Hassibi (2013) showed that properly portrayed via showed the Generalized Normal Hassibi (2013) showed that such such a process process can can be be propcommunication resources, which motivates the2012; research communication resources, which motivates the on the event-based scheduling (Heemels et al., Han and erly portrayed via the Generalized Closed Skew Normal communication resources, which motivates the research research erly portrayed via the Generalized Closed Skew Normal (GCSN) distribution (Genton, 2004).Closed Skew Normal erly portrayed via the Generalized on the event-based scheduling (Heemels et al., 2012; Han on the event-based scheduling (Heemels et al., 2012; Han et al., 2015; Cheng et al., 2016). distribution (Genton, 2004). on the event-based scheduling (Heemels et al., 2012; Han (GCSN) (GCSN) (GCSN) distribution distribution (Genton, (Genton, 2004). 2004). et al., al., 2015; 2015; Cheng Cheng et et al., al., 2016). et et al., 2015; Cheng et al., 2016). 2016).innovation plays the dom- In our view, there exist some inherent connections between In many studies, measurement In our view, there exist some inherent connections between In our view, there exist some inherent connections between theour quantization andinherent event-based scheduling, due view, thereproblem exist some connections between In many many studies, studies, measurement innovation schemes. plays the the domdomIn measurement innovation inating in the design of event-based es- In the quantization problem and event-based scheduling, due In manyrole studies, measurement innovation plays plays the In domthe quantization problem and event-based scheduling, due to the fact that both lead to the hidden-truncation process. the quantization problem and event-based scheduling, due inating role in the design of event-based schemes. In esinating in design In timation theory, innovationschemes. is defined as to the fact that both lead to the hidden-truncation process. inating role role in the themeasurement design of of event-based event-based schemes. In esesto the fact that both lead to the hidden-truncation process. Such a connection motivates us to recheck the Gaussian to the fact that both lead to the hidden-truncation process. timation theory, measurement innovation is defined as timation theory, measurement innovation is defined as the difference between the measurement collected via the motivates us recheck Gaussian timation theory, measurement innovation is defined as Such Such aaa connection connection motivates us to to recheck the the Gaussian approximation assumption in the event-based scheme and connection motivates us to recheck the Gaussian the difference difference between the the measurement collected via via the the Such the between the collected sensor of current time, and one-step measurement assumption in the event-based scheme and the difference between the measurement measurement collected viaprethe approximation approximation assumption in the event-based scheme and develop a novelassumption estimation in algorithm via GCSN distribuapproximation the event-based scheme and sensor of current time, and the one-step measurement presensor current the measurement prediction at theand remote estimator of the last step. develop aa novel estimation algorithm via GCSN distribusensor of ofcomputed current time, time, and the one-step one-step measurement pre- tion. develop novel estimation algorithm via GCSN distribuThe main contribution of this paper can be stated develop a novel estimation algorithm via GCSN distribudiction computed at the remote estimator of the last step. diction the remote of last step. Innovation revealsat the potential contribution of the The main contribution of this paper can be stated diction computed computed at the remote estimator estimator of the the lastnewly step. tion. tion. The contribution of paper can as follows. By introducing GCSN distribution into The main main contributionthe of this this paper can be be stated stated Innovation reveals the the potential contribution of the the newly Innovation reveals contribution of obtained measurement for minimizing the estimation error tion. as follows. By introducing the GCSN distribution into Innovation reveals the potential potential contribution of the newly newly as follows. By introducing the GCSN distribution into the event-based scheduling, an exact MMSE estimation as follows. By introducing the GCSN distribution into obtained measurement measurement for minimizing the estimation estimation error the event-based scheduling, an exact MMSE estimation obtained minimizing the error covariance. In a type offor widely used event-based schedule, obtained measurement for minimizing the estimation error the event-based scheduling, an exact MMSE estimation algorithm which scheduling, does not resort to theMMSE Gaussian approxthe event-based an exact estimation covariance. In a type of widely used event-based schedule, covariance. aa type widely used the packet isIn onlyof its innovation exceedsschedule, the pre- algorithm which doesis not resort to the covariance. Insent type ofwhen widely used event-based event-based schedule, algorithm which resort to Gaussian approxassumption proposed the Gaussian first timeapproxfor the algorithm which does does not not resort for to the the Gaussian approxthe packet packet is sent only only when when its its innovation exceeds the prepre- imation the sent exceeds the designated threshold; not since small innovation assumption is proposed for the first time for the the packet is is sent onlyOtherwise when its innovation innovation exceeds the pre- imation imation assumption is proposed for the first time multidimensional linear Gaussian Markov process. assumption is proposed for the first time for for the the designated threshold; Otherwise prediction not since since small innovation designated threshold; Otherwise innovation implies that the measurement at the remote imation multidimensional linear Gaussian Markov process. designated threshold; Otherwise not not since small small innovation multidimensional linear Gaussian Markov process. multidimensional linear Gaussian Markov process. implies that the measurement prediction at the remote n b implies that the measurement prediction at the remote estimator side has already been closed to the real data implies that the measurement prediction at the remote Notations: Rnn is the n-dimensional Euclidean space. yabb estimator side has You already been closed to theprior real state data Notations: the n-dimensional Euclidean space. yyaab estimator side has the R is the Euclidean space. (Wu et al., 2013; and been Xie, closed 2013). to The col{ya , . . . , R ybn} is the column vector with entries estimator side has already already been closed to the real real data data Notations: Notations: R isrepresents the n-dimensional n-dimensional Euclidean space. ya (Wu et al., 2013; You and Xie, 2013). The prior state ′ col{y , . . . , y } represents the column vector with entries (Wu et al., 2013; You and Xie, 2013). The prior state a b , . . . , y } represents the column vector with entries distribution scheduling is usually , . . . y x denotes the transpose of a vector or a matrix. a b (Wu et al., under 2013; innovation-based You and Xie, 2013). The prior state ycol{y a b col{y , . . . , y } represents the column vector with entries a b ′ distribution under under innovation-based innovation-based scheduling scheduling is is usually usually , . . . , yybb .. x′′ E[(x denotes the transpose or ′ a distribution a the of vector or a a matrix. matrix. a ,, .. .. .. ,, y) ] aisvector the covariance of − E[x])(y − E[y])of distribution innovation-based scheduling is usually yyyCov(x, ⋆ yb . x x denotes denotes the transpose transpose of a The work of under L. He was supported by Natural Science Foundation ′′ a vector or a matrix. ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ ′ ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ random variable (r.v.) x and y; Cov(x, x) is abbreviated as The work of L. He was supported by Natural Science Foundation of China (NSFC) under Grant No. 61503337 and China Postdoctoral The work of L. He was supported by Natural Science Foundation ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ The work of L. He was supported by Natural Science Foundation random variable (r.v.) x and y; Cov(x, x) is abbreviated as random variable (r.v.) x and y; Cov(x, x) is abbreviated as of China (NSFC) under Grant No. 61503337 and China Postdoctoral Var(x). Similarly define the conditional variance Var(x|y). Science Foundation under Grant No. 2015M571870; The work of J. of China (NSFC) under Grant No. 61503337 and China Postdoctoral random variable (r.v.) x and y; Cov(x, x) is abbreviated as of China (NSFC) under Grant No. 61503337 and China Postdoctoral Var(x). Similarly define the conditional Science Foundation Grant No. 2015M571870; The work of J. Var(x). Similarly define the variance Var(x|y). p(x) represents probability density variance functionVar(x|y). (pdf) of Chen was supportedunder by NSFC under No. U1401253. Science Foundation under Grant No. 2015M571870; The Var(x). Similarlythe define the conditional conditional variance Var(x|y). Science Foundation under Grant No. Grant 2015M571870; The work work of of J. J. p(x) represents the probability density function (pdf) of Chen was supported by NSFC under Grant No. U1401253. p(x) represents the probability density function Chen was supported by NSFC under Grant No. U1401253. p(x) represents the probability density function (pdf) (pdf) of of Chen was supported by NSFC under Grant No. U1401253. Copyright © 2016, 2016 IFAC 157Hosting by Elsevier Ltd. All rights reserved. 2405-8963 © IFAC (International Federation of Automatic Control) Copyright 2016 IFAC 157 Copyright © 2016 IFAC 157 Peer review© of International Federation of Automatic Copyright ©under 2016 responsibility IFAC 157Control. 10.1016/j.ifacol.2016.10.389
2016 IFAC NECSYS 158 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
Fig. 1. Event-based sensor data scheduling a r.v. P(a ≤ x ≤ b) denotes the probability with the limited range [a, b]. Nn (µ, Σ) is the n-dimensional Gaussian distribution with mean µ and variance Σ, corresponding pdf is denoted as φn (x; µ, Σ). Φn ([a, b]; µ, Σ) represents the integration of φn (x; µ, Σ) within the limited range [a, b], where, a, b ∈ Rn and the integration is component-wise. 2. PROBLEM SETUP We consider the sensor data scheduling for remote state estimation problem as depicted in Fig. 1. The discrete linear time-invariant process is described as follows: xk+1 = Axk + ωk , (1) (2) yk = Cxk + vk , where, xk ∈ Rn is the system state, yk ∈ Rm is the measurement collected via local sensor at time k. Gaussian white noises ωk and vk are assumed to be zero-mean with variances Q ≥ 0 and R > 0, respectively. The initial state x0 is also Gaussian with distribution Nn (µ0 , Σ0 ). Two sequences {ωk } and {vk } are independent, which implies that E[ωi vj′ ] = 0, (∀i, j). The pair (A, C) is detectable and √ (A, Q) is stabilizable. To save scarce communication resource which may be used by other applications, the scheduler equipped at the sensor side decides whether to transmit the current measurement yk over the network to the remote estimator or not at each time step. Denote γk = 1 as the packet is sent and γk = 0 otherwise, which will be shortly specified in Eq.(5). To concentrate on analyzing the senor scheduling scheme, we ignore other imperfect effects such as packet loss, time delay and quantization, etc. Denote the information set available at the estimator side of time k ∈ Z+ as Ik {γ1 y1 , . . . , γk yk }. (3) Assume that there is a feedback channel from the remote estimator to the local scheduler for reliably transmitting the measurement predict yˆk− E[yk |Ik−1 ] at any time instant k ∈ Z+ , where, I0 {x0 }. 1 On that basis, the scheduler can compute the measurement innovation zk specified as Eq. (4) for indicating how accurately the remote estimator can predict the current measurement yk before it sends, and thus how to transmit yk may help to reduce such a prediction error covariance. The measurement innovation is defined as zk yk − yˆk− = yk − C x ˆ− (4) k. We have E[zk ] = 0 due to projection law (Kailath et al., 2000). The scheduling sequence {γk } is specified via the following event-based scheme for a given threshold δ ∈ Rm : 0, if −δ ≤ zk ≤ δ, (5) γk = 1, otherwise. 1 Such a feedback is based on the IEEE 802.15.4/ZigBee protocol. Details can reference Wu et al. (2013).
158
Here, the inequality holds for the component-wise sense. In such a scheme, saving one step of transmission (γk = 0) does not mean the corresponding error covariance increases as that of open-loop prediction. On the contrary, such an event implies that the measurement predict has already been quite reliable, which helps to reduce the communication cost while keeping the estimation error covariance under a pre-designed level. The primary concern of this paper is to propose an exact MMSE estimation algorithm for system (1) ∼ (2) under scheduling scheme (5). We assume that no a prior probability information of state predict is available. 3. PRELIMINARIES This paper chooses the Bayesian method as depicted in Ho and Lee (1964) to obtain the MMSE estimation algorithm. To ease analysis, the equalities frequently used in subsequent sections are provided in this section. 1) Compute the prior pdf via Chapman-Kolmogorov equation p(xk+1 |Ik ) = p(xk+1 |xk )p(xk |Ik )dxk . (6)
2) Obtain the posterior pdf via Bayesian rule on account of the event really happens: • In the case of γk+1 = 1, i.e., Ik+1 = Ik ∪ {yk+1 }. p(xk+1 |Ik+1 ) =
p(yk+1 |Ik , xk+1 )p(xk+1 |Ik ) . p(yk+1 |Ik ) (7)
where, p(yk+1 |Ik ) = p(yk+1 |Ik , xk+1 )p(xk+1 |Ik )dxk+1 .
(8)
• In the case of γk+1 = 0, i.e., Ik+1 = Ik ∪ {−δ ≤ zk+1 ≤ δ}. p(xk+1 |Ik+1 ) P(−δ ≤ zk+1 ≤ δ|Ik , xk+1 ) . (9) =p(xk+1 |Ik ) P(−δ ≤ zk+1 ≤ δ|Ik )
Different from classic Kalman filter, the integrations in (6) and (8) are computationally prohibitive, since γk cannot be expressed as the linear combination of the elements of Ik , which clearly reveals the inherent nonlinear property of the event-based mechanism. For facilitate of analysis, Wu et al. (2013); You and Xie (2013) assumed that the prior state distribution is Gaussian, i.e., − ˆ− p(xk+1 |Ik ) = φ(xk+1 ; x k+1 , Pk+1 ),
(10)
− ′ ˆ− where, xˆ− k+1 E[xk+1 |Ik ], Pk+1 E[(xk+1 − x k+1 )(·) ]. Such an assumption substantially simplifies the computation and approximates the corresponding filter to be linear. Though nice estimation algorithm is obtained, the rationality of such a Gaussianality approximation lacks rigorous verification.
In next two sections, the Gaussian approximation assumption is shown to break down and an exact MMSE estimation algorithm, which does not resort to any approximation, is also provided.
2016 IFAC NECSYS September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
4. GCSN DISTRIBUTION FOR GAUSSIAN HIDDEN-TRUNCATION PROCESS Suppose {xk } and {yk } are jointly Gaussian. For notational convenience, we define three coefficients as follows: Kk Cov(xk , y1k ), Lk Var(y1k ), Mk Var(xk ), which implies that Kk ∈ Rn×(km) , Lk ∈ R(km)×(km) and Mk ∈ Rn×n . We can directly obtain ′ Var(xk |y1k ) = Mk − Kk L−1 k Kk ,
Var(y1k |xk ) = Lk − Kk′ Mk−1 Kk . With a slight confusion of notation, we use Ik as sets when referring to information pattern, and as vectors when referring to the information sets. Then we have the following result. The proof uses some properties of moment generating function of r.v. and is skipped due to space limitation. Lemma 1. The state xk conditioned on the information set Ik can be decomposed into two independent terms: k (11) xk |Ik = ξk + Kk L−1 k y1 |Ik , where, ξk ∼ Nn (0, Var(xk |y1k )). The corresponding estimation error covariance is given by −1 ′ k Var(xk |Ik ) = Var(xk |y1k ) + Kk L−1 k Var(y1 |Ik )Lk Kk . (12)
It is shown in Eq. (12) that the posterior pdf under eventbased scheme is generally non-Gaussian. The Gaussian approximation assumption is valid only when Var(y1k |Ik ) approaches 0, which implies δ → 0 and the measurements are sent at any time instant. Under this extreme situation, k k xk |Ik → xk |y1k ∼ Nn (Kk L−1 (13) k y1 , Var(xk |y1 )), which reduces to the classic MMSE estimator (Kailath et al. (2000), pp.81). The last term of Eq. (12) explicitly quantifies the performance deterioration incurred from reducing the communication cost. As expected, the problem becomes substantially challenging when dropping the Gaussian approximation assumption and the techniques employed in the aforementioned works (Wu et al. (2013); You and Xie (2013)) no longer apply. To accurately characterize the Gaussian hidden-truncation process from the event-trigger scheme (5), we develop a new framework by introducing the closed skew normal (CSN) distribution. Definition 2. (Genton (2004), pp.26). Consider m ≥ 1, n ≥ 1, µ ∈ Rn , ν ∈ Rm , D is an arbitrary m × n matrix, Σ and ∆ are positive definite matrices of dimension n×n and m × m, respectively. Then the pdf of a CSN distribution is given by CSNn,m (x; µ, Σ, D, ν, ∆) Φm (−∞, D(x − µ); ν, ∆) . (14) =φn (x; µ, Σ) Φm (−∞, 0; ν, ∆ + DΣD′ ) Here, µ, Σ, D can be deemed as location, scale and skewness parameters, respectively. When taking D = 0, Eq.(14) reduces to the pdf of the Gaussian distribution with mean µ and variance Σ. It should be noted that for the general CSN distribution, µ and Σ cannot be deemed as mean and variance, which are in fact involved to compute and can be found in Genton (2004). Intuitively, for independent r.v. X ∼ Nn (µ, Σ) and Z ∼ Nm (ν, ∆), Eq. (14) explicitly characterizes the pdf of the conditional event (Genton (2004), pp.28) 159
159
X|{Z − D(X − µ) < 0}, where, D is the weighted matrix. For the bilateral truncation scheme (5), we need the modified version of Eq. (14), which is called the generalized closed skew normal (GCSN) distribution. Definition 3. (Sukhavasi and Hassibi (2013)). The pdf of the GCSN distribution is given by: GCSNn,m (x; µ, Σ, D, Sm , ∆) Φm (Sm ; D(x − µ), ∆) . (15) =φn (x; µ, Σ) Φm (Sm ; 0, ∆ + DΣD′ ) Here, Sm is the region determined via the specified realization. With such a tool for characterizing the bilateral truncationprocess, an exact MMSE estimation algorithm for system under event-based scheme is proposed in next section. 5. MMSE ESTIMATION ALGORITHM UNDER EVENT-BASED SCHEME WITH GCSN DISTRIBUTION Before presenting the main result of the paper, this section starts with some lemmas. Specifically, Lemma 4 and Lemma 5 play the dominating role in computing the multiplication of two Gaussian pdfs. Lemma 4 was put forward in Ho and Lee (1964), corresponding proof uses the mathematical trick known as “completion of squares” for multivariate argument. Lemma 5 uses Schur complement lemma (Kailath et al. (2000), pp.725). Most of the proofs are omitted in consideration of limited space. Lemma 4. Given x, a ∈ Rn , y ∈ Rm and symmetric matrices A ∈ Rn×n , B ∈ Rm×m , which are positively defined such that CAC ′ + B is non-singular (C ∈ Rm×n ). Then (x − a)′ A−1 (x − a) + (y − Cx)′ B −1 (y − Cx) =(x − µ)′ Σ−1 (x − µ) + (y − Ca)′ (CAC ′ + B)−1 (y − Ca), where, µ a + AC ′ (CAC ′ + B)−1 (y − Ca), Σ (A−1 + C ′ B −1 C)−1 .
Lemma 5. For X ∈ Rn×n , C ∈ Rm×n and R ∈ Rm×m , the following equality holds true: |Y ||CXC ′ + R| = |X||R|, (16)
where, Y (X −1 + C ′ R−1 C)−1 .
As will be shown in the proof of Theorem 8, Lemma 6 plays the dominating role in computing P(−δ ≤ zk+1 ≤ δ|Ik ) of Eq. (9). Proposition 7 is needed in calculating the integration of Eq. (6) and Eq. (8). Lemma 6. Let a, b ∈ Rm , p, q ∈ Rn and B ∈ Rm×n . If V ∼ Nn (µ1 , Σ), we have q Φm ([a + BV, b + BV ]; µ2 , Ω)φn (V ; µ1 , Σ)dV p a + Bµ1 − µ2 b + Bµ1 − µ2 =Φm+n , ; 0, ΣW , p − µ1 q − µ1 where,
ΣW
BΣB ′ + Ω −ΣB ′
− BΣ . Σ
(17)
2016 IFAC NECSYS 160 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
On the basis of Lemma 6, the following proposition can be easily obtained. Proposition 7. Let a, b ∈ Rm and B ∈ Rm×n . If V ∼ Nn (µ1 , Σ), we have EV [Φm ([a + BV, b + BV ]; µ2 , Ω)] =Φm ([a + Bµ1 − µ2 , b + Bµ1 − µ2 ]; 0, BΣB ′ + Ω). Now we shall tie together the ideas of the previous lemmas to state the main result of this paper. Theorem 8. The pdfs of the MMSE estimation algorithm for system (1) ∼ (2) under scheme (5) is given as follows: p(xk+1 |Ik ) − − − =GCSN (xk+1 ; µ− k+1 , Σk+1 , Dk+1 , [Ak , Bk ], ∆k+1 ), p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , [Ak+1 , Bk+1 ], ∆k+1 ). Corresponding coefficients are updated as: 1) Predict µ− k+1 = Aµk , ′ Σ− k+1 = AΣk A + Q, − −1 = Dk Σk A′ (Σ− , Dk+1 k+1 )
= ∆k +
Dk Σk Dk′
2) Correction • If γk+1 = 1 ′ Kk+1 = Σ− k+1 C , Lk+1 =
−
′ CΣ− k+1 C
We begin by verifying the recursion of the GCSN distributions, which does not resort to the event-based scheme. Consider a stable system with A = 0.9, C = 1, Q = 0.5, R = 0.5 and the initial state x0 ∼ N (5, 1). Let the packet of the first step fail to send and the remaining transmissions all be successful. Fig. 3 shows the correctness of the algorithm (18) ∼ (20). µ=0,Σ=5 µ=0,Σ=5,a=−0.5,b=10,D=5,∆=5 µ=0,Σ=5,a=−50,b=500,D=50,∆=0.1
0.45
+ R,
−1 ˜k+1 , µk+1 = µ− k+1 + Kk+1 Lk+1 y
Σk+1 =
This section leverages some numerical examples to show the effectiveness of the proposed estimation algorithm. Fig. 2 gives the pdf of one Gaussian distribution and two GCSN distributions. The blue curve shows that µ is no longer the expectation of a GCSN distribution. The red curve shows the truncation effect.
0.5
y˜k+1 = yk+1 − Cµ− k+1 , Σ− k+1 − − Dk+1 ,
(18)
− − ′ Dk+1 Σ− k+1 (Dk+1 ) .
6. SIMULATION
0.4 0.35
(19)
′ Kk+1 L−1 k+1 (Kk+1 ) , − ∆k+1 = ∆k+1 .
0.3 pdf
∆− k+1
Remark 10. We notice from Eq. (21) ∼ (22) that the dimensions of some coefficients linearly increase when the packets are not sent due to the event-based scheme. Such a problem arises from the bilateral hidden-truncation process (5) and can be explained via Lemma 6, where state-augmentation approach is adopted for dealing with the integration of limited range. One possible way to deal with this problem is to transform mechanism (5) in a way that turns ΣW of Eq. (17) into a diagonal matrix, which will be left as our future work.
Dk+1 = The region is updated as Ak+1 = Ak + Dk+1 (µ− k+1 − µk+1 ),
0.25 0.2 0.15 0.1
(20) Bk+1 = Bk + Dk+1 (µ− k+1 − µk+1 ). • If γk+1 = 0 µk+1 = µ− , , Σk+1 = Σ− k+1− k+1 0 Dk+1 ∆− k+1 , (21) Dk+1 = , ∆k+1 = 0 R C
0.05 0 −6
160
−2
0
2
4
6
8
x
Fig. 2. PDF of a GCSN distribution
x− − µ− ), Θk = C(ˆ k+1 k+1 Ak Bk , Bk+1 = . (22) Ak+1 = Θk − δ Θk + δ
8 Real State x
k
Estimated State
6
4
2 State
Both (20) and (22) are updated with the initial conditions Ak0 −δ, Bk0 δ, (23) where, k0 denotes the first time instant when the packet is not sent. Remark 9. To the best of our knowledge, Theorem 8 first proposes an exact MMSE estimation algorithm for the multidimensional linear Gaussian Markov process with the threshold-type event-based scheme. The introduced GCSN distribution avoids the Gaussian approximation assumption widely used in the event-based schedules, i.e., Wu et al. (2013); You and Xie (2013), which is proved to be inaccurate in Lemma 1.
−4
0
−2
−4
−6
0
20
40
60
80
100
Time
Fig. 3. State versus estimate under GCSN distribution
2016 IFAC NECSYS September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
As for the event-based schedule, consider an unstable system with A = 1.2, C = 1, Q = 0.25, R = 0.25 and x0 ∼ N (1, 1). The threshold for the event-based scheme (5) is δ = 0.2. The comparison between the real state and the corresponding estimate is given in Fig. 4. Fig. 5 shows the packet transmission sequence under event-based scheme. Here, 1 denotes the packet is not sent, which is different from Eq. (5). The estimation error and the trace of error covariance are provided in Fig. 6 and Fig. 7, respectively. The effectiveness of the proposed MMSE algorithm under event-based scheme via GCSN distribution is thus verified. 7. CONCLUSION This paper comes up with an event-based scheme for remote state estimation over a bandwidth-limited channel. Our main contribution lies in the design of an exact MMSE estimation algorithm that accurately exploits the additional information of the event and does not resort to the Gaussian approximation assumption. There are many interesting extension directions to explore for this work: finding more appropriate event-based scheme to reduce computational complexity; establishing the quantitative relationship among estimation performance, communication rate and event-triggered threshold; comparing system performance of the present algorithm with those in Wu et al. (2013) and Han et al. (2015), etc.
161
REFERENCES ˚ Astr¨om, K. and Bernhardsson, B. (2002). Comparison of Riemann and Lebesque sampling for first order stochastic systems. In Proceedings of the 41st IEEE Conference on Decision and Control, volume 2, 2011–2016. IEEE. Chen, J., Yu, Q., Cheng, P., Sun, Y., Fan, Y., and Shen, X. (2011). Game theoretical approach for channel allocation in wireless sensor and actuator networks. IEEE Transactions on Automatic Control, 56(10), 2332– 2344. Cheng, P., Qi, Y., Xin, K., Chen, J., and Xie, L. (2016). Energy-efficient data forwarding for state estimation in multi-hop wireless sensor networks. IEEE Transactions on Automatic Control, 61(5), 1322–1327. Genton, M. (2004). Skew-elliptical Distributions and Their Applications: A Journey Beyond Normality. CRC Press. Han, D., Mo, Y., Wu, J., Sinopoli, B., and Shi, L. (2015). Stochastic event-triggered sensor schedule for remote state estimation. IEEE Transactions on Automatic Control, 60(10), 2661–2675. Heemels, W., Johansson, K., and Tabuada, P. (2012). An introduction to event-triggered and self-triggered control. In Proceddings of the IEEE 51st Annual Conference on Decision and Control, 3270–3285. IEEE. Ho, Y. and Lee, R. (1964). A Bayesian approach to problems in stochastic estimation and control. IEEE Transactions on Automatic Control, 9(4), 333–339.
20000 Real State xk
1.5
Estimated State 15000
1
0.5
State
Estimation error
10000
5000
0
−0.5 0 −1 −5000
0
10
20
30 Time
40
50
60
−1.5
Fig. 4. State versus estimate under event-based schedule
0
10
20
30 Time
40
50
60
30 Time
40
50
60
Fig. 6. Estimation error
1 0.22
0.9
Trace of the estimation error covariance
0.8
Packet dropout
0.7 0.6 0.5 0.4 0.3 0.2 0.1 0
0
10
20
30 Time
40
50
0.21
0.2
0.19
0.18
0.17
60 0.16
Fig. 5. Scheduling sequence: 1 represents the packet is not sent and 0 represents successful transmission 161
0
10
20
Fig. 7. Trace of the estimation error covariance
2016 IFAC NECSYS 162 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
Johansson, K., Pappas, G., Tabuada, P., and Tomlin, C. (2014). Guest editorial special issue on control of cyberphysical systems. IEEE Transaction on Automatic Control, 59(12), 3120–3121. Kailath, T., Sayed, A., and Hassibi, B. (2000). Linear Estimation. Prentice Hall. Sukhavasi, R. and Hassibi, B. (2013). The Kalman-like particle filter: Optimal estimation with quantized innovations/measurements. IEEE Transactions on Signal Processing, 61(1), 131–136. Wu, J., Jia, Q., Johansson, K., and Shi, L. (2013). Eventbased sensor data scheduling: Trade-off between communication rate and estimation quality. IEEE Transactions on Automatic Control, 58(4), 1041–1046. You, K. and Xie, L. (2013). Kalman filtering with scheduled measurements. IEEE Transactions on Signal Processing, 61(6), 1520–1530. Appendix A. PROOF OF THEOREM 8 Proof. The proof is conducted via mathematical induction. Here, dimensions are omitted for sake of brevity. 1) Prediction of step 1 Using C-K equation (6) for system (1) with initial − distribution p(x0 |I0 ) = φ(x0 ; µ0 , Σ0 ), we have µ− 1 , Σ1 as in Eq. (18) when taking k = 0. 2) Correction of step 1 • If the packet is sent, i.e., I1 = I0 ∪ {γ1 = 1} From Bayesian rule, we have p(x1 |I1 ) = φ(x1 ; µ1 , Σ1 ), where, µ1 , Σ1 are given as Eq. (19). • If the packet is not sent, i.e., I1 = I0 ∪ {γ1 = 0} Eq. (9) can be computed as p(x1 |I1 ) =
− − Φ([A1 , B1 ]; C(x1 − µ1 ), R) . φ(x1 ; µ− 1 , Σ1 ) − ′ Φ([A1 , B1 ]; 0, CΣ1 C + R)
where, A1 , B1 are as in Eq. (23). When taking − µ1 = µ − 1 , Σ1 = Σ1 , D1 = C, ∆1 = R. we have p(x1 |I1 ) = GCSN (x1 ; µ1 , Σ1 , D1 , [A1 , B1 ], ∆1 ). 3) Prediction of step k+1 Take Sk [Ak , Bk ]. Without loss of generality, assume p(xk |Ik ) = GCSN (xk ; µk , Σk , Dk , Sk , ∆k ). From Eq. (6), we have ∞ p(xk+1 |Ik ) = φ(xk+1 ; Axk , Q)φ(xk ; µk , Σk ) −∞
Φ(Sk ; Dk (xk − µk ), ∆k ) dxk . Φ(Sk ; 0, ∆k + Dk Σk Dk′ ) Using both Lemma 4 and Lemma 5, as well as pdf of the Gaussian distribution, φ(xk+1 ; Axk , Q)φ(xk ; µk , Σk ) 12 ♯ |Σk ||Σ− k+1 | − = φ(xk ; µ♯k , Σ♯k )φ(xk+1 ; µ− k+1 , Σk+1 ), |Σk ||Q| ×
− where, µ− k+1 , Σk+1 are given as Eq. (18), and −1 (xk+1 − µ− µ♯k µk + Σk A′ (Σ− k+1 ) k+1 ),
Σ♯k [(Σk )−1 + A′ Q−1 A]−1 . In view of both Lemma 5 and Proposition 7, as well − as taking Dk+1 and ∆− k+1 as Eq. (18), we have 162
p(xk+1 |Ik ) − − − =GCSN (xk+1 ; µ− k+1 , Σk+1 , Dk+1 , Sk , ∆k+1 ). 4) Correction of step k+1 • If Ik+1 = Ik ∪ {γk+1 = 1} Computing p(yk+1 |Ik ) can conduct the similar way as computing p(xk+1 |Ik ). Further using the Gaussian property of p(yk+1 |xk+1 ), yielding p(xk+1 |Ik+1 ) =φ(xk+1 ; µk+1 , Σk+1 ) − − Φ(Sk ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) × . − − ′ Φ(Sk+1 ; 0, ∆− k+1 + Dk+1 Σk+1 (Dk+1 ) ) where, µk+1 and Pk+1 are given in Eq. (19), Ak+1 and Bk+1 are as in Eq. (20). Further take Dk+1 , ∆k+1 as in Eq. (19), we have p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , Sk+1 , ∆k+1 ). • If Ik+1 = Ik ∪ {γk+1 = 0} The crucial step in Eq. (9) is to compute P(−δ ≤ zk+1 ≤ δ|Ik ) Ωk+1 = . (A.1) − − − ′ Φ(Sk ; 0, ∆k+1 + Dk+1 Σ− k+1 (Dk+1 ) ) where, Ωk+1 yˆ− +δ k+1 − ′ φ(yk+1 ; Cµ− k+1 , CΣk+1 C + R) − −δ yˆk+1
× Φ([Mk+1 , Nk+1 ]; 0, Λk+1 )dyk+1 ,
− − Mk+1 Ak − Dk+1 Kk+1 L−1 k+1 (yk+1 − Cµk+1 ), − − Nk+1 Bk − Dk+1 Kk+1 L−1 k+1 (yk+1 − Cµk+1 ),
− − ′ Λk+1 ∆− k+1 + Dk+1 Σk+1 (Dk+1 ) . By virtue of Lemma 6, we have Ωk+1 = Φ([Ak+1 , Bk+1 ]; 0, ΣW ), where, Ak+1 , Bk+1 are given as Eq. (22), and − − − − ′ Dk+1 Σ− k+1 (Dk+1 ) + ∆k+1 Dk+1 Kk+1 . ΣW = − ′ Kk+1 (Dk+1 )′ Lk+1 Augment the product of the following two cdf− − s Φ(Sk ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) as well as − − Φ([ˆ yk+1 − δ, yˆk+1 + δ]; Cxk+1 , R) into another Gaussian cdf with dimension m+n, then combine it with Eq. (A.1), yielding p(xk+1 |Ik+1 ) − =φ(xk+1 ; µ− k+1 , Σk+1 )
Φ(Sk+1 ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) , Φ(Sk+1 ; 0, ΣW ) where, Dk+1 , ∆k+1 are given as Eq. (21). When − taking µk+1 = µ− k+1 , Σk+1 = Σk+1 , we have ×
ΣW = ∆k+1 + Dk+1 Σk+1 (Dk+1 )′ . In a word, p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , Sk+1 , ∆k+1 ). The induction step is thus complete.
ScienceDirect
IFAC-PapersOnLine 49-22 (2016) 157–162
Optimal Estimation Algorithm Design Optimal Estimation Algorithm Design Optimal Algorithm Design underEstimation Event-based Sensor Data under Event-based Sensor Data ⋆ under Event-based Sensor Data Scheduling ⋆⋆ Scheduling Scheduling Lidong He ∗∗ Yifei Qi ∗∗ Jiming Chen ∗∗ Lidong He ∗∗ Yifei Qi ∗∗ Jiming Chen ∗∗ Lidong Lidong He He Yifei Yifei Qi Qi Jiming Jiming Chen Chen ∗ State Key Laboratory of Industrial Control Technology, College of ∗ ∗ Key Laboratory of Industrial Control Technology, College of ∗ State State Key of Control Technology, College Science and Engineering, Zhejiang University State Control Key Laboratory Laboratory of Industrial Industrial Control Technology, College of of Control Science and Engineering, Zhejiang University Control Science and Engineering, Zhejiang University (E-mail: [email protected], [email protected], Control Science and Engineering, Zhejiang University (E-mail: [email protected], [email protected], (E-mail: [email protected], [email protected]). (E-mail: [email protected], [email protected], [email protected], [email protected]). [email protected]). [email protected]). Abstract: This paper investigates a remote state estimation problem over a communication Abstract: This investigates a state over Abstract: Thisis paper paper investigates a remote remotesignificance state estimation estimation problem over aaa communication communication channel, which of theoretical and practical in the problem study of cyber-physical systems Abstract: This paper investigates a remote state estimation problem over communication channel, which is of theoretical and practical significance in the study of cyber-physical systems channel, which is of theoretical and practical significance in the study of cyber-physical systems (CPS). consideration of bothand bandwidth estimation an event-based channel,Inwhich is of theoretical practical limitation significanceand in the study ofquality, cyber-physical systems (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based schedule is proposed to decide whether to transmit the current measurement to the remote (CPS). In consideration of both bandwidth limitation and estimation quality, an event-based schedule is proposed to decide whether to transmit the current measurement to the remote schedule is proposed to decide whether to transmit the current measurement to the estimator or not. To overcome the difficulties induced by the bilaterally hidden-truncation schedule is proposed to decide whether to transmit the current measurement to the remote remote estimator or not. To overcome the difficulties induced by the bilaterally hidden-truncation estimator or event-based not. To To overcome overcome the difficulties induced by the the bilaterally hidden-truncation process from scheme,the we difficulties introduce the generalized closed skew normal distribution estimator or not. induced by bilaterally hidden-truncation process from event-based scheme, we introduce the generalized closed skew normal distribution process from scheme, we the closed skew normal to accurately portray the state distribution. Furthermore, an exact MMSE algorithm process from event-based event-based scheme, we introduce introduce the generalized generalized closed skewestimation normal distribution distribution to accurately portray the state distribution. Furthermore, an exact MMSE estimation algorithm to accurately portray the state distribution. Furthermore, an exact MMSE estimation is designed, which does not rely on the Gaussian approximation assumption. Numerical examples to accurately portray the state distribution. Furthermore, an exact MMSE estimation algorithm algorithm is designed, which does not rely on the Gaussian approximation assumption. Numerical examples is designed, which does not rely on the Gaussian approximation assumption. Numerical are provided to show the effectiveness of the proposed algorithm. is designed, which does not rely on the Gaussian approximation assumption. Numerical examples examples are provided to show the effectiveness of the proposed algorithm. are provided provided to to show show the the effectiveness effectiveness of of the the proposed proposed algorithm. algorithm. are © 2016, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: CPS, remote state estimation, event-based scheduling, MMSE estimation Keywords: CPS, CPS, remote state state estimation, event-based event-based scheduling, MMSE MMSE estimation Keywords: Keywords: CPS, remote remote state estimation, estimation, event-based scheduling, scheduling, MMSE estimation estimation 1. INTRODUCTION assumed to be Gaussian, which benefits of obtaining the 1. INTRODUCTION INTRODUCTION assumed to be Gaussian, benefits of obtaining the 1. assumed Gaussian, which benefits obtaining the closed-form expression of which the estimation error covariance 1. INTRODUCTION assumed to to be be Gaussian, which benefits of of obtaining the closed-form expression of the estimation error covariance closed-form expression of the estimation error covariance the packet is kept saving due to its small innovation. With the integration of communication, computation and when closed-form expression of the estimation error covariance With the the integration of communication, communication, computation and when the the packet packet is is kept kept saving saving due due to to its its small small innovation. innovation. With integration of and control elements, cyber-physical systemscomputation (CPS) enable a when when packet is kept saving due to its small innovation. With the integration of communication, computation and In a the slightly different problem setting, Sukhavasi and control elements, cyber-physical systems (CPS) enable a control elements, cyber-physical systems (CPS) wealth new applications, including environment moni-aa In a slightly different problem setting, Sukhavasi and control of elements, cyber-physical systems (CPS) enable enable In a slightly different problem setting, Sukhavasi Hassibi (2013) have indicated that the hidden-truncation a slightly different problem setting, Sukhavasi and and wealth and of new new applications, including environment moni- In wealth applications, including environment monitoring control, smart home, unmanned aerial vehicle (2013) have indicated that the hidden-truncation wealth of of new applications, including environment moni- Hassibi Hassibi (2013) have indicated that the hidden-truncation process caused quantization is generally non-Gaussian. Hassibi (2013) by have indicated that the hidden-truncation toring and control, control, smart(Chen home, et unmanned aerial vehicle toring and smart home, unmanned aerial vehicle and modern battlefields al., 2011). Despite of process caused by quantization is generally non-Gaussian. toring and control, smart home, unmanned aerial vehicle Such caused by is non-Gaussian. a finding us the remark of ˚ Astr¨ om and Bernprocess caused recalls by quantization quantization is generally generally non-Gaussian. and exciting modern application battlefields (Chen (Chen et al., al., 2011). Despite of process and et Despite of ˚ the prospect, some2011). technical issues aa finding recalls us the remark of A str¨ omisand Bern˚ and modern modern battlefields battlefields (Chen et al., 2011). Despite of Such Such finding recalls us the remark of A str¨ Bernhardsson (2002): “An interesting conclusion that the ˚ a finding recalls us the remark of Astr¨oom m and and Bernthe exciting exciting application prospect, some technical issues Such the prospect, technical issues must be fullyapplication analyzed before CPS some become commonplace (2002): “An interesting conclusion is that the the exciting application prospect, some technical issues hardsson hardsson (2002): “An interesting conclusion is that steady state probability distribution of the control error (2002): “An interesting conclusion is that the the must be be fullyetanalyzed analyzed before CPSthem, become commonplace must before CPS commonplace (Johansson al., 2014). Among much effort has hardsson steady state probability distribution of the control error must be fully fully analyzed before CPS become become commonplace steady state probability distribution of the control error (for event-based sample) distribution is non-Gaussian even if theerror disstate probability of the control (Johansson et to al.,examining 2014). Among Among them, much estimation effort has has steady (Johansson et al., 2014). them, much effort been devoted how to improve (for sample) isInstead non-Gaussian even if the dis(Johansson et al., 2014). Among them, much effort has turbances (for event-based event-based sample) even if are Gaussian.” of Gaussian, event-based sample) is is non-Gaussian non-Gaussian even Sukhavasi if the the disdisbeen devoted devoted to to examining examining how to to improve estimation been how estimation performance, the consideration limited sensing and (for turbances are Gaussian.” Instead of Gaussian, Sukhavasi been devotedwith to examining how toofimprove improve estimation turbances are Gaussian.” Instead of Gaussian, Sukhavasi and Hassibi (2013) showed that such a process can be propturbances are Gaussian.” Instead of Gaussian, Sukhavasi performance, with the consideration of limited sensing and performance, with the consideration of limited sensing and communication resources, which motivates research (2013) that such aaClosed processSkew can be propperformance, with the consideration of limitedthe sensing and and and Hassibi Hassibi (2013) showed that properly portrayed via showed the Generalized Normal Hassibi (2013) showed that such such a process process can can be be propcommunication resources, which motivates the2012; research communication resources, which motivates the on the event-based scheduling (Heemels et al., Han and erly portrayed via the Generalized Closed Skew Normal communication resources, which motivates the research research erly portrayed via the Generalized Closed Skew Normal (GCSN) distribution (Genton, 2004).Closed Skew Normal erly portrayed via the Generalized on the event-based scheduling (Heemels et al., 2012; Han on the event-based scheduling (Heemels et al., 2012; Han et al., 2015; Cheng et al., 2016). distribution (Genton, 2004). on the event-based scheduling (Heemels et al., 2012; Han (GCSN) (GCSN) (GCSN) distribution distribution (Genton, (Genton, 2004). 2004). et al., al., 2015; 2015; Cheng Cheng et et al., al., 2016). et et al., 2015; Cheng et al., 2016). 2016).innovation plays the dom- In our view, there exist some inherent connections between In many studies, measurement In our view, there exist some inherent connections between In our view, there exist some inherent connections between theour quantization andinherent event-based scheduling, due view, thereproblem exist some connections between In many many studies, studies, measurement innovation schemes. plays the the domdomIn measurement innovation inating in the design of event-based es- In the quantization problem and event-based scheduling, due In manyrole studies, measurement innovation plays plays the In domthe quantization problem and event-based scheduling, due to the fact that both lead to the hidden-truncation process. the quantization problem and event-based scheduling, due inating role in the design of event-based schemes. In esinating in design In timation theory, innovationschemes. is defined as to the fact that both lead to the hidden-truncation process. inating role role in the themeasurement design of of event-based event-based schemes. In esesto the fact that both lead to the hidden-truncation process. Such a connection motivates us to recheck the Gaussian to the fact that both lead to the hidden-truncation process. timation theory, measurement innovation is defined as timation theory, measurement innovation is defined as the difference between the measurement collected via the motivates us recheck Gaussian timation theory, measurement innovation is defined as Such Such aaa connection connection motivates us to to recheck the the Gaussian approximation assumption in the event-based scheme and connection motivates us to recheck the Gaussian the difference difference between the the measurement collected via via the the Such the between the collected sensor of current time, and one-step measurement assumption in the event-based scheme and the difference between the measurement measurement collected viaprethe approximation approximation assumption in the event-based scheme and develop a novelassumption estimation in algorithm via GCSN distribuapproximation the event-based scheme and sensor of current time, and the one-step measurement presensor current the measurement prediction at theand remote estimator of the last step. develop aa novel estimation algorithm via GCSN distribusensor of ofcomputed current time, time, and the one-step one-step measurement pre- tion. develop novel estimation algorithm via GCSN distribuThe main contribution of this paper can be stated develop a novel estimation algorithm via GCSN distribudiction computed at the remote estimator of the last step. diction the remote of last step. Innovation revealsat the potential contribution of the The main contribution of this paper can be stated diction computed computed at the remote estimator estimator of the the lastnewly step. tion. tion. The contribution of paper can as follows. By introducing GCSN distribution into The main main contributionthe of this this paper can be be stated stated Innovation reveals the the potential contribution of the the newly Innovation reveals contribution of obtained measurement for minimizing the estimation error tion. as follows. By introducing the GCSN distribution into Innovation reveals the potential potential contribution of the newly newly as follows. By introducing the GCSN distribution into the event-based scheduling, an exact MMSE estimation as follows. By introducing the GCSN distribution into obtained measurement measurement for minimizing the estimation estimation error the event-based scheduling, an exact MMSE estimation obtained minimizing the error covariance. In a type offor widely used event-based schedule, obtained measurement for minimizing the estimation error the event-based scheduling, an exact MMSE estimation algorithm which scheduling, does not resort to theMMSE Gaussian approxthe event-based an exact estimation covariance. In a type of widely used event-based schedule, covariance. aa type widely used the packet isIn onlyof its innovation exceedsschedule, the pre- algorithm which doesis not resort to the covariance. Insent type ofwhen widely used event-based event-based schedule, algorithm which resort to Gaussian approxassumption proposed the Gaussian first timeapproxfor the algorithm which does does not not resort for to the the Gaussian approxthe packet packet is sent only only when when its its innovation exceeds the prepre- imation the sent exceeds the designated threshold; not since small innovation assumption is proposed for the first time for the the packet is is sent onlyOtherwise when its innovation innovation exceeds the pre- imation imation assumption is proposed for the first time multidimensional linear Gaussian Markov process. assumption is proposed for the first time for for the the designated threshold; Otherwise prediction not since since small innovation designated threshold; Otherwise innovation implies that the measurement at the remote imation multidimensional linear Gaussian Markov process. designated threshold; Otherwise not not since small small innovation multidimensional linear Gaussian Markov process. multidimensional linear Gaussian Markov process. implies that the measurement prediction at the remote n b implies that the measurement prediction at the remote estimator side has already been closed to the real data implies that the measurement prediction at the remote Notations: Rnn is the n-dimensional Euclidean space. yabb estimator side has You already been closed to theprior real state data Notations: the n-dimensional Euclidean space. yyaab estimator side has the R is the Euclidean space. (Wu et al., 2013; and been Xie, closed 2013). to The col{ya , . . . , R ybn} is the column vector with entries estimator side has already already been closed to the real real data data Notations: Notations: R isrepresents the n-dimensional n-dimensional Euclidean space. ya (Wu et al., 2013; You and Xie, 2013). The prior state ′ col{y , . . . , y } represents the column vector with entries (Wu et al., 2013; You and Xie, 2013). The prior state a b , . . . , y } represents the column vector with entries distribution scheduling is usually , . . . y x denotes the transpose of a vector or a matrix. a b (Wu et al., under 2013; innovation-based You and Xie, 2013). The prior state ycol{y a b col{y , . . . , y } represents the column vector with entries a b ′ distribution under under innovation-based innovation-based scheduling scheduling is is usually usually , . . . , yybb .. x′′ E[(x denotes the transpose or ′ a distribution a the of vector or a a matrix. matrix. a ,, .. .. .. ,, y) ] aisvector the covariance of − E[x])(y − E[y])of distribution innovation-based scheduling is usually yyyCov(x, ⋆ yb . x x denotes denotes the transpose transpose of a The work of under L. He was supported by Natural Science Foundation ′′ a vector or a matrix. ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ ′ ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ random variable (r.v.) x and y; Cov(x, x) is abbreviated as The work of L. He was supported by Natural Science Foundation of China (NSFC) under Grant No. 61503337 and China Postdoctoral The work of L. He was supported by Natural Science Foundation ] is the covariance of Cov(x, y) E[(x − E[x])(y − E[y]) ⋆ The work of L. He was supported by Natural Science Foundation random variable (r.v.) x and y; Cov(x, x) is abbreviated as random variable (r.v.) x and y; Cov(x, x) is abbreviated as of China (NSFC) under Grant No. 61503337 and China Postdoctoral Var(x). Similarly define the conditional variance Var(x|y). Science Foundation under Grant No. 2015M571870; The work of J. of China (NSFC) under Grant No. 61503337 and China Postdoctoral random variable (r.v.) x and y; Cov(x, x) is abbreviated as of China (NSFC) under Grant No. 61503337 and China Postdoctoral Var(x). Similarly define the conditional Science Foundation Grant No. 2015M571870; The work of J. Var(x). Similarly define the variance Var(x|y). p(x) represents probability density variance functionVar(x|y). (pdf) of Chen was supportedunder by NSFC under No. U1401253. Science Foundation under Grant No. 2015M571870; The Var(x). Similarlythe define the conditional conditional variance Var(x|y). Science Foundation under Grant No. Grant 2015M571870; The work work of of J. J. p(x) represents the probability density function (pdf) of Chen was supported by NSFC under Grant No. U1401253. p(x) represents the probability density function Chen was supported by NSFC under Grant No. U1401253. p(x) represents the probability density function (pdf) (pdf) of of Chen was supported by NSFC under Grant No. U1401253. Copyright © 2016, 2016 IFAC 157Hosting by Elsevier Ltd. All rights reserved. 2405-8963 © IFAC (International Federation of Automatic Control) Copyright 2016 IFAC 157 Copyright © 2016 IFAC 157 Peer review© of International Federation of Automatic Copyright ©under 2016 responsibility IFAC 157Control. 10.1016/j.ifacol.2016.10.389
2016 IFAC NECSYS 158 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
Fig. 1. Event-based sensor data scheduling a r.v. P(a ≤ x ≤ b) denotes the probability with the limited range [a, b]. Nn (µ, Σ) is the n-dimensional Gaussian distribution with mean µ and variance Σ, corresponding pdf is denoted as φn (x; µ, Σ). Φn ([a, b]; µ, Σ) represents the integration of φn (x; µ, Σ) within the limited range [a, b], where, a, b ∈ Rn and the integration is component-wise. 2. PROBLEM SETUP We consider the sensor data scheduling for remote state estimation problem as depicted in Fig. 1. The discrete linear time-invariant process is described as follows: xk+1 = Axk + ωk , (1) (2) yk = Cxk + vk , where, xk ∈ Rn is the system state, yk ∈ Rm is the measurement collected via local sensor at time k. Gaussian white noises ωk and vk are assumed to be zero-mean with variances Q ≥ 0 and R > 0, respectively. The initial state x0 is also Gaussian with distribution Nn (µ0 , Σ0 ). Two sequences {ωk } and {vk } are independent, which implies that E[ωi vj′ ] = 0, (∀i, j). The pair (A, C) is detectable and √ (A, Q) is stabilizable. To save scarce communication resource which may be used by other applications, the scheduler equipped at the sensor side decides whether to transmit the current measurement yk over the network to the remote estimator or not at each time step. Denote γk = 1 as the packet is sent and γk = 0 otherwise, which will be shortly specified in Eq.(5). To concentrate on analyzing the senor scheduling scheme, we ignore other imperfect effects such as packet loss, time delay and quantization, etc. Denote the information set available at the estimator side of time k ∈ Z+ as Ik {γ1 y1 , . . . , γk yk }. (3) Assume that there is a feedback channel from the remote estimator to the local scheduler for reliably transmitting the measurement predict yˆk− E[yk |Ik−1 ] at any time instant k ∈ Z+ , where, I0 {x0 }. 1 On that basis, the scheduler can compute the measurement innovation zk specified as Eq. (4) for indicating how accurately the remote estimator can predict the current measurement yk before it sends, and thus how to transmit yk may help to reduce such a prediction error covariance. The measurement innovation is defined as zk yk − yˆk− = yk − C x ˆ− (4) k. We have E[zk ] = 0 due to projection law (Kailath et al., 2000). The scheduling sequence {γk } is specified via the following event-based scheme for a given threshold δ ∈ Rm : 0, if −δ ≤ zk ≤ δ, (5) γk = 1, otherwise. 1 Such a feedback is based on the IEEE 802.15.4/ZigBee protocol. Details can reference Wu et al. (2013).
158
Here, the inequality holds for the component-wise sense. In such a scheme, saving one step of transmission (γk = 0) does not mean the corresponding error covariance increases as that of open-loop prediction. On the contrary, such an event implies that the measurement predict has already been quite reliable, which helps to reduce the communication cost while keeping the estimation error covariance under a pre-designed level. The primary concern of this paper is to propose an exact MMSE estimation algorithm for system (1) ∼ (2) under scheduling scheme (5). We assume that no a prior probability information of state predict is available. 3. PRELIMINARIES This paper chooses the Bayesian method as depicted in Ho and Lee (1964) to obtain the MMSE estimation algorithm. To ease analysis, the equalities frequently used in subsequent sections are provided in this section. 1) Compute the prior pdf via Chapman-Kolmogorov equation p(xk+1 |Ik ) = p(xk+1 |xk )p(xk |Ik )dxk . (6)
2) Obtain the posterior pdf via Bayesian rule on account of the event really happens: • In the case of γk+1 = 1, i.e., Ik+1 = Ik ∪ {yk+1 }. p(xk+1 |Ik+1 ) =
p(yk+1 |Ik , xk+1 )p(xk+1 |Ik ) . p(yk+1 |Ik ) (7)
where, p(yk+1 |Ik ) = p(yk+1 |Ik , xk+1 )p(xk+1 |Ik )dxk+1 .
(8)
• In the case of γk+1 = 0, i.e., Ik+1 = Ik ∪ {−δ ≤ zk+1 ≤ δ}. p(xk+1 |Ik+1 ) P(−δ ≤ zk+1 ≤ δ|Ik , xk+1 ) . (9) =p(xk+1 |Ik ) P(−δ ≤ zk+1 ≤ δ|Ik )
Different from classic Kalman filter, the integrations in (6) and (8) are computationally prohibitive, since γk cannot be expressed as the linear combination of the elements of Ik , which clearly reveals the inherent nonlinear property of the event-based mechanism. For facilitate of analysis, Wu et al. (2013); You and Xie (2013) assumed that the prior state distribution is Gaussian, i.e., − ˆ− p(xk+1 |Ik ) = φ(xk+1 ; x k+1 , Pk+1 ),
(10)
− ′ ˆ− where, xˆ− k+1 E[xk+1 |Ik ], Pk+1 E[(xk+1 − x k+1 )(·) ]. Such an assumption substantially simplifies the computation and approximates the corresponding filter to be linear. Though nice estimation algorithm is obtained, the rationality of such a Gaussianality approximation lacks rigorous verification.
In next two sections, the Gaussian approximation assumption is shown to break down and an exact MMSE estimation algorithm, which does not resort to any approximation, is also provided.
2016 IFAC NECSYS September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
4. GCSN DISTRIBUTION FOR GAUSSIAN HIDDEN-TRUNCATION PROCESS Suppose {xk } and {yk } are jointly Gaussian. For notational convenience, we define three coefficients as follows: Kk Cov(xk , y1k ), Lk Var(y1k ), Mk Var(xk ), which implies that Kk ∈ Rn×(km) , Lk ∈ R(km)×(km) and Mk ∈ Rn×n . We can directly obtain ′ Var(xk |y1k ) = Mk − Kk L−1 k Kk ,
Var(y1k |xk ) = Lk − Kk′ Mk−1 Kk . With a slight confusion of notation, we use Ik as sets when referring to information pattern, and as vectors when referring to the information sets. Then we have the following result. The proof uses some properties of moment generating function of r.v. and is skipped due to space limitation. Lemma 1. The state xk conditioned on the information set Ik can be decomposed into two independent terms: k (11) xk |Ik = ξk + Kk L−1 k y1 |Ik , where, ξk ∼ Nn (0, Var(xk |y1k )). The corresponding estimation error covariance is given by −1 ′ k Var(xk |Ik ) = Var(xk |y1k ) + Kk L−1 k Var(y1 |Ik )Lk Kk . (12)
It is shown in Eq. (12) that the posterior pdf under eventbased scheme is generally non-Gaussian. The Gaussian approximation assumption is valid only when Var(y1k |Ik ) approaches 0, which implies δ → 0 and the measurements are sent at any time instant. Under this extreme situation, k k xk |Ik → xk |y1k ∼ Nn (Kk L−1 (13) k y1 , Var(xk |y1 )), which reduces to the classic MMSE estimator (Kailath et al. (2000), pp.81). The last term of Eq. (12) explicitly quantifies the performance deterioration incurred from reducing the communication cost. As expected, the problem becomes substantially challenging when dropping the Gaussian approximation assumption and the techniques employed in the aforementioned works (Wu et al. (2013); You and Xie (2013)) no longer apply. To accurately characterize the Gaussian hidden-truncation process from the event-trigger scheme (5), we develop a new framework by introducing the closed skew normal (CSN) distribution. Definition 2. (Genton (2004), pp.26). Consider m ≥ 1, n ≥ 1, µ ∈ Rn , ν ∈ Rm , D is an arbitrary m × n matrix, Σ and ∆ are positive definite matrices of dimension n×n and m × m, respectively. Then the pdf of a CSN distribution is given by CSNn,m (x; µ, Σ, D, ν, ∆) Φm (−∞, D(x − µ); ν, ∆) . (14) =φn (x; µ, Σ) Φm (−∞, 0; ν, ∆ + DΣD′ ) Here, µ, Σ, D can be deemed as location, scale and skewness parameters, respectively. When taking D = 0, Eq.(14) reduces to the pdf of the Gaussian distribution with mean µ and variance Σ. It should be noted that for the general CSN distribution, µ and Σ cannot be deemed as mean and variance, which are in fact involved to compute and can be found in Genton (2004). Intuitively, for independent r.v. X ∼ Nn (µ, Σ) and Z ∼ Nm (ν, ∆), Eq. (14) explicitly characterizes the pdf of the conditional event (Genton (2004), pp.28) 159
159
X|{Z − D(X − µ) < 0}, where, D is the weighted matrix. For the bilateral truncation scheme (5), we need the modified version of Eq. (14), which is called the generalized closed skew normal (GCSN) distribution. Definition 3. (Sukhavasi and Hassibi (2013)). The pdf of the GCSN distribution is given by: GCSNn,m (x; µ, Σ, D, Sm , ∆) Φm (Sm ; D(x − µ), ∆) . (15) =φn (x; µ, Σ) Φm (Sm ; 0, ∆ + DΣD′ ) Here, Sm is the region determined via the specified realization. With such a tool for characterizing the bilateral truncationprocess, an exact MMSE estimation algorithm for system under event-based scheme is proposed in next section. 5. MMSE ESTIMATION ALGORITHM UNDER EVENT-BASED SCHEME WITH GCSN DISTRIBUTION Before presenting the main result of the paper, this section starts with some lemmas. Specifically, Lemma 4 and Lemma 5 play the dominating role in computing the multiplication of two Gaussian pdfs. Lemma 4 was put forward in Ho and Lee (1964), corresponding proof uses the mathematical trick known as “completion of squares” for multivariate argument. Lemma 5 uses Schur complement lemma (Kailath et al. (2000), pp.725). Most of the proofs are omitted in consideration of limited space. Lemma 4. Given x, a ∈ Rn , y ∈ Rm and symmetric matrices A ∈ Rn×n , B ∈ Rm×m , which are positively defined such that CAC ′ + B is non-singular (C ∈ Rm×n ). Then (x − a)′ A−1 (x − a) + (y − Cx)′ B −1 (y − Cx) =(x − µ)′ Σ−1 (x − µ) + (y − Ca)′ (CAC ′ + B)−1 (y − Ca), where, µ a + AC ′ (CAC ′ + B)−1 (y − Ca), Σ (A−1 + C ′ B −1 C)−1 .
Lemma 5. For X ∈ Rn×n , C ∈ Rm×n and R ∈ Rm×m , the following equality holds true: |Y ||CXC ′ + R| = |X||R|, (16)
where, Y (X −1 + C ′ R−1 C)−1 .
As will be shown in the proof of Theorem 8, Lemma 6 plays the dominating role in computing P(−δ ≤ zk+1 ≤ δ|Ik ) of Eq. (9). Proposition 7 is needed in calculating the integration of Eq. (6) and Eq. (8). Lemma 6. Let a, b ∈ Rm , p, q ∈ Rn and B ∈ Rm×n . If V ∼ Nn (µ1 , Σ), we have q Φm ([a + BV, b + BV ]; µ2 , Ω)φn (V ; µ1 , Σ)dV p a + Bµ1 − µ2 b + Bµ1 − µ2 =Φm+n , ; 0, ΣW , p − µ1 q − µ1 where,
ΣW
BΣB ′ + Ω −ΣB ′
− BΣ . Σ
(17)
2016 IFAC NECSYS 160 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
On the basis of Lemma 6, the following proposition can be easily obtained. Proposition 7. Let a, b ∈ Rm and B ∈ Rm×n . If V ∼ Nn (µ1 , Σ), we have EV [Φm ([a + BV, b + BV ]; µ2 , Ω)] =Φm ([a + Bµ1 − µ2 , b + Bµ1 − µ2 ]; 0, BΣB ′ + Ω). Now we shall tie together the ideas of the previous lemmas to state the main result of this paper. Theorem 8. The pdfs of the MMSE estimation algorithm for system (1) ∼ (2) under scheme (5) is given as follows: p(xk+1 |Ik ) − − − =GCSN (xk+1 ; µ− k+1 , Σk+1 , Dk+1 , [Ak , Bk ], ∆k+1 ), p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , [Ak+1 , Bk+1 ], ∆k+1 ). Corresponding coefficients are updated as: 1) Predict µ− k+1 = Aµk , ′ Σ− k+1 = AΣk A + Q, − −1 = Dk Σk A′ (Σ− , Dk+1 k+1 )
= ∆k +
Dk Σk Dk′
2) Correction • If γk+1 = 1 ′ Kk+1 = Σ− k+1 C , Lk+1 =
−
′ CΣ− k+1 C
We begin by verifying the recursion of the GCSN distributions, which does not resort to the event-based scheme. Consider a stable system with A = 0.9, C = 1, Q = 0.5, R = 0.5 and the initial state x0 ∼ N (5, 1). Let the packet of the first step fail to send and the remaining transmissions all be successful. Fig. 3 shows the correctness of the algorithm (18) ∼ (20). µ=0,Σ=5 µ=0,Σ=5,a=−0.5,b=10,D=5,∆=5 µ=0,Σ=5,a=−50,b=500,D=50,∆=0.1
0.45
+ R,
−1 ˜k+1 , µk+1 = µ− k+1 + Kk+1 Lk+1 y
Σk+1 =
This section leverages some numerical examples to show the effectiveness of the proposed estimation algorithm. Fig. 2 gives the pdf of one Gaussian distribution and two GCSN distributions. The blue curve shows that µ is no longer the expectation of a GCSN distribution. The red curve shows the truncation effect.
0.5
y˜k+1 = yk+1 − Cµ− k+1 , Σ− k+1 − − Dk+1 ,
(18)
− − ′ Dk+1 Σ− k+1 (Dk+1 ) .
6. SIMULATION
0.4 0.35
(19)
′ Kk+1 L−1 k+1 (Kk+1 ) , − ∆k+1 = ∆k+1 .
0.3 pdf
∆− k+1
Remark 10. We notice from Eq. (21) ∼ (22) that the dimensions of some coefficients linearly increase when the packets are not sent due to the event-based scheme. Such a problem arises from the bilateral hidden-truncation process (5) and can be explained via Lemma 6, where state-augmentation approach is adopted for dealing with the integration of limited range. One possible way to deal with this problem is to transform mechanism (5) in a way that turns ΣW of Eq. (17) into a diagonal matrix, which will be left as our future work.
Dk+1 = The region is updated as Ak+1 = Ak + Dk+1 (µ− k+1 − µk+1 ),
0.25 0.2 0.15 0.1
(20) Bk+1 = Bk + Dk+1 (µ− k+1 − µk+1 ). • If γk+1 = 0 µk+1 = µ− , , Σk+1 = Σ− k+1− k+1 0 Dk+1 ∆− k+1 , (21) Dk+1 = , ∆k+1 = 0 R C
0.05 0 −6
160
−2
0
2
4
6
8
x
Fig. 2. PDF of a GCSN distribution
x− − µ− ), Θk = C(ˆ k+1 k+1 Ak Bk , Bk+1 = . (22) Ak+1 = Θk − δ Θk + δ
8 Real State x
k
Estimated State
6
4
2 State
Both (20) and (22) are updated with the initial conditions Ak0 −δ, Bk0 δ, (23) where, k0 denotes the first time instant when the packet is not sent. Remark 9. To the best of our knowledge, Theorem 8 first proposes an exact MMSE estimation algorithm for the multidimensional linear Gaussian Markov process with the threshold-type event-based scheme. The introduced GCSN distribution avoids the Gaussian approximation assumption widely used in the event-based schedules, i.e., Wu et al. (2013); You and Xie (2013), which is proved to be inaccurate in Lemma 1.
−4
0
−2
−4
−6
0
20
40
60
80
100
Time
Fig. 3. State versus estimate under GCSN distribution
2016 IFAC NECSYS September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
As for the event-based schedule, consider an unstable system with A = 1.2, C = 1, Q = 0.25, R = 0.25 and x0 ∼ N (1, 1). The threshold for the event-based scheme (5) is δ = 0.2. The comparison between the real state and the corresponding estimate is given in Fig. 4. Fig. 5 shows the packet transmission sequence under event-based scheme. Here, 1 denotes the packet is not sent, which is different from Eq. (5). The estimation error and the trace of error covariance are provided in Fig. 6 and Fig. 7, respectively. The effectiveness of the proposed MMSE algorithm under event-based scheme via GCSN distribution is thus verified. 7. CONCLUSION This paper comes up with an event-based scheme for remote state estimation over a bandwidth-limited channel. Our main contribution lies in the design of an exact MMSE estimation algorithm that accurately exploits the additional information of the event and does not resort to the Gaussian approximation assumption. There are many interesting extension directions to explore for this work: finding more appropriate event-based scheme to reduce computational complexity; establishing the quantitative relationship among estimation performance, communication rate and event-triggered threshold; comparing system performance of the present algorithm with those in Wu et al. (2013) and Han et al. (2015), etc.
161
REFERENCES ˚ Astr¨om, K. and Bernhardsson, B. (2002). Comparison of Riemann and Lebesque sampling for first order stochastic systems. In Proceedings of the 41st IEEE Conference on Decision and Control, volume 2, 2011–2016. IEEE. Chen, J., Yu, Q., Cheng, P., Sun, Y., Fan, Y., and Shen, X. (2011). Game theoretical approach for channel allocation in wireless sensor and actuator networks. IEEE Transactions on Automatic Control, 56(10), 2332– 2344. Cheng, P., Qi, Y., Xin, K., Chen, J., and Xie, L. (2016). Energy-efficient data forwarding for state estimation in multi-hop wireless sensor networks. IEEE Transactions on Automatic Control, 61(5), 1322–1327. Genton, M. (2004). Skew-elliptical Distributions and Their Applications: A Journey Beyond Normality. CRC Press. Han, D., Mo, Y., Wu, J., Sinopoli, B., and Shi, L. (2015). Stochastic event-triggered sensor schedule for remote state estimation. IEEE Transactions on Automatic Control, 60(10), 2661–2675. Heemels, W., Johansson, K., and Tabuada, P. (2012). An introduction to event-triggered and self-triggered control. In Proceddings of the IEEE 51st Annual Conference on Decision and Control, 3270–3285. IEEE. Ho, Y. and Lee, R. (1964). A Bayesian approach to problems in stochastic estimation and control. IEEE Transactions on Automatic Control, 9(4), 333–339.
20000 Real State xk
1.5
Estimated State 15000
1
0.5
State
Estimation error
10000
5000
0
−0.5 0 −1 −5000
0
10
20
30 Time
40
50
60
−1.5
Fig. 4. State versus estimate under event-based schedule
0
10
20
30 Time
40
50
60
30 Time
40
50
60
Fig. 6. Estimation error
1 0.22
0.9
Trace of the estimation error covariance
0.8
Packet dropout
0.7 0.6 0.5 0.4 0.3 0.2 0.1 0
0
10
20
30 Time
40
50
0.21
0.2
0.19
0.18
0.17
60 0.16
Fig. 5. Scheduling sequence: 1 represents the packet is not sent and 0 represents successful transmission 161
0
10
20
Fig. 7. Trace of the estimation error covariance
2016 IFAC NECSYS 162 September 8-9, 2016. Tokyo, Japan
Lidong He et al. / IFAC-PapersOnLine 49-22 (2016) 157–162
Johansson, K., Pappas, G., Tabuada, P., and Tomlin, C. (2014). Guest editorial special issue on control of cyberphysical systems. IEEE Transaction on Automatic Control, 59(12), 3120–3121. Kailath, T., Sayed, A., and Hassibi, B. (2000). Linear Estimation. Prentice Hall. Sukhavasi, R. and Hassibi, B. (2013). The Kalman-like particle filter: Optimal estimation with quantized innovations/measurements. IEEE Transactions on Signal Processing, 61(1), 131–136. Wu, J., Jia, Q., Johansson, K., and Shi, L. (2013). Eventbased sensor data scheduling: Trade-off between communication rate and estimation quality. IEEE Transactions on Automatic Control, 58(4), 1041–1046. You, K. and Xie, L. (2013). Kalman filtering with scheduled measurements. IEEE Transactions on Signal Processing, 61(6), 1520–1530. Appendix A. PROOF OF THEOREM 8 Proof. The proof is conducted via mathematical induction. Here, dimensions are omitted for sake of brevity. 1) Prediction of step 1 Using C-K equation (6) for system (1) with initial − distribution p(x0 |I0 ) = φ(x0 ; µ0 , Σ0 ), we have µ− 1 , Σ1 as in Eq. (18) when taking k = 0. 2) Correction of step 1 • If the packet is sent, i.e., I1 = I0 ∪ {γ1 = 1} From Bayesian rule, we have p(x1 |I1 ) = φ(x1 ; µ1 , Σ1 ), where, µ1 , Σ1 are given as Eq. (19). • If the packet is not sent, i.e., I1 = I0 ∪ {γ1 = 0} Eq. (9) can be computed as p(x1 |I1 ) =
− − Φ([A1 , B1 ]; C(x1 − µ1 ), R) . φ(x1 ; µ− 1 , Σ1 ) − ′ Φ([A1 , B1 ]; 0, CΣ1 C + R)
where, A1 , B1 are as in Eq. (23). When taking − µ1 = µ − 1 , Σ1 = Σ1 , D1 = C, ∆1 = R. we have p(x1 |I1 ) = GCSN (x1 ; µ1 , Σ1 , D1 , [A1 , B1 ], ∆1 ). 3) Prediction of step k+1 Take Sk [Ak , Bk ]. Without loss of generality, assume p(xk |Ik ) = GCSN (xk ; µk , Σk , Dk , Sk , ∆k ). From Eq. (6), we have ∞ p(xk+1 |Ik ) = φ(xk+1 ; Axk , Q)φ(xk ; µk , Σk ) −∞
Φ(Sk ; Dk (xk − µk ), ∆k ) dxk . Φ(Sk ; 0, ∆k + Dk Σk Dk′ ) Using both Lemma 4 and Lemma 5, as well as pdf of the Gaussian distribution, φ(xk+1 ; Axk , Q)φ(xk ; µk , Σk ) 12 ♯ |Σk ||Σ− k+1 | − = φ(xk ; µ♯k , Σ♯k )φ(xk+1 ; µ− k+1 , Σk+1 ), |Σk ||Q| ×
− where, µ− k+1 , Σk+1 are given as Eq. (18), and −1 (xk+1 − µ− µ♯k µk + Σk A′ (Σ− k+1 ) k+1 ),
Σ♯k [(Σk )−1 + A′ Q−1 A]−1 . In view of both Lemma 5 and Proposition 7, as well − as taking Dk+1 and ∆− k+1 as Eq. (18), we have 162
p(xk+1 |Ik ) − − − =GCSN (xk+1 ; µ− k+1 , Σk+1 , Dk+1 , Sk , ∆k+1 ). 4) Correction of step k+1 • If Ik+1 = Ik ∪ {γk+1 = 1} Computing p(yk+1 |Ik ) can conduct the similar way as computing p(xk+1 |Ik ). Further using the Gaussian property of p(yk+1 |xk+1 ), yielding p(xk+1 |Ik+1 ) =φ(xk+1 ; µk+1 , Σk+1 ) − − Φ(Sk ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) × . − − ′ Φ(Sk+1 ; 0, ∆− k+1 + Dk+1 Σk+1 (Dk+1 ) ) where, µk+1 and Pk+1 are given in Eq. (19), Ak+1 and Bk+1 are as in Eq. (20). Further take Dk+1 , ∆k+1 as in Eq. (19), we have p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , Sk+1 , ∆k+1 ). • If Ik+1 = Ik ∪ {γk+1 = 0} The crucial step in Eq. (9) is to compute P(−δ ≤ zk+1 ≤ δ|Ik ) Ωk+1 = . (A.1) − − − ′ Φ(Sk ; 0, ∆k+1 + Dk+1 Σ− k+1 (Dk+1 ) ) where, Ωk+1 yˆ− +δ k+1 − ′ φ(yk+1 ; Cµ− k+1 , CΣk+1 C + R) − −δ yˆk+1
× Φ([Mk+1 , Nk+1 ]; 0, Λk+1 )dyk+1 ,
− − Mk+1 Ak − Dk+1 Kk+1 L−1 k+1 (yk+1 − Cµk+1 ), − − Nk+1 Bk − Dk+1 Kk+1 L−1 k+1 (yk+1 − Cµk+1 ),
− − ′ Λk+1 ∆− k+1 + Dk+1 Σk+1 (Dk+1 ) . By virtue of Lemma 6, we have Ωk+1 = Φ([Ak+1 , Bk+1 ]; 0, ΣW ), where, Ak+1 , Bk+1 are given as Eq. (22), and − − − − ′ Dk+1 Σ− k+1 (Dk+1 ) + ∆k+1 Dk+1 Kk+1 . ΣW = − ′ Kk+1 (Dk+1 )′ Lk+1 Augment the product of the following two cdf− − s Φ(Sk ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) as well as − − Φ([ˆ yk+1 − δ, yˆk+1 + δ]; Cxk+1 , R) into another Gaussian cdf with dimension m+n, then combine it with Eq. (A.1), yielding p(xk+1 |Ik+1 ) − =φ(xk+1 ; µ− k+1 , Σk+1 )
Φ(Sk+1 ; Dk+1 (xk+1 − µ− k+1 ), ∆k+1 ) , Φ(Sk+1 ; 0, ΣW ) where, Dk+1 , ∆k+1 are given as Eq. (21). When − taking µk+1 = µ− k+1 , Σk+1 = Σk+1 , we have ×
ΣW = ∆k+1 + Dk+1 Σk+1 (Dk+1 )′ . In a word, p(xk+1 |Ik+1 ) =GCSN (xk+1 ; µk+1 , Σk+1 , Dk+1 , Sk+1 , ∆k+1 ). The induction step is thus complete.