Resilient State Estimation for Discrete-Time Linear Systems

Resilient State Estimation for Discrete-T ime Linear Systems Alexandre Kircher 1 , La urent Bako 1 , Eric Blanco 1 , Moha med Benallouch 2 Abstract — This paper proposes a resilient state estimator for L TI discre te-time systems. The dyn amic equ ation of the system is assumed to be affected by a bounded process noise. As to the a vailable measurements, they ar e potentially corrupted by a noise of both dense and impul siv e natu res. In th is settin g, we construct the estimator as the m ap which associates to the measureme nts, the minimizin g set of an appropriate (conv ex) perfo rmance function . It is then shown that the proposed estimator enjoys the property of resilience, that is, it ind uces an estimation error which, un der certain condi tions, is inde- pendent of the extreme values of the (impulsive) measurement noise. Theref ore , the estimation error may be bounded while the measurement noise is virtually unbound ed. Moreov er , th e expression of the b ou n d depend s explicit ly on the degree of observability of the system being observed and on the considered perf ormance function. Finally , a few simulation results are p ro vided to illustrate the resilience property . Index terms —Secure state estimation, sensor attacks, outli ers, resilient estimators, Cyber -physical systems. I . I N T RO D U C T I O N W e consider in this work the pr oblem of designing state estimators which would be resilient against an (unknown) sparse noise sequence affecting th e measurements. By sparse noise we re f er here to a sign al sequen c e which is of i mpulsive nature, that is, a seq uence which is most of the time eq ual to zero, except at a few instants wh ere it can take on ar bitrarily large values. Th e problem is relevant for example, in the su- pervision of Cyber-Physical Systems [5]. In th is application, the supervisory d ata may b e collected by spatially distributed sensors and then sent to a distant pro cessing unit through some com municatio n network. Durin g the transmission, the data m ay in cur inter mittent packet losses or adversarial attacks con sisting in e. g., the injection o f arbitrary signals. This estimation problem was investigated through many different approach es. Since the measuremen ts are assumed to b e affected b y a sequence o f o utliers which is sparse in time, a natu ral scheme of solutio n to the state estimatio n problem may be to first de tect the occurr ences of the nonzero instances of that sparse noise, remove the co r rupted data and then pr oceed with classical estimation method s such as the Kalman filter or Luenb e rger ty pe of ob server [13], [15]. Another category of app roaches, w h ich are inspired by some recen t results in compre ssive s ampling [4] , [8], rely on sparsity-indu cing optimization tec hniques. A striking feature of these metho ds is that they do no t treat separately the tasks of d etection, data clean ing and estimation . Instead, an 1 A. Kircher , L. Bako and E. Blanco are with Uni ver sité de L yon, Laboratoi re Ampère (Ecole Centrale L yon, CNRS UMR 5005), Ecully F- 69134. E-mail: alexand re.kircher@ec -lyon.fr 2 M. Benallouch is with E CAM L yon, 40 Montée Saint-Barthél émy , 69321 L yon, France . implicit discrimination of the wrong data is induced by some specific prop erties o f the to-be-minim ized cost fu nction. One of the first work s that puts forward this ap proach for the resilient state estimation prob lem is the one reporte d in [7]. There, it is assumed that only a fixed number of sensors are subject to attacks (sparse but otherwise arb itrary disturbanc e s). Th e challeng e then resides in the fact that at each time instant, one doe s not know which senso r is compro mised. Note howev er that the assumptions in [7] were qu ite restricti ve as no process no ise or measurem ent noise (o ther than the sparse attack signal) was con sidered. These limitations ope n ways for later extensions in many directions. For example, [18] s uggests a ref ormulatio n which reduces c omputatio nal c o st by using the concep t of ev ent- triggered up date; [14] c onsiders an observation model which includes d ense no ise alon g with the sparse attack signal. In [6], the assum ption of a fixed nu mber of attacked sensors is relaxed. Finally , th e recent pap er [11] propo ses a unified framework fo r analyz ing resilienc e capabilities o f m ost o f these optimization- based estimators. Although a b ound on the estimation er ror was derived in this p aper, it is not quantitatively related to the p roperties (e.g., obser vability) of the dy namic system being ob ser ved. The contribution of the curr ent paper is the design of a (conve x) o ptimization - based resilient estima to r for L TI discrete-time systems. The available model of the system assumes b ounde d n oise in both the dy namics an d the ob- servation equation with th e latter being possibly affected by an un known but sparse attack sign al. Contrary to th e settings in so m e existing works, we did not imp ose h ere a ny restriction on the number of sensors which are subject to attacks, that is, any sensor can be compr omised at any time. Our main theor etical result concerns th e r e silien c e an a lysis of the pro posed estimator . W e show that th e estimatio n error associated with the n ew estimator can be ma d e, under certain conditions, insensitive to th e amplitude of the attack signal. Our bound , although necessarily conservati ve, has the impo r tant advantage o f being explicitly expressible in function of the pro perties of the considered dy namic system. This m akes it a valuable qualitative tool for assessing the impact of the estimator’ s d esign par ameters and th a t of the system matrices on th e qu ality of the estimation. For example, it reflects the in tuition that the more observable the system is, the larger th e num ber of instances o f gro ss values (of the output noise) it can han dle and the smaller the error boun d. Outline. The rest of the paper is organ ized as follows. The estimation setting is defined in Section I I. In Section III we elab orate o n the pro posed optimization -based estimator: Necessary te c hnical tools are introd uced in Section III-A for the statement and the proof o f the main result in Section III-B. Section I V illustrates the p erform ance of the estimation method in simulation; Sec tio n V provides some con cluding remarks. Notatio ns. Throu g hout th is p aper, R ≥ 0 (respectively R > 0 ) designates the set of no nnegative (r espectively positiv e) reals. W e note R a the set of (column) vectors with a real elemen ts and for any vector z in R a , z i with i in { 1 , ..., a } is the i -th com p onent of z . Mo reover , R a × b is the set of r e al matrices with a rows an d b c o lumns. If M ∈ R a × b , then M ⊤ will designate the transposed matrix of M . No tation k·k will repr esent a given norm over a giv en set (which will be specified when necessary ). k·k 2 is the Euclid ean no rm, defined by k z k 2 = √ z ⊤ z fo r all z in R a . k·k 1 will designate the ℓ 1 -norm , defined by k z k 1 = P a i =1 | z i | for z ∈ R a . For a finite set S , th e notation |S | will r e f er to the cardinality of S . I I . T H E E S T I M AT I O N P R O B L E M Consider the following discrete-time L inear T ime- In variant (L TI) system Σ :  x t +1 = Ax t + w t y t = C x t + f t (1) where x t ∈ R n is the state vector at tim e t , y t ∈ R n y is the outp ut vector a t time t ; A ∈ R n × n the dyn amic matr ix of the system and C ∈ R n y × n y is the ob servation matrix . w t ∈ R n and f t ∈ R n y model respectively the process noise and the outp ut noise both of wh ich are unkn own. W e shall howe ver make the inform al assumption s th at { w t } is bo unded with a r elativ ely small amplitude. As to the sequence { f t } , it can take on p otentially ar b itrarily large values, that is, no e xplicit bound is imp osed on its amplitude . This type of n oise can model for example, ordinary measure- ment noise (of ‘m oderate amount’) together with interm itten t faulty measurem e nts, attack signals o r packet losses o n d ata transmitted ov er a com m unication ne twork. For conv enience, one can also view f t as the sum of two n o ise co mpone n ts, a dense no ise, representing a b o unded noise induced by the sensors, and a sparse noise sequenc e, i.e. , a noise who se instances are equal to zero mo st of the time but who se non- zero elements can take o n arbitrary values. Problem. The problem con sidered in this paper is one o f estimating the states x 0 , . . . , x T − 1 of the system (1) on a time period T giv en T m e asurements y 0 , ..., y T − 1 of the system output. W e shall seek an acc u rate est imate o f th e state despite the uncertainties in the system equations (1) modeled by w t and f t the chara cteristics o f which are described above. In par ticular, we would like the to-be-d esigned estimato r to produ ce an estimate such th at the estimation error is, when possible, ind epende n t of the maximum am plitude of { f t } . Such an estimato r will then be called resilient. I I I . R E S I L I E N T O P T I M I Z A T I O N - B A S E D E S T I M ATO R W e pro p ose a co n vex optimization -based so lution to th e state estimation prob lem d efined above. Giv en the system matrices A and C and T outp u t measurem e nts y 0 , ..., y T − 1 , consider a p erform ance f unction F : R n × T → R ≥ 0 defined by F ( Z ) = λ X t ∈T ′ k z t +1 − Az t k 2 2 + X t ∈T k y t − C z t k 1 , (2) where T = { 0 , . . . , T − 1 } , T ′ = { 0 , . . . , T − 2 } and Z =  z 0 · · · z T − 1  , i.e., the vectors z t ∈ R n are the columns o f the matrix Z . Here, λ > 0 is a user-defined parameter which aims at balancing th e co ntributions of the two terms inv o lved in the expression of the pe r forman ce index F . This idea o f weigh ting the ter ms contain ed in F could also be don e differently depen d ing on the time ind ex, for example by takin g terms of the form k W t ( z t +1 − Az t ) k 2 2 and k V t ( y t − C z t ) k 1 , where W t and V t would b e positive- definite weighting matrices. Let P ( R n × T ) den ote the collection of subsets o f R n × T . Then th e propo sed estimator is defined as the set-valued map Ψ : R n y × T → P ( R n × T ) which maps the av ailable measuremen ts Y ,  y 0 · · · y T − 1  to the subset Ψ ( Y ) of R n × T defined by Ψ( Y ) = arg min Z ∈ R n × T F ( Z ) . (3) By assuming that the pair ( A, C ) is obser vable, it can be checked that F is coer civ e, i.e. , it satisfies lim k Z k → + ∞ F ( Z ) = + ∞ for any norm k·k on R n × T . It follows that the estimator Ψ expressed in (3 ) is well-defined in the sense that the u nderly ing optimization pr o blem in (3) admits a solution [16]. Note howe ver th at the minimizer need not be un ique. M o reover , since the objective fun c tion F is conv ex, the elem e nts of the so- defined state estimato r Ψ( Y ) can be d etermined efficiently for a giv en Y . Many numerical solvers ca n be used fo r this pur p ose, see e.g . [10], [1], [1 9] for the comp utational asp e c ts. The rest of the paper will focus on assessing the resilien c e proper ties of the estimator (3). For this purp ose we need some prelimina r y tech nical results. A. Pr e liminaries T o begin with the an alysis, we introduce some u sef ul technical too ls, the first of which is th e class of K ∞ function s (see, e.g., [12]). T h is class of fun ctions will be used to measure the incre a sing ra te of the estimation error . Definition 1 (class- K ∞ function s) . A function α : R ≥ 0 → R ≥ 0 is said to be of class- K ∞ if it is co ntinuo u s, zer o at zer o, strictly in cr easing a nd satisfies lim s → + ∞ α ( s ) = + ∞ . Using this defin ition we can state a technical lem ma wh ich will play an importan t role in the analysis. Lemma 1. Let G : R n × m → R ≥ 0 be a nonnegative continuo us function sa tisfying the following pr operties: • P ositive defi n iteness: G ( S ) = 0 if an d only if S = 0 • Relaxed ho mogeneity: Ther e exists a K ∞ function σ such th at G ( S ) ≥ σ ( 1 λ ) G ( λS ) for a ll λ ∈ R > 0 . Then for an y norm k·k o n R n × m , ther e exists d > 0 such that for all S ∈ R n × m , G ( S ) ≥ dσ ( k S k ) . Pr oo f. W e start by observin g that the u nit hyper sphere D = { S ∈ R n × m : k S k = 1 } is a com pact set in the topolog y induced by th e norm k·k . By the extreme v alue theorem, G b eing continu ous, admits nece ssarily a min imum value on D , i.e., there is S ⋆ ∈ D such that G ( S ) ≥ d , G ( S ⋆ ) > 0 for all S ∈ D . For a ny nonzero S ∈ R n × m , S k S k ∈ D so that G ( S k S k ) ≥ d . On the other hand, by the r elaxed homog e n eity of G , G ( S ) ≥ σ ( k S k ) G ( S k S k ) ≥ dσ ( k S k ) . Moreover , this inequa lity holds for S = 0 . It ther e fore holds true for any S ∈ R n × m . For future uses in the paper, consid e r n ow the function H : R n × T → R ≥ 0 defined by H ( Z ) = λ 2 X t ∈T ′ k z t +1 − Az t k 2 2 + X t ∈T k C z t k 1 (4) Note the r esemblance between F ( Z ) and H ( Z ) . Th ey only d iffer by the ab sen ce of y t in the second term of H and the factor of the first term wh ich is λ in the first case and λ/ 2 in th e second . Lemma 2 (Lower Boun d on H ) . Let k·k b e a norm on R n × T . Con sid er the function H defined in (4) un der th e assumption that ( A, C ) is observable. Th en H ( Z ) ≥ Dq ( k Z k ) ∀ Z ∈ R n × T (5) wher e q : R ≥ 0 → R ≥ 0 is the function d efined b y ∀ α ∈ R ≥ 0 , q ( α ) = min( α, α 2 ) (6) and D = min k Z k =1 H ( Z ) > 0 . (7) Pr oo f. Th e idea of th e proo f is to ch eck that H satisfies the condition s of Le mma 1 and then ap ply it to conclu de. First, note that continu ity and nonn egati vity of H are obvious. As to the relaxed homo geneity prope rty , it can be checked straightfor ward ly that it ho lds with σ = q . Finally , setting H ( Z ) = 0 im plies that z t +1 = Az t and C z t = 0 for all t = 0 , . . . , T − 1 . It immediately f o llows th at C A t z 0 = 0 and so, O z 0 = 0 where O =  C ⊤ · · · ( C A n − 1 ) ⊤  ⊤ is the ob ser vability matr ix of the system. By the obser vability assumption, we g et that z 0 = 0 an d consequen tly , that Z = 0 . Therefo re H is positive-definite. The statement of the lemma now follows b y applyin g Le m ma 1 . T o proceed fur ther, let u s intro duce a f ew no tations. W e use the no ta tio n I = { 1 , . . . , n y } to denote a label set for the sensors describ ed by the ob servation equa tion in (1) an d T = { 0 , . . . , T − 1 } to the set of time ind exes. For i ∈ I , c ⊤ i denotes the i -th r ow of the observation matrix C . The next defin itio n intr oduces a p arameter to gau ge th e resilience properties of an estimator of the form de fin ed in (3). Definition 2 ( r - Resilience ind ex p r ) . Let r b e a nonnegative inte ger . Assume tha t the system Σ in (1) is o bservable. W e define th e r -Resilience index of the estimator Ψ in (3) (when applied to Σ ) as the real nu mber p r given by p r = sup Z 6 =0 Z ∈ R n × T sup Λ r ⊂I ×T | Λ r | = r P ( i,t ) ∈ Λ r   c ⊤ i z t   H ( Z ) (8) wher e H is as defined in (4) . The supr emum is taken her e over a ll nonzer o Z in R n × T and over all subsets Λ r of I × T with car dinality r . The index p r can be interpreted a s a quantitativ e measu re o f the obser vability of the system Σ . The observability is needed here to ensure that the den ominator H ( Z ) o f (8) is different from zero when ever Z 6 = 0 (see the positi ve definiteness proof o f H in the pr oof Lemma 2 above). Furthermo re, it should b e r e m arked that P ( i,t ) ∈ Λ r   c ⊤ i z t   ≤ H ( Z ) for any Λ r ⊂ I × T , whic h implies that the defining suprema o f p r are well-defined . What the r -Resilience parame ter assesses is h ow much the estimator can ha n dle d ata corru ption as it represents the worst ratio between the weight of r c orrup ted estimates (which take any value and be p otentially placed anywhere in time) an d the weight o f th e who le estimated trajectory . As a result, the lower p r is, the m ore resilient the estimator is expected to be. The next section gives more backg roun d to the introd uction of p r and which role it exactly p lays in the resilience analysis of the estimato r . From a comp u tational viewpoint we observe that the parameter p r is hard to comp ute in general. In effect, obtaining p r numerically would req uire solvin g a nonconve x and com binatoria l optimization p roblem. Th is is in deed a common ch aracteristic of the con c epts which a re usually used to assess resilienc e; for example the p opular Restricted Isometry Property (RIP) co n stant [3] is comparatively as hard to ev aluate. Nevertheless, if we restrict attentio n to estimation problems wher e the process noise { w t } would be identically equal to zer o, then by ad ding in (8) the additional constraint that z t +1 = Az t , p r can be e xactly computed using the me th od in [17] or mor e cheaply overestimated using the one in [2]. B. Characterization of the r esilience pr o perty The main result o f this paper con sists in the character iza- tion of the resilience property of the state estimator (3 ). More specifically , our result states that the estimation error, i.e. , the difference between the real state tra je c tory and th e estimated one, is up per bounde d by a bound which d oes not depen d on the amplitud e of the ou tlier s containe d in { f t } provid ed that the nu mber of such outlier s is below some th reshold. Before stating the main th eorem, let u s introd u ce a last notation to be used in the analysis. Le t ε ≥ 0 be a g iv en number . For any ad missible sequen c e { f t } t ∈T in (1), we can split the in dex set I × T into two disjoint label sets, J ε = { ( i, t ) ∈ I × T : | f it | ≤ ε } , (9) indexing those 1 f it which are boun ded by ε and J c ε = { ( i, t ) ∈ I × T : | f it | > ε } indexing those f it which are possibly un boun ded. I t is importan t to keep in min d that ε is ju st a param eter for deco m posing the noise sequence in two p arts in v iew of the analy sis (and not a bo und on f it ). The pa r ticular situation where ε = 0 reflects the appro ach where one would view any n onzero f it as an outlier . Theorem 1 (Upper bo und o n the estima tio n error) . Con sider the system Σ defined b y (1) with o utput measur ement Y and consider the estimator (3) . Let ε ∈ R ≥ 0 and r = |J c ε | . If Σ is observable and p r < 1 / 2 , then fo r all ˆ X =  ˆ x 0 · · · ˆ x T − 1  ∈ Ψ( Y ) , k E k ≤ h  2 β Σ ( ε ) D (1 − 2 p r )  (10) wher e E =  ˆ x 0 − x 0 · · · ˆ x T − 1 − x T − 1  , k·k is any given norm on R n × T , β Σ ( ε ) is d efined by β Σ ( ε ) = λ X t ∈T ′ k w t k 2 2 + X ( i,t ) ∈J ε | f it | , (11) the func tion h : R ≥ 0 → R ≥ 0 is defi ned by ∀ α ∈ R ≥ 0 , h ( α ) = max  α, √ α  (12) and D is g iven as in (7) fr om the norm k·k . Pr oo f. By d efinition (3) of the estimator Ψ , it hold s th at for all ˆ X ∈ Ψ ( Y ) , F ( ˆ X ) ≤ F ( X ) , that is, λ X t ∈T ′ k ˆ x t +1 − A ˆ x t k 2 2 + X t ∈T k y t − C ˆ x t k 1 ≤ λ X t ∈T ′ k x t +1 − Ax t k 2 2 + X t ∈T k y t − C x t k 1 = λ X t ∈T ′ k w t k 2 2 + X t ∈T k f t k 1 . (13) Next, w e der iv e a lower bo und o n the left hand side o f (13). For every t in T , let e t = ˆ x t − x t . Then k ˆ x t +1 − A ˆ x t k 2 2 = k ˆ x t +1 − x t +1 − A ( ˆ x t − x t ) + w t k 2 2 ≥ k e t +1 − Ae t + w t k 2 2 ≥ 1 2 k e t +1 − Ae t k 2 2 − k w t k 2 2 . (14) The last inequ ality uses the identity (see Lem ma 3 in Append ix A for a proof ) k z 1 − z 2 k 2 2 ≥ 1 2 k z 1 k 2 2 − k z 2 k 2 2 ∀ ( z 1 , z 2 ) ∈ R n × R n . (15 ) Similarly , we can write k y t − C ˆ x t k 1 = k y t − C x t − C ( ˆ x t − x t ) k 1 = k f t + C e t k 1 As a co nsequenc e , the seco nd term of th e left-han d-side of (13) is expre ssible as X t ∈T k y t − C ˆ x t k 1 = X ( i,t ) ∈I ×T   f it + c ⊤ i e t   . 1 f it denotes the i -th entry of the ve ctor f t . Now , depend ing on if th e couple ( i, t ) belon gs to J ε or not, we apply the triangle inequality prop erty of the absolute value differently , the two cases being ∀ ( i, t ) ∈ J ε ,   f it + c ⊤ i e t   ≥ | c ⊤ i e t | − | f it | ∀ ( i, t ) ∈ J c ε ,   f it + c ⊤ i e t   ≥ | f it | − | c ⊤ i e t | It follows th at X t ∈T k y t − C ˆ x t k 1 ≥ X ( i,t ) ∈J ε | c ⊤ i e t | − X ( i,t ) ∈J c ε | c ⊤ i e t | − X ( i,t ) ∈J ε | f it | + X ( i,t ) ∈J c ε | f it | . Combining this with (1 3) and (1 4) an d re-arr a nging, yield s λ 2 X t ∈T ′ k e t +1 − Ae t k 2 2 + X ( i,t ) ∈J ε | c ⊤ i e t | − X ( i,t ) ∈J c ε | c ⊤ i e t | ≤ 2  λ X t ∈T ′ k w t k 2 2 + X ( i,t ) ∈J ε | f it |  (16) On the r ight hand side o f (16), w e recogn iz e 2 β Σ ( ε ) as in (10). As to the term o n the lef t hand side, it is equal to H ( E ) − 2 P ( i,t ) ∈J c ε | c ⊤ i e t | . Indepe n dently , |J c ε | = r so by d efinition (8) o f the index p r , X ( i,t ) ∈J c ε | c ⊤ i e t | ≤ p r H ( E ) (17) Consequently , it follows fro m (16) and (1 7) th at (1 − 2 p r ) H ( E ) ≤ H ( E ) − 2 X ( i,t ) ∈J c ε | c ⊤ i e t | ≤ 2 β Σ ( ε ) . Since p r is assumed to be smaller than 1 / 2 , 1 − 2 p r > 0 . Therefo re, we can wr ite H ( E ) ≤ 2 β Σ ( ε ) 1 − 2 p r (18) Thanks to L emma 2, w e h ave H ( E ) ≥ D q ( k E k ) fo r any giv en n orm k·k on R n × T . This implies th at q ( k E k ) ≤ 2 β Σ ( ε ) D (1 − 2 p r ) Now ob serve that th e functio n h defined in (12 ), is the in verse function of q , meaning that for e very λ ∈ R ≥ 0 , h ( q ( λ )) = λ . Moreover , h is an increasin g fun ction. Ap plying h to both members of the previous in equality gives the desired result. The resilienc e p r operty of th e estimator (3) lies here in the fact that, u n der the con d itions of Theo rem 1 , the bo und in (10) on th e estimation erro r does not depend o n the magnitud es of the extreme values o f the noise sequence { f it } ( i,t ) ∈I ×T . Considering in particular the fu n ction β Σ ( ε ) , we remar k tha t it can be overestimated as follows β Σ ( ε ) ≤ λ X t ∈T ′ k w t k 2 2 + |J ε | ε. (19) W e r e cognize two terms in th e upp er bo und of β Σ ( ε ) : (i) the first o ne is a sum which simply represents the uncertain ty broug ht b y the dense noise w t over th e whole state trajectory and which does not depen d on ε ; (ii) the seco n d on e is a bound on th e sum of tho se instances of f it whose magnitud e is smaller that ε . Because β Σ is a fun ction of ε , the bo und in (10) represen ts indeed a family o f bou nds parameter ized b y ε . Since ε is a mere analysis device, a qu estion would b e h ow to select it for the analy sis to a c h ieve the smallest bou nd. Such values, say ε ⋆ , satisfy ε ⋆ ∈ arg min ε ≥ 0  h  2 β Σ ( ε ) D (1 − 2 p r )  : r = |J c ε | , p r < 1 / 2  . Another interesting p o int is that the inequ a lity stated by Theorem 1 h olds fo r any no rm on R n × T . No te th ough that the value of the boun d depend s (throu gh the p arameter D defined in (7)) on the specific norm used to measu r e the esti- mation error . Moreover , different choices of the p erform a nce- measuring norm w ill result in different geo metric forms for the uncertain set, that is, the ball (in the ch osen n orm) centered at the true state with rad ius equal to the u pper bou nd displayed in (1 0). W e also ob serve that the smaller the p arameter p r is, the tigh ter the err o r boun d will be, wh ich sug g ests that the estima to r is mo r e re silient whe n p r is lower . A similar reasoning ap plies to the n umber D which is desired to be large her e. T h ese two parame ter s (i.e. , p r and D ) reflect proper ties of the s ystem whose state is being estimated. The y can be in terpreted, to some extent, as me asures of the d egre e of ob servability of the system. In co nclusion, the estimator inherits partially its resilience proper ty from characteristics of the system being ob served. This is co nsistent with the well-known fact that the mo re obser vable a system is, the mor e ro bustly its state can be estimated from outp ut measuremen ts. Finally , an interesting property of the estimator can b e stated in the a bsence of den se n oise: Corollary 1. Consider the system Σ defi ned by ( 1) and let r = |J c 0 | (which mea ns that we consider every n o nzer o occurr ence o f f it as an outlier). If p r < 1 / 2 , an d if w t = 0 for all t , th en the estimator d e fined by (3) r etrieves exactly the state trajectory of the system. Pr oo f. Th is follows directly fr o m the fact that β Σ (0) = 0 in the case where the re is no d ense noise w t and ε = 0 . Therefo re, we have the exact recoverability of every state of the system (1) by the e stimator wh en there is no pro cess noise. Accord ing to o ur analysis, the num ber of outliers that can be h andled by the estimator in this case can be undere stima ted by max  r : p r < 1 / 2  . (20) I V . S I M U L A T I O N R E S U LT S In this section, we p r esent th e simulation results of a system desgin ed as (1) with A =  − 0 . 11 − 0 . 34 − 0 . 34 0 . 46  , C =  1 . 4 − 0 . 94  20 40 60 80 100 − 2 0 2 sample t x 1 t x t ˆ x t ˆ x L t 20 40 60 80 100 − 2 0 2 sample t x 2 t x t ˆ x t ˆ x L t Fig. 1: State of the system a n d its estimates (resilient estimator and smoo ther) in absence of sparse noise 20 40 60 80 100 − 20 0 20 sample t x 1 t x t ˆ x t ˆ x L t 20 40 60 80 100 − 10 0 10 sample t x 2 t x t ˆ x t ˆ x L t 20 40 60 80 100 − 50 0 50 sample t y t s t y wt Fig. 2: State, estima te d states (thr o ugh r esilient estimation and smooth ing) and outpu t of th e system in pr esence o f sparse no ise w t is a gaussian white n o ise of un it variance. T he den se compon ent of f t , wh ich will b e called v t in this section, is a gaussian white noise of signal-to-noise r atio e qual to 30dB, while the sparse c o mpon ent of f t , which will be called s t , is a sparse vector wh ose non-zero elemen ts are random ly selected an d given a r andom value: as a result of this structu re, we note y wt = C x t + v t the u ncorru pted output of the system. Th e estimated states were the n obtained b y directly solv ing the optimisation prob lem defined in (3) with λ = 1 / 5 throug h CVX [10]. T o give a basis for comparison , we also estimated the state o f the system throu gh a Rauch - T un g-Striebel smoother which is an extension of the Kalman filter to o ffline estimation [9]. Figure 1 presen ts the classic case where there is n o sp a rse noise cor r upting the output of the system. This is the scenario handled by classic estimato r s such as the Kalman Filter o r in our case the Rauch-T ung-Striebel s moother . W e can howe ver notice that ou r estimator giv es satisfying results, fitting th e trajectory of the real state an d g iving very similar resu lts to the smooth er . I t is all the m ore interesting as ou r estimator does not take into acc ount th e statistical prop e rties of the noises inv o lved in th e system, co ntrary to th e smooth e r which requires a tunin g to appro a c h the variance of tho se noises. Figure 2 now pre sen ts the case where twenty corr u pted values were added to the outpu t of the system. The smoother tries to comp ensate the attacks, as it can b e n oted that the estimate di verges when a cor r uption occu rs, but it is entirely normal giv en th at the Kalman filter theor y is de sig n ed around noises in the fo rm of wh ite g a ussian pr ocesses on ly . Figure 3 compa r es the trajectory o f the real state and the estimated state obtained thr ough our resilient estimator . Even in the p resence of corru pted measuremen ts of arb itrarily large magnitud e, the estimator still manages to e fficiently track the trajectory of the r e a l states, sho wing that its perfo rmance are not really degrad e d in tha t case. 20 40 60 80 100 − 2 0 2 sample t x 1 t x t ˆ x t 20 40 60 80 100 − 2 0 2 sample t x 2 t x t ˆ x t Fig. 3: State of the system a n d its estimate (r esilient estimator) in presen ce of sparse noise V . C O N C L U S I O N In this pa p er , we considered the p r oblem of estimating the state o f line a r d iscr e te-time systems in the face of uncertainties mo deled as process an d measuremen t noise in the system eq uations. The mea su rement noise sequence assumes v alues o f po ssibly ar bitrarily large amplitude w h ich occur in termittently in time. For this prob lem we p r oposed an e stima to r b ased on the resolution of a conv ex op timization problem . In p articular, we proved a resilienc e prop erty for the proposed estimator, that is, the resulting estimation error is bou n ded by a b ound which is indepen dent of the extrem e values of the measur e m ent no ise provided that the num ber of o ccurren ces (over time and over th e wh o le set of sensors) of such extreme v alues is limited . Future works will aim at genera lizin g the resilient pro perties to a wider class of estimators and applying the estimation framew ork to relev ant practical cases. A P P E N D I X A. Additiona l elements to the pr oo f of Theorem 1 Lemma 3. Let G : R n × m → R ≥ 0 be a con vex function satisfying the pr o perties of positive defi niteness and r e laxed homogeneity ( for a given K ∞ function σ ) a s both defined in Lemma 1. Then , for all ( S 1 , S 2 ) ∈ R n × m × R n × m , G ( S 1 − S 2 ) ≥ 2 σ (1 / 2) G ( S 1 ) − G ( S 2 ) (21) Pr oo f. As G is convex, G  1 2 ( S 1 − S 2 ) + S 2 2  ≤ 1 2 G ( S 1 − S 2 ) + 1 2 G ( S 2 ) (22) which, by m ultiplying the who le ine quality by 2 , ca n b e rewritten as G ( S 1 − S 2 ) ≥ 2 G ( S 1 / 2) − G ( S 2 ) (23) Moreover , by assump tion, G verifies th e relaxed homogen e- ity pro perty with a K ∞ function σ : it entails th at ∀ S 1 ∈ R n × m , G ( S 1 / 2) ≥ σ (1 / 2 ) G ( S 1 ) (24) which, when in je c te d in (23), gives the desire d result. In the case w h ere G = k·k 2 2 , as nor ms are ho m ogeneo us, for every λ ∈ R > 0 and z ∈ R n , G ( z ) = G ( λz ) /λ 2 . It follows that Lemma 3 can be ap plied to G for σ such that ∀ α ∈ R ≥ 0 , σ ( α ) = α 2 , yielding k z 1 − z 2 k 2 2 ≥ 1 2 k z 1 k 2 2 − k z 2 k 2 2 ∀ ( z 1 , z 2 ) ∈ R n × R n . (25 ) R E F E R E N C E S [1] M. ApS. The MOSEK optimization tool box for MA T LAB. [2] L. Bako. On a class of optimization-b ased robust estimat ors. IEEE T ransactions on Automatic Contr ol , 62(11):5990– 5997, 2017. [3] E. J. Candes. T he restricte d isometry property a nd its impl icati ons for compressed sensing. Comptes r endus mathe matique , 346(9-10):589 – 592, 2008. [4] E. J. Candès and M. B. W akin. An introductio n to compressi ve sampling. IEEE Signal Proc essing Society , 25:21–30, 2008. [5] A. Cardenas, S. Am in, and S. Sastry . Se cure control: T o ward s survi v able cyber -physical systems. In International Confer ence on Distrib uted Computing Systems W orkshops, Beijing, China , pages 495–500, 2008. [6] Y . H. Chang, Q. Hu, and C. J. T omlin. Sec ure estimation based kalman filter for cybe r–physica l systems against sensor attac ks. Automatica , 95:399–41 2, 2018. [7] H. Fawz i, P . T abuada , and S. Diggavi . Secure estimatio n and control for cybe r-physi cal s ystems under adversari al atta cks. IEEE T ransactions on Automatic Contr ol , 59(6):1454–1 467, 2014. [8] S. Foucart and H. Rauhut. A mathemati cal intr oduct ion to compressive sensing . Birkhäuser , 2013. [9] A. Gelb . Applied optimal estimation . MIT press, 1974. [10] M. C. Grant and S. P . Boyd. CVX: Matlab software for discipline d con vex programming, version 2.1. 2017. [11] D. Han, Y . Mo, and L. Xie. Con ve x optimizat ion based state estimat ion against sparse inte grity attacks. IEEE T ransac tion on Automatic Contr ol (DOI: 10.1109/T AC.2019.2891458) , 2019. [12] C. M. Kell ett. A compendi um of comparison function results. Mathemat ics of Contr ol, Signal s, and Syste ms , 26:339–374, 2014. [13] S. Mishra, Y . Shoukry , N. Karamchandani, S. N. Digga vi, and P . T ab uada. Secure state estimation against sensor attacks in the presence of noise. IEEE T ransact ions on Contr ol of Network Syste ms , 4(1):49–59 , 2017. [14] M. Pajic , I. Lee , and G. J. Pappas. Attack-resil ient state estimati on for noisy dynami cal systems. IEE E T ransactio ns on Contr ol of Netwo rk Systems , 4(1):82–92, 2017. [15] F . Pasqual etti, F . Dorfler , and F . Bullo. Attack detec tion and identi- fication in cyber -physica l systems. IEE E T ransaction s on Automati c Contr ol , 58(11):271 5–2729, 2013. [16] R. T . Rockafella r . Conv ex Analysis . Princeton Univ ersity Press. [17] Y . Sharon, J. Wright, and Y . Ma. Minimum s um of distances estimator : Robust ness and stabi lity . In American Contr ol Confer ence , St. Louis, MO, USA , pages 524–530, 2009. [18] Y . Shoukry and P . T abuada. Event-t riggere d s tate observe rs for sparse sensor noise /atta cks. IEEE T ransact ions on Automatic Contr ol , 61(8):2079 –2091, 2016. [19] J. F . Sturm. Using SeDuMi 1.02, a MA TLAB toolbox for optimiz ation ov er symmetric cones. Optimiza tion methods and software , 11(1- 4):625–653 , 1999.