Robust stability of moving horizon estimation for nonlinear systems with bounded disturbances using adaptive arrival cost

IET Research Journals Rob ust stability of mo ving horizon estimation f or nonlinear systems with bounded disturbances using adaptive arriv al cost ISSN 1751-8644 doi: 0000000000 www .ietdl.org Nestor N. Deniz 1 , Marina H. Mur illo 1 , Guido Sanchez 1 , Lucas M. Genzelis 1 , Leonardo L. Giov anini 1 1 Research institute for signals, systems and computational intelligence, Ciudad Universitaria UNL, Ruta Nac. No 168, km 472.4, FICH, 4to Piso (3000) Santa Fe - Argentina * E-mail: ndeniz@sinc.unl.edu.ar Abstract: In this paper , the robust stability and conv ergence to the true state of moving horizon estimator based on an adaptiv e arrival cost are established for nonlinear detectable systems . Robust global asymptotic stability is shown for the case of non- vanishing bounded disturbances whereas the convergence to the true state is prov ed for the case of vanishing disturbances. Se veral simulations were made in order to show the estimator behaviour under different operational conditions and to compare it with the state of the ar t estimation methods. 1 Introduction State estimation plays a fundamental role in feedback control, sys- tem monitoring and system optimization because noisy measure- ments is the only information av ailable from the system. Sev eral methods have been developed for accomplishing such task (see [1]; [2]; among others). All these methods ha ve been dev eloped upon the assumption on the knowledge of noises and model of the system, as well as, the absence of constraints. In practice, these assumptions are not easily satisfied and research efforts were focused on approaches that do not relay on such require- ments (see [3], [4], [5], among others). For example, an H ∞ filter is designed minimizing the H ∞ norm of the mapping between dis- turbances and estimation error. In [3], [6] and [7] an approach that solves a least-square estimation problem is introduced. Both meth- ods are based on the adequate selection of the uncertainty model instead of relying on statistical assumptions on noises. In these approaches, uncertainty models are formulated based on the avail- able information of the system. In the same way , robust estimation algorithms based on as min-max robust filtering, set-valued estima- tion and guaranteed cost paradigm, hav e attracted the attention of the research community (see [4], [8]). Building on the success on moving horizon control, moving hori- zon estimation (MHE) has attracted attention of researchers since the pioneering w ork of [9] (see also [10], [11] and [12]). The interest in such estimation methods stems from the possibility of dealing with limited amount of data, instead of using all the infor- mation av ailable from the beginning, and the ability to incorporate constraints. In recent years, both theoretical properties of various MHE schemes as well as efficient computational methods for real- time implementation hav e been studied (see [13], [14], [15], [16], [17], [18]). In particular , it is of interest to establish robust stability and estimate conv ergence properties. In recent years sev eral results hav e been obtained for different algorithms, advancing from ide- alistic assumptions (observability and no disturbances) to realistic situations (detectability and bounded disturbances). For nonlinear observable systems, [12] established the asymp- totic stability of the estimation error for the standard cost function. Furthermore, if the disturbances are asymptotically vanishing the estimation error is robust asymptotically stable and it asymptot- ically conv erges to zero ([19]-[20]). [14] and [21] proposed an estimation scheme, based on least-square cost function of the esti- mation residuals, that guaranteed the boundedness of estimation error for observable systems subject to bounded additiv e distur- bances. Finally , for the general case of nonlinear detectable systems subject to bounded disturbances, [22] and [23] showed the robust global asymptotic stability (RGAS) and conv ergence of estimation error in case of bounded or vanishing disturbances, respectively . In these works, the least-square objective function was modified by adding a max-term. [22] established RGAS for the full informa- tion estimator while [23] established RGAS and con vergence for the moving horizon estimator . Furthermore, for a particular choice of the weights of the objectiv e function, [23] established these results for the least-squares type objectiv e function. This paper introduces the RGAS and con vergence analysis for the moving horizon estimator based on adaptive arri val cost proposed in [18] in the practical case of nonlinear detectable systems subject to bounded disturbances. T o establish robust stability properties for MHE it is crucial that the prior weighting in the cost function is chosen properly . In various schemes the necessary assumptions in the prior weighting are difficult to verify ([12], [19]), while in others can be verified a prior [23]. In the MHE scheme analysed in this work, the assumption on the prior weighting can be verified a prior by design. Furthermore, the disturbances gains become uniform (i.e., they are valid independent of N ), allo wing to extend the stability analysis to full information estimators with least-square type cost functions. The rest of the paper is or ganized as follo ws: Section 2 introduces the notation, definitions and properties that will be used through the paper . Section 3 presents the main result and shows its connections with previous stability analysis. Section 4 discusses simple exam- ples, pre viously used in the literature, with the purpose of illustrating the concepts and also in order to show the difference with others MHE algorithms. Finally , Section 5 presents conclusions. 2 Preliminaries and setup 2.1 Notation Let Z [ a,b ] denotes the set of integers in the interval [ a, b ] ⊆ R , and Z ≥ a denotes the set of integers greater or equal to a . Bold- face symbols denote sequences of finite or infinite length, i.e., w : = { w k 1 , . . . , w k 2 } for some k 1 , k 2 ∈ Z ≥ 0 and k 1 < k 2 , respectively . W e denote x j | k as the finite sequence x giv en at time k ∈ Z ≥ 0 and j ∈ [ k 1 , k 2 ] . By | x | we denote the Euclidean norm of IET Research Journals, pp. 1–10 c  The Institution of Engineering and T echnology 2015 1 a vector x ∈ R n . Let k x k : = sup k ∈ Z ≥ 0 | x k | denote the supreme norm of the sequence x and k x k [ a,b ] : = sup k ∈ Z [ a,b ] | x k | . A func- tion γ : R ≥ 0 → R ≥ 0 is of class K if γ is continuous, strictly increasing and γ (0) = 0 . If γ is also unbounded, it is of class K ∞ . A function ζ : R ≥ 0 → R ≥ 0 is of class L if ζ ( k ) is non increas- ing and lim k →∞ ζ ( k ) = 0 . A function β : R ≥ 0 × Z ≥ 0 → R ≥ 0 is of class KL if β ( · , k ) is of class K for each fixed k ∈ Z ≥ 0 , and β ( r, · ) of class L for each fixed r ∈ R ≥ 0 . The follo wing inequalities hold for all β ∈ KL , γ ∈ K and a j ∈ R ≥ 0 with j ∈ Z [1 ,n ] γ ( a 1 + a 2 + . . . + a n ) ≤ γ ( na 1 ) + . . . + γ ( na n ) , β ( a 1 + a 2 + . . . + a n , k ) ≤ β ( na 1 , k ) + . . . + β ( na n , k ) . (1) The preceding inequalities hold since max { a j } is included in the sequence { a 1 , a 2 , . . . , a n } and K functions are non-negati ve strictly increasing functions. Bounded sequences: A sequence w is bounded if k w k is finite. The set of bounded sequences w is denoted as W ( w max ) : = { w : k w k ≤ w max } for some w max ∈ R ≥ 0 Con vergent sequences: A bounded infinite sequence w is con- ver gent if | w k | → 0 as k → ∞ . Let denote the set of conv ergent sequences C : C w : = { w ∈ W ( w max ) | w is con vergent } Analogously , C v is defined for the sequence v . 2.2 Problem statement Let us consider the state estimation problem for nonlinear discrete time systems of the form x k +1 = f ( x k , w k ) , x 0 = x 0 y k = h ( x k ) + v k , (2) where x k ∈ X ⊆ R n , w k ∈ W ⊆ R p , y k ∈ Y ⊆ R m , v k ∈ V ⊆ R m are the state, process noise, measurement and estimation resid- uals vectors, respectiv ely . The process disturbance w k and estima- tion residuals v k are unkno wn but assumed to be bounded, i.e, w ∈ W ( w max ) and v ∈ V ( v max ) for some w max , v max ∈ R ≥ 0 . X , Y , W and V are compact and conv ex sets with the null vector 0 belongs to them. In the following we assume that f : R n × R p → R n is continuous, locally Lipschitz on x k and h : R n → R m is continuous. The solution to the system (2) at time k is denoted by x ( k ; x 0 , w ) , with initial condition x 0 and process disturbance sequence w . Furthermore, the initial condition x 0 is unknown, but a prior knowledge ¯ x 0 is assumed to be available and its error is assumed to be bounded, i.e., ¯ x 0 ∈ X 0 : = { ¯ x 0 : | x 0 − ¯ x 0 | ≤ e max } , X 0 ⊆ X . The solution of the estimation problem aims to find at time k an estimate ˆ x k | k of the current state x k minimizing a performance metric using by the MHE. At each sampling time k , given the previous N measurements y := { y k − N , . . . , y k − 1 } , the following optimization problem is solved min ˆ x k − N | k , ˆ w j | k Ψ : = Γ k − N | k  ˆ x k − N | k  + k − 1 P j = k − N `  ˆ w j | k , ˆ v j | k  s.t.        ˆ x j +1 | k = f  ˆ x j | k , ˆ w j | k  , j ∈ Z [ k − N ,k − 1] y j = h  ˆ x j | k  + ˆ v j | k , ˆ x j | k ∈ X , ˆ w j | k ∈ W , ˆ v j | k ∈ V , (3) where ˆ x k − j | k is the optimal estimated and ˆ w j | k is the optimal process noise estimate at sample k − j j = 0 , 1 , . . . , N based on measurements y k − j av ailable at time k . The process noise ˆ w j | k : = { ˆ w k − N | k , . . . , ˆ w k − 1 | k } and ˆ x k − N | k are the optimiza- tion v ariables. The stage cost `  w j | k , v j | k  penalizes the estimated process noise sequence ˆ w j | k and the estimation residuals ˆ v j | k = y j − h  ˆ x j | k  , while Γ k − N  ˆ x k − N | k  penalizes the prior esti- mated ˆ x k − N | k . The adequate choice of ` ( · ) and Γ k − N ( · ) , and their parameters, allows to ensure the robust stability of the estimator [23]. While the estimation window is not full, k ≤ N , problem (3) can be reformulated and solved as a full information problem min ˆ x k − N | k , ˆ w Ψ : = Γ 0 | k  ˆ x 0 | k  + k − 1 P j =0 `  ˆ w j | k , ˆ v j | k  s.t.        ˆ x j +1 | k = f  ˆ x j | k , ˆ w j | k  , j ∈ Z [0 ,k − 1] y j = h  ˆ x j | k  + ˆ v j | k , ˆ x j | k ∈ X , ˆ w j | k ∈ W , ˆ v j | k ∈ V , as k increases this problem becomes (3) for all k ≥ N . In pre vious works, the rob ust stability of MHE has been achie ved by modifying the standard least-square cost function through the inclusion of a max –term ([22]; [23]) or by a suitable choice of the cost’ s function parameters ([23]). Another mechanism to solve this problem is combining a suitable choice of the stage cost `  ˆ w j | k , ˆ v j | k  with a time–varying prior weight of the form Γ k − N | k  ˆ x k − N | k  = k ˆ x k − N | k − ¯ x k − N k P − 1 k − N | k , (4) whose parameters  P − 1 k − N | k , ¯ x k − N  are recursiv ely updated using the information available at time k ([18], 2017). The prior weighting is defined in this way to av oid the introduction of artificial cycling in the estimation process (see [19]). In this approach, the prior weight matrix P k − N | k is giv en by  k − N = y k − N − ˆ y k − N | k , N k = h 1 + ˆ x T k − N | k − 1 P k − N − 1 ˆ x k − N | k − 1 i σ |  k − N | 2 2 α k = 1 − 1 N k , W k = " I − P k − N − 1 ˆ x k − N | k − 1 ˆ x T k − N | k − 1 1 + ˆ x T k − N | k − 1 P k − N − 1 ˆ x k − N | k − 1 # P k − N − 1 , P k − N =  1 α k W k if 1 α k T r ( W k ) ≤ c, W k otherwise , (5) where σ, σ w , c, λ ∈ R > 0 , c > λ, P 0 = λI n × n and σ  σ w , where σ w denotes the process noise variance. The prior knowledge of the window ¯ x k − N is updated using a smoothed estimate ([24]) ¯ x k − N = ˆ x k − N | k − 1 . (6) The optimization problem (3) can be reformulated in terms of the ini- tial condition ˆ x 0 and the estimated process noises and the residuals along the entire trajectory as follows min ˆ x 0 | k , ˆ w Ψ : = k − 1 P j = k − N `  ˆ w j | k , ˆ v j | k  + k − N − 1 P j =1 α k − N − j k `  ˆ w j | k , ˆ v j | k  + α k − N k Γ 0 | k  ˆ x 0 | k  s.t.        ˆ x j +1 | k = f  ˆ x j | k , ˆ w j | k  , j ∈ Z [0 ,k − 1] , α ∈ (0 , 1] y j = h  ˆ x j | k  + ˆ v j | k , ˆ x j | k ∈ X , ˆ w j | k ∈ W , ˆ v j | k ∈ V , This formulation of problem (3) allo ws to explicitly see the ef fect of past data on the current state estimate ˆ x k | k . In this formulation IET Research Journals, pp. 1–10 2 c  The Institution of Engineering and T echnology 2015 it is easy to see the e xponential averaging of these data. Allowing α change in time, the past data has different affects on the current estimates depending on ˆ x k | k . Before proceeding to the development of the main results, we state the main properties and assumptions about the prior weighting Γ k − N . The updating mechanism (5) is a time-v arying filter whose inputs are ˆ x k − N | k − 1 ˆ x T k − N | k − 1 and the initial condition P 0 . It generates recursively a real-time estimation of P k − N | k by updating P k − N − 1 | k − 1 with an exponential time-av eraging of ˆ x k − N | k − 1 ˆ x T k − N | k − 1 . The updating mechanism (5) only use data and it does not rely on a model of the system. The sequence P k | k k ≥ 0 is positiv e definite, it is decreasing in norm and it is bounded. The proof of these properties follo ws similar steps as in [18]. Assumption 1. The prior weighting Γ k − N is a continuous function Γ k − N : R n → R lower bounded by γ p ∈ K ∞ and upper bounded by ¯ γ p ∈ K ∞ such that: γ p  | ˆ x k − N | k − ¯ x k − N |  ≤ Γ k − N  ˆ x k − N | k  Γ k − N  ˆ x k − N | k  ≤ ¯ γ p  | ˆ x k − N | k − ¯ x k − N |  (7) for all ˆ x ∈ X and γ p ( r ) ≥ c p r a , ¯ γ p ( r ) ≤ ¯ c p r a . (8) wher e 0 ≤ c p ≤ ¯ c p and a ∈ R ≥ 1 . Giv en prior weighting updating scheme (5) inequality (7) satisfies [18] | P − 1 0 | r a ≤ Γ k − N  ˆ x k − N | k  ≤ | P − 1 ∞ | r a . (9) Definition 1. The system (2) is incrementally input/output-to-state stable if ther e exist functions β ∈ K L and γ 1 , γ 2 ∈ K such that for every two initial states z 1 , z 2 ∈ R n , and any two disturbances sequences w 1 , w 2 the following holds for all k ∈ Z ≥ 0 : | x ( k , z 1 , w 1 ) − x ( k , z 2 , w 2 ) | ≤ max { β ( | z 1 − z 2 | , k ) , γ 1 ( k w 1 − w 2 k ) , γ 2 ( k h ( x 1 ) − h ( x 2 ) k ) } ≤ β ( | z 1 − z 2 | , k )+ γ 1 ( k w 1 − w 2 k ) + γ 2 ( k v 1 − v 2 k ) (10) This definition combines the concepts of output-to-state-stability (OSS) and input-to-state-stability (ISS). As stated in [25], the notion of IOSS represents a natural combination of the ideas of str ong observability and ISS, and it was called detectability in [26] and str ong unboundedness observability in [27]. In addition, the exis- tence of an observ er for the system (2), which is incrementally input- output-to-state stable ( i -IOSS) instead of IOSS (see Remark 24 in [25]), is assumed. Note that k h ( x 1 ) − h ( x 2 ) k = k v 1 − v 2 k , since y k = h ( x k ) + v k These assumptions will help us to bound the functions in volv ed in the definition of i -IOSS and to relate them with the terms of the MHE cost function (stage cost and prior weight). In the following sections the updating mechanism (5) and the assumption of i -IOSS [28] will be used to prove robust stabil- ity of the proposed MHE in the presence of bounded disturbances and con vergence to the true state in the case of conv ergent distur- bances. Some assumptions about functions related to system (2) and Definition 1 will be helpful in the sequel. Assumption 2. The function β ( r , s ) ∈ K L and satisfies the follow- ing inequality β ( r, s ) ≤ c β r p s − q (11) for some c β ∈ R ≥ 0 , p ∈ R ≥ 0 and q ∈ R ≥ 0 and q ≥ p . Assumption 3. The stage cost ` ( · , · ) : R p × R m → R is a con- tinuous function bounded by γ w , γ v , ¯ γ w , ¯ γ v ∈ K ∞ such that the following inequalities ar e satisfied ∀ w ∈ W and v ∈ V γ w ( w ) + γ v ( v ) ≤ ` ( w, v ) ≤ ¯ γ w ( w ) + ¯ γ v ( v ) . (12) Functions γ 1 and γ 2 from Definition 1 are related with the bounds of stage cost ¯ γ w , γ w , ¯ γ v and γ v through the following inequalities γ 1  3 γ − 1 w ( r )  ≤ c 1 r α 1 , γ 2  3 γ − 1 v ( r )  ≤ c 2 r α 2 (13) for c 1 , c 2 , α 1 , α 2 > 0 . Inequalities (11) to (13) were used in previ- ous works ([22]; [23]). In this work, we claim that the proposed estimator holds the prop- erty of being robust global asymptotic stable, which is defined as follows. Definition 2. Consider the system described (2) subject to dis- turbances w ∈ W ( w max ) and v ∈ V ( v max ) for w max ∈ R ≥ 0 , v max ∈ R ≥ 0 with prior estimate ¯ x 0 ∈ X ( e max ) for e max ∈ R ≥ 0 . The moving horizon state estimator given by equation (3) with adap- tive prior weight is r obustly globally asymptotically stable (RGAS) if ther e exists functions Φ ∈ KL and π w , π v ∈ K such that for all x 0 ∈ X , all ¯ x 0 ∈ X 0 , the following is satisfied for all k ∈ Z ≥ 0 | x k − ˆ x k | ≤ Φ ( | x 0 − ¯ x 0 | , k ) + π w  k w k [0 ,k − 1]  + π v  k v k [0 ,k − 1]  . (14) W e want to show that if system (2) is i-IOSS, then Assumptions 1, 2 and 3 are fulfilled and the proposed MHE estimator with adap- tiv e arriv al cost weight matrix is RGAS. Furthermore, if the process disturbance and measurement noise sequences are conv ergent (i.e., w ∈ C w , v ∈ C v ), then ˆ x k | k → x k as k → ∞ . 3 Robust stability of moving horizon estimation under bounded disturbances W e are ready to deri ve the main result: RGAS of the proposed mo v- ing horizon estimator with a large enough estimation horizon N for nonlinear detectable systems under bounded disturbances. Further- more, a KL function e xist such that (14) is v alid with this Φ , π w and π v for all estimation horizon N ≥ N . Theorem 1. Consider an i-IOSS system (2) with disturbances w ∈ W ( w max ) , v ∈ V ( v max ) . Assume that the arrival cost weight matrix of the MHE pr oblem Γ k − N is updated using the adaptive algorithm (5) . Mor eover , Assumptions 1 , 2 and 3 are fulfilled and initial condition x 0 is unknown, but a prior estimate ¯ x 0 ∈ X 0 is available. Then, the MHE estimator (3) is RGAS . Proof . The optimal cost of problem (3) is giv en by Ψ ∗ N = Ψ  ˆ x ∗ k − N | k , ˆ w ∗ [ k − N , k − 1]  = Γ k − N  ˆ x ∗ k − N | k  + k − 1 X j = k − N `  ˆ w ∗ j | k , ˆ v ∗ j | k  , IET Research Journals, pp. 1–10 c  The Institution of Engineering and T echnology 2015 3 which is bounded (Assumptions 1 and 3) ∀ | ˆ w j | k | and ∀ | ˆ v j | k | for all j ∈ Z [ k − N ,k − 1] by Ψ ∗ N ≤ γ p  | ˆ x ∗ k − N | k − ¯ x k − N |  + N γ w  | ˆ w ∗ j | k |  + N γ v  | ˆ v ∗ j | k |  , Ψ ∗ N ≥ γ p  | ˆ x ∗ k − N | k − ¯ x k − N |  + N γ w  | ˆ w ∗ j | k |  + N γ v  | ˆ v ∗ j | k |  . Due optimality , the following inequalities hold ∀ k ∈ [ k − N , k − 1] Ψ  ˆ x ∗ k − N | k , ˆ w ∗  ≤ Ψ ( x k − N , w ) , ≤ ¯ γ p ( | x k − N − ¯ x k − N | ) + N ¯ γ w ( k w k ) + N ¯ γ v ( k v k ) , (15) then, taking into account the lower and upper bounds we ha ve | ˆ x k − N | k − ¯ x k − N | ≤ γ − 1 p ( ¯ γ p ( | x k − N − ¯ x k − N | ) + N ¯ γ w ( k w k ) + N ¯ γ v ( k v k )) . By mean of Assumptions 1 and 3 the last inequality can be written as follows | ˆ x k − N | k − ¯ x k − N | ≤ γ − 1 p (3 ¯ γ p ( | x k − N − ¯ x k − N | )) + γ − 1 p (3 N ¯ γ w ( k w k )) + γ − 1 p (3 N ¯ γ v ( k v k )) , ≤ 3 1 a | P − 1 0 |  | P − 1 ∞ | 1 a | x k − N − ¯ x k − N | + N 1 a ¯ γ 1 a w ( k w k ) + N 1 a ¯ γ 1 a v ( k v k )  . Analogously , bounds for | ˆ v j | k | and | ˆ w j | k | can be found | ˆ w j | k | ≤ γ − 1 w  3 N ¯ γ p ( | x k − N − ¯ x k − N | )  + γ − 1 w (3 ¯ γ w ( k w k )) + γ − 1 w (3 ¯ γ v ( k v k )) , | ˆ v j | k | ≤ γ − 1 v  3 N ¯ γ p ( | x k − N − ¯ x k − N | )  + γ − 1 v (3 ¯ γ w ( k w k )) + γ − 1 v (3 ¯ γ v ( k v k )) . (16) Next, let us consider some sample k ∈ Z ≥ N . Assuming that system (2) is i-IOSS with z 1 = x k − N , z 2 = ˆ x k − N | k , w 1 = { w j } , w 2 = { ˆ w j | k } , v 1 = { v j } and v 2 = { ˆ v j | k } for all j ∈ Z [ k − N ,k − 1] . Since x ( k ) = x ( N , z 1 , w 1 ) , ˆ x ( k ) = ˆ x k | k = x ( N , z 2 , w 2 ) we obtain | x k − ˆ x k | k | ≤ β  | x k − N − ˆ x k − N | k | , N  + γ 1     w j − ˆ w j | k     + γ 2     v j − ˆ v j | k     . (17) In order to get a finite upper bound for the estimation error, the three terms in the right hand side of equation (17) must be upper bounded. The first term can be written β  | x k − N − ˆ x k − N | k | , N  ≤ β (2 | x k − N − ¯ x k − N | , N ) + β  2 | ˆ x k − N | k − ¯ x k − N | , N  ≤ β (2 | x k − N − ¯ x k − N | , N ) + β 2 3 1 a | P − 1 ∞ | 1 a | P − 1 0 | | x k − N − ¯ x k − N | + 2 3 1 a N 1 a | P − 1 0 | ¯ γ 1 a w ( k w k ) + 2 3 1 a N 1 a | P − 1 0 | ¯ γ 1 a v ( k v k ) , N ! ≤ β (2 | x k − N − ¯ x k − N | , N ) + β 6 3 1 a | P − 1 ∞ | 1 a | P − 1 0 | | x k − N − ¯ x k − N | , N ! + β 6 3 1 a N 1 a | P − 1 0 | ¯ γ 1 a w ( k w k ) , N ! + β 6 3 1 a N 1 a | P − 1 0 | ¯ γ 1 a v ( k v k ) , N ! . Using Assumptions 1 and 2, function β ( · ) is bounded by β  | x k − N − ˆ x k − N | k | , N  ≤ c β 2 p N q | x k − N − ¯ x k − N | p + c β 6 p 3 p a | P − 1 ∞ | p a | P − 1 0 | p N q | x k − N − ¯ x k − N | p + c β 6 p 3 p a N p a | P − 1 0 | ¯ γ p a w ( k w k ) + c β 6 p 3 p a N p a | P − 1 0 | ¯ γ p a v ( k v k ) ≤ c β 2 p N q | x k − N − ¯ x k − N | p + | P − 1 ∞ | | P − 1 0 | ! p c β 6 p 3 p a N q | x k − N − ¯ x k − N | p + c β 6 p 3 p a N p a − q | P − 1 0 | ¯ γ p a w ( k w k ) + c β 6 p 3 p a N p a − q | P − 1 0 | ¯ γ p a v ( k v k ) . T aking in account that P − 1 k is a symmetric positive definite matrix for all k ∈ Z [0 , ∞ ) , then | P − 1 k | ≤ λ max  P − 1 k  , where λ max  P − 1 k  denotes the maximal eigenv alue of matrix P − 1 k . Denoting λ min  P − 1 k  as the minimal eigen value of matrix P − 1 k and taking in account that | P − 1 k | ≤ | P − 1 k +1 | , the maximum con- ditioning number of matrix P − 1 k can be defined as C P − 1 : = λ max  P − 1 ∞  /λ min  P − 1 0  , then β  | x k − N − ˆ x k − N | k | , N  can be bounded by β  | x k − N − ˆ x k − N | k | , N  ≤ c β 18 p | P − 1 0 |  ¯ γ p a w ( k w k ) + ¯ γ p a v ( k v k )  + c β N q  2 p + C p P − 1 18 p  | x k − N − ¯ x k − N | p . (18) The first term in the right side of this equation is bounded due the assumption that | x k − N − ¯ x k − N | ∈ X 0 ( e max ) , while the second term are finite constants. T o extend the validness of (18) to the full estimation horizon, an extension of the function β at the beginning of the estimation, N = 0 , is required. The second term in the right hand side of equation (17), can be bounded by the following inequality ∀ j ∈ Z [ k − N ,k − 1] IET Research Journals, pp. 1–10 4 c  The Institution of Engineering and T echnology 2015 γ 1     w j − ˆ w j | k     ≤ γ 1  k w k +    ˆ w j | k     ≤ γ 1  k w k + γ − 1 w  3 N ¯ γ p ( | x k − N − ¯ x k − N | )  + γ − 1 w (3 ¯ γ w ( k w k )) + γ − 1 w (3 ¯ γ v ( k v k ))  . Recalling Assumption 3, the reader can verify the following inequality γ 1     w j − ˆ w j | k     ≤ c 1 3 α 1 | P − 1 ∞ | α 1 N α 1 | x k − N − ¯ x k − N | aα 1 + c 1 3 α 1 ¯ γ α 1 v ( k v k ) + γ 1  3  k w k + γ − 1 w (3 ¯ γ w ( k w k ))  . (19) In an equivalent manner , a bound for the third term in the right hand side of equation (17) can be found γ 2     v j − ˆ v j | k     ≤ c 2 3 α 2 | P − 1 ∞ | α 2 N α 2 | x k − N − ¯ x k − N | aα 2 + c 2 3 α 2 ¯ γ α 2 w ( k w k ) + γ 2  3  k v k + γ − 1 v (3 ¯ γ v ( k v k ))  . (20) Once an upper bound for the three terms of equation (17) were found, defining ζ : = max { p, aα 1 , aα 2 } , η : = min { q , α 1 , α 2 } and ρ : = max { p, α 1 , α 2 } , equation (17) can be rewritten as follo ws | x k − ˆ x k | k | ≤ | x k − N − ¯ x k − N | ζ N η  C ρ P − 1  c β 18 p + c 1 3 α 1 λ α 1 min  P − 1 0  + c 2 3 α 2 λ α 1 min  P − 1 0  + c β 2 p  +   c β 18 p ¯ γ p a w ( k w k ) | P − 1 0 | + γ 1  3  k w k + γ − 1 w (3 ¯ γ w ( k w k ))  + c 2 3 α 2 ¯ γ α 2 w ( k w k )  +   c β 18 p ¯ γ p a v ( k v k ) | P − 1 0 | + γ 2  3  k v k + γ − 1 v (3 ¯ γ v ( k v k ))  + c 1 3 α 1 ¯ γ α 1 v ( k v k )  . (21) Defining the functions ¯ β ( r, s ) , φ w ( r ) and φ v ( r ) for all r ≥ 0 and s ∈ Z ≥ 1 as follows ¯ β ( r, s ) : = r ζ s η  C ρ P − 1  c β 18 p + λ α 1 min  P − 1 0   c 1 3 α 1 + c 2 3 α 2   + c β 2 p  , (22) φ w ( r ) : = c β 18 p ¯ γ p a w ( k r k ) | P − 1 0 | + γ 1  3  k r k + γ − 1 w (3 ¯ γ w ( k r k ))  + c 2 3 α 2 ¯ γ α 2 w ( k r k ) , (23) φ v ( r ) : = c β 18 p ¯ γ p a v ( k r k ) | P − 1 0 | + γ 2  3  k r k + γ − 1 v (3 ¯ γ v ( k r k ))  + c 1 3 α 1 ¯ γ α 1 v ( k r k ) . (24) equation (21) can be written ∀ k ∈ Z [1 ,N − 1] as follows | x k − ˆ x k | k | ≤ ¯ β ( | x k − N − ¯ x k − N | , N ) + φ w ( k w k ) + φ v ( k v k ) . (25) T o guarantee the validity of previous results on the entire time horizon we must extend the definition of β ( r, s ) to s = 0 . Because of ¯ β ( r, s ) ∈ K L , ¯ β ( r, 0) ∈ K L and ¯ β ( r, 0) ≥ ¯ β ( r, k ) for k ∈ Z ≥ 1 , it is sufficient to define ¯ β ( r, 0) ≥ k β ¯ β ( r, 1) for some k β ∈ R > 1 to extend the definition of ¯ β ( r, s ) for all k ∈ Z ≥ 0 . W e would like to determinate the decreasing rate for the function ¯ β ( r, s ) N samplings time in the future. In order to do that, let define the constants µ ∈ R > 0 , δ > 2 + µ 1 + µ and r max : = max { 1 δ  ¯ β ( e max , 0) + φ w ( k w k ) + φ v ( k v k )  , δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) } The minimum horizon length required to accomplish a decreasing rate δ will be giv en by N ≥  δ ζ r ζ − 1 max C ρ P − 1  c β 18 p + λ α 1 min  P − 1 0  ( c 1 3 α 1 + c 2 3 α 2 ) + c β 2 p  1 η (26) Adopting an estimator with a window length greater or equal to N such that ¯ β ( δ r, N ) ≤  N N  η r , (27) the ef fects of the initial conditions will vanish with a decreasing rate δ . As k → ∞ , the estimation will entry to the bounded set X ( w , v ) ∈ X defined by the noises of the system X ( w , v ) : = {| x k + j − ˆ x k + j | k + j | ≤ δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) } . (28) This set define the minimum size region of error space X that the error can achie ve by removing the effect of errors in initial condi- tions ( e max ). Equation (27) establish a trade off between speed of con ver gence and window length, which is related with the size of X ( w , v ) . For any MHE with adaptiv e arriv al cost and window length N ≥ N two situations can be considered • The estimator removed the effects of x 0 on ˆ x k + j | k + j such that x k + j − ˆ x k + j | k + j ∈ X ( w , v ) , and • The estimator has not removed the effects of x 0 on ˆ x k + j | k + j such that x k + j − ˆ x k + j | k + j / ∈ X ( w , v ) , Assuming the first situation and recalling equations (25) and (27), the following inequalities hold | x k + N − ˆ x k + N | k + N | ≤ ¯ β ( | x k − ¯ x k | , k ) + φ w ( k w k ) + φ v ( k v k ) , ≤ | x k − ¯ x k | δ  N N  η + φ w ( k w k ) + φ v ( k v k ) , ≤ (2 + µ ) ( φ w ( k w k ) + φ v ( k v k )) , ≤ δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) . (29) This equation implies the fact that the estimation error x k + j − ˆ x k + j | k + j ∈ X ( w , v ) ∀ j ∈ Z ≥ 0 . IET Research Journals, pp. 1–10 c  The Institution of Engineering and T echnology 2015 5 In the other case, when the estimation error is outside of X ( w , v ) , equations (25) and (27) are recalled again and the follo wing inequal- ities hold | x k + N − ˆ x k + N | k + N | ≤ | x k − ¯ x k | δ  N N  η + φ w ( k w k ) + φ v ( k v k ) , ≤ | x k − ¯ x k | δ  N N  η + | x k − ¯ x k | δ (1 + µ )  N N  η , ≤| x k − ¯ x k |  N N  η  2 + µ δ (1 + µ )  . (30) Since δ > 2+ µ 1+ µ , then ∀ N ≥ N we hav e θ : =  N N  η  2 + µ δ (1 + µ )  < 1 . (31) Equations (30) and (31) rev eal a contractiv e behaviour of the esti- mation error with θ as contraction factor . For some finite time k ∗ the estimation error will decrease until x k ∗ + j − ˆ x k ∗ + j | k ∗ + j ∈ X ( w , v ) . In an equiv alent formulation, equations (29) and (30) put in evi- dence the existence of a positi ve inv ariant set and a L yapunov like function for the proposed estimator . From equation (30), one can see that for the case that the estimation error belong to the set X ( w, v ) C ∩ X , the estimation error decreases in a factor of θ every N sampling time. T aking in account the general case in which | x k − ˆ x k | k | ∈ X for k ∈ Z ≥N , follo wing the same proce- dure as in [23], we could define i : = b k N c (where b·c denotes the floor function) and j : = k mod N , therefore k = iN + j . Combin- ing equations (29) and (30) and the fact that | x j − ˆ x j | ≤ δ r max for j ∈ Z [0 ,N − 1] one can obtain | x k − ˆ x k | k | ≤ max {| x j − ¯ x j | θ i , δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) } , ≤ max { θ i  ¯ β ( | x 0 − ¯ x 0 | , j ) + φ w ( k w k ) + φ v ( k v k )) , δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) } , ≤ θ i ¯ β ( | x 0 − ¯ x 0 | , j ) + δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) , ≤ Φ ( | x 0 − ¯ x 0 | , k ) + δ (1 + µ ) ( φ w ( k w k ) + φ v ( k v k )) . (32) where Φ ( | x 0 − ¯ x 0 | , k ) : = θ i ¯ β ( | x 0 − ¯ x 0 | , j ) j ∈ [0 , N − 1] . Since ¯ β ( r, s ) ∈ K L , function Φ ( | x 0 − ¯ x 0 | , k ) could increase in the steps Z [ iN − 1 ,iN ] for i ≥ 1 (recall definition in Equation (22)). Therefore, define ¯ Φ ( | x 0 − ¯ x 0 | , k ) which is an upper bound for Φ ( | x 0 − ¯ x 0 | , k ) . T aking in account that noises at time ≥ k do not affect the estimation at time k , equation (32) can be rewritten as | x k − ˆ x k | k | ≤ ¯ Φ ( | x 0 − ¯ x 0 | , k ) + δ (1 + µ )  φ w ( k w k ) [0 ,k − 1] + φ v ( k v k ) [0 ,k − 1]  . (33) This equation is just equation (14) with Φ ( | x 0 − ¯ x 0 | , k ) = ¯ Φ ( | x 0 − ¯ x 0 | , k ) , (34) π w  k w k [0 ,k − 1]  = δ (1 + µ ) φ w  k w k [0 ,k − 1]  , (35) π v  k v k [0 ,k − 1]  = δ (1 + µ ) φ v  k v k [0 ,k − 1]  , (36) therefore the estimator proposed in equations (3) is RGAS. Finally , in order to prov e that the estimation error | x k − ˆ x k | k | → 0 when | w k | ∈ C w , | v k | ∈ C v , we must note that equation (25) holds for k w k [ k − N ,k − 1] instead of k w k and k v k [ k − N ,k − 1] instead k v k (it can be done omitting last step in equation (15)). From a qualitative point of view , taking in account that function ¯ Φ ( | x 0 − ¯ x 0 | , k ) ∈ KL and sequences w and v are con vergent, the right hand side of equation (33) tends to zero as k → ∞ .  The proof of Theorem 1 is constructive and provides an esti- mate of the estimation horizon N required to guarantee RGAS of the MHE proposed in this work. The estimates N and functions Φ , π w and π v can be quite conserv ative, since their deriv ation inv olved conservati ve estimates of noises, errors, stage costs and arriv al cost. Note that the minimum horizon necessary to guarantee RGAS N depends on r max , which depends on the class of disturbances con- sidered (upper bounds of noises and error), the initial value of the prior weighting matrix P 0 and the bounds of the stage cost. The minimum horizon length is independent of k w k , k v k , and the same N ensures RGAS for all bounded disturbances and bounded prior error , like the result obtained by [23]). This implies that we can prove the RGAS property for full information estimator with least–square objectiv e function. Remark 1. Functions φ w and φ v in equations (23) and (24) , and hence π w and π v in equations (35) and (36) , do not depend on the estimation horizon N which means that the moving horizon estima- tor with adaptive arrival cost is RGAS with uniform gains given by (23) and (24) . 4 Examples The following examples will be used to illustrate the results pre- sented in the previous sections and compare the performance of the estimators. The examples considered in this work are taken from [23] for a direct comparison of the results. 4.1 Example 1 The first example considers the system x ( t + 1) =  0 . 8 x 0 ( t ) + 0 . 2 x 1 ( t ) + 0 . 5 w ( t ) − 0 . 3 x 0 ( t ) + 0 . 5 cos( x 1 ( t ))  (37) y ( t ) = x 1 ( t ) + v ( t ) The stage cost is chosen as ` ( w, v ) = 10 w 2 + 10 v 2 and the hori- zon length is N = 10 . The prior weighting is chosen as Γ( χ ) = 0 . 1( χ − ˆ x ( t | t )) T ( χ − ˆ x ( t | t )) for the MAX estimator ([23]) and Γ t ( χ ) = ( χ − ˆ x ( t | t )) T Π − 1 k ( χ − ˆ x ( t | t )) for the AD AP estimator (our method), where Π 0 = 10 I 2 and Π k is obtained using equa- tions (5) with σ = 0 . 2 and c = 1 e 6 . The MAX estimator uses δ = 1 , δ 1 = κ N with κ = 0 . 89 2 and δ 2 = 1 / N (see equation (3) of [23]). The full information estimator (FIE MAX, see [22]) is configured with the same parameter used by [23], maintaining the stage cost and prior weighting Γ 0 , and δ = 1 , δ 1 = κ t and δ 2 = 1 /t . T able 1 Example 1 averaged MSE over 300 tr ials. FIE MAX AD AP MAX EKF x 0 0.02040 0.02176 0.02206 0.02296 x 1 0.00135 0.00151 0.00156 0.00154 IET Research Journals, pp. 1–10 6 c  The Institution of Engineering and T echnology 2015 0 10 20 30 40 50 60 70 80 k 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 x 0 R E A L F I E M A X A D A P M A X E K F Fig. 1 : Comparison between ADAP (red dash dotted), MAX (blue dashed), FIEMAX (green dotted), EKF (magenta) estimators, and real system state (black solid). 0 10 20 30 40 50 60 70 80 k 0.0 0.1 0.2 0.3 0.4 0.5 x 1 R E A L F I E M A X A D A P M A X E K F Fig. 2 : Comparison between ADAP (red dash dotted), MAX (blue dashed), FIEMAX (green dotted), EKF (magenta) estimators, and real system state (black solid). T able 1 sho ws the mean square estimation error of each estimator av eraged ov er 300 trials. It can be seen that the proposed estima- tor average mean square estimation error is smaller than MAX ones and closer to FIE MAX. The main performance difference between AD AP and FIEMAX estimators is the inclusion of the max term in the last one, which allo ws to follow the sudden changes (see Figures 1 and 2). Figures 1 and 2 sho ws simulation results with initial condition x 0 = [0 . 5 , 0] T and prior estimate ¯ x 0 = [0 , 0] T . The process and measurement disturbances w and v are sampled from an uniform distribution o ver the intervals [ − 0 . 3 , 0 . 3] and [ − 0 . 2 , 0 . 2] , respec- tiv ely . This figure shows that the estimators that use the max term are able of following the sudden changes, howe ver in the remaining of the signal the MAX estimator is mo ving away of the FIEMAX while AD AP remains closer . 4.1.1 MHE in the presence of variable measurement noise: Now the MHE estimator is e valuated in the presence of time- varying measurement noise. The variance of the measurement noise is changed from 0 . 2 to 1 . 0 between times 20 and 40, then it returns to 0 . 2 . T able 2 Example 1 aver . MSE over 300 tr ials with variable measurement noise. AD AP MAX FIE MAX x 0 0.02068 0.03067 0.00761 x 1 0.00290 0.00335 0.00068 T able 2 shows the average mean square error in the presence of variable measurement noise. In this case we can see that the behaviour of the proposed estimator is marginally affected by the variations of the measurement, while the mean square error of x 0 of other estimators increase significantly . These behaviours are due to the adaptation capabilities of the prior weighting updating mech- anism, which is able of tracking the changes of noises, in the case of AD AP estimator, and the effect of the max term in MAX and FIEMAX estimators. Figure 3 shows the ev olution of the trace of P − 1 k − N used in the prior weight of ADAP estimator in both examples. It can be seen that the trace of both matrices grow in similar way , howe ver when the measurement noise changes its v ariance from 0 . 2 to 1 . 0 the trace of P − 1 k − N increases its value (from 12 . 5 to 22 . 5 ) and them both traces hav e the same behaviour . 4.2 Example 2 As a second example, we consider a second order gas-phase irre- versible reaction of the form 2 A → B . This example has been considered in the context of moving horizon estimation in [29], [22] and [23]. Assuming an isothermal reaction and that the ideal gas law holds, the system dynamics ˙ x =  − 2 k x 2 0 k x 2 0  h ( x ) = x 0 + x 1 (38) where x = [ x 0 , x 1 ] , x 0 is the partial pressure of the reactant A , x 1 is the partial pressure of the product B , and k = 0 . 16 is the reac- tion rate constant. The measured output of the system is the total pressure. The system is affected with additiv e process and measure- ment noise w and v drawn from normal distrib utions with zero mean and cov ariance Q w = 0 . 001 2 I 2 and R v = 0 . 1 2 , respectiv ely . The stage cost and prior weighting are chosen as ` ( w, v ) = w T Q − 1 w w + R − 1 v v 2 and Γ t ( χ ) = ( χ − ˆ x ( t | t )) T Π − 1 k ( χ − ˆ x ( t | t )) with Π 0 = (1 / 36) I 2 , where Π k is determined by an extended Kalman filtering recursion in the case of the MAX estimator and the adapti ve method in the case of the AD AP estimator with σ = 0 . 1 and c = 1 e 6 . For the MAX estimator we use δ 1 = 1 / N , δ 2 = 1 and δ = 0 . In the case of the AD AP estimator, the stage cost weight matrices are chosen as Q w = 0 . 001 I 2 and R v = 0 . 1 . W e use a multiple shooting strategy with a sampling time of ∆ = 0 . 1 and we add the restrictions x 0 ≥ 0 and x 1 ≥ 0 . T able 3 shows the values of the mean squared error computed from the time 10 (in order to neglect the initial transient error) up to the simulation end time and averaged over 300 trials for horizon sizes of N = 5 and N = 10 . 0 10 20 30 40 50 60 70 80 k 0 5 10 15 20 25 30 35 T r a c e ( P − 1 k ) Fig. 3 : Comparison of the ev olution of trace ( P − 1 k − N ) used by AD AP estimator for time–varying (red dash dotted) and constant (blue dashed) measurement noise parameters. IET Research Journals, pp. 1–10 c  The Institution of Engineering and T echnology 2015 7 T able 3 Example 2 averaged MSE over 300 tr ials and diff erent horizon size. N=2 N=5 N=10 AD AP MAX ADAP MAX AD AP MAX FIE x 0 0.18808 0.58652 0.03367 0.04615 0.00171 0.00772 0.00024 x 1 0.23037 0.66768 0.04074 0.05077 0.00285 0.00951 0.00120 0 10 20 30 40 50 60 70 80 k 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 x 0 A D A P M A X F I E R E A L 0 10 20 30 40 50 60 70 80 k 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 x 0 A D A P M A X F I E R E A L 0 10 20 30 40 50 60 70 80 k 0.0 0.5 1.0 1.5 2.0 2.5 3.0 x 0 A D A P M A X F I E R E A L Fig. 4 : Comparison between ADAP (red dash dotted), MAX (blue dashed), FIEMAX (green dotted) estimators and real system state (black solid) for different horizon length ( N = 2 , 5 and 10 ). Figures 4 and 5 show simulation results with x 0 = [3 , 1] T and ¯ x 0 = [0 . 1 , 4 . 5] T and horizons of sizes N = 2 , 5 and 10 , along with results for a full information estimator using the same parameters: the same stage cost ` ( · ) , prior weighting Γ 0 , δ = 0 , δ 1 = 1 /t and δ 2 = 1 . These figures show that the behaviour of AD AP estimator hardly change with horizon length (only the startup beha viours show differences) and no offset in the estimates, while the beha viour of the MAX estimator changes significantly . In addition to the cycling effect caused by the use of the filtered estimate to update ¯ x k − N [24], the MAX estimator also exhibits offset in the estimate that depends on the estimation horizon length. 0 10 20 30 40 50 60 70 80 k 0.0 0.5 1.0 1.5 2.0 2.5 x 1 A D A P M A X F I E R E A L 0 10 20 30 40 50 60 70 80 k 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 x 1 A D A P M A X F I E R E A L 0 10 20 30 40 50 60 70 80 k 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 x 1 A D A P M A X F I E R E A L Fig. 5 : Comparison between ADAP (red dash dotted), MAX (blue dashed), FIEMAX (green dotted) estimators and real system state (black solid) for different horizon length ( N = 2 , 5 and 10 ). 5 Conclusions In this paper we established robust global asymptotic stability for moving horizon estimator with a least-square type cost function for nonlinear detectable (i-IOSS) systems in presence of bounded distur- bances. It was also sho wn that the estimation error con verges to zero in case that disturbances con verge to zero. This was done for an esti- mator which uses a least-square type cost function whose arrival cost us updated using adaptiv e estimation methods. An advantage of this updating mechanism is that the required conditions on prior weight- ing are such that it can be chosen off-line. Furthermore, it introduces a feedback mechanism between the arriv al cost weight and the esti- mation errors that automatically controls the amount of information used to compute it, which allows to shorten the estimation horizon. The standard least-square type cost function is typically used in practical applications and RGAS has been proved in [23]. Howe ver , for this formulation, the disturbances gains depend on the estimation horizon. Hence, this result does not allow to establish robust global asymptotic stability for a full information estimator . W e sho wed that IET Research Journals, pp. 1–10 8 c  The Institution of Engineering and T echnology 2015 changing the updating mechanism of arriv al cost weight the dis- turbances gains becomes uniform, allowing to extend the stability analysis to full information estimators with least-square type cost functions. Ackno wledgment The authors wish to thank the Consejo Nacional de In vestiga- ciones Cientificas y T ecnicas (CONICET) from Argentina, for their support. IET Research Journals, pp. 1–10 c  The Institution of Engineering and T echnology 2015 9 6 References 1 Andrew H Jazwinski. Stochastic processes and filtering theory . Courier Corpora- tion, 2007. 2 John L Crassidis and John L Junkins. Optimal estimation of dynamic systems . Chapman and Hall/CRC, 2004. 3 Huaizhong Li and Minyue Fu. A linear matrix inequality approach to robust h/sub/spl infin//filtering. IEEE Tr ansactions on Signal Pr ocessing , 45(9):2338– 2350, 1997. 4 Ali H Sayed. A framework for state-space estimation with uncertain models. IEEE T ransactions on A utomatic Control , 46(7):998–1013, 2001. 5 Franco Blanchini and Stefano Miani. Set-theoretic methods in control. systems & control: Foundations & applications. Birkhäuser . Boston, MA , 2008. 6 Laurent El Ghaoui and Hervé Lebret. Rob ust solutions to least-squares prob- lems with uncertain data. SIAM Journal on matrix analysis and applications , 18(4):1035–1064, 1997. 7 K Hu and J Y uan. Improved robust h inf. filtering for uncertain discrete-time switched systems. IET Contr ol Theory & Applications , 3(3):315–324, 2009. 8 Xing Zhu, Y eng Chai Soh, and Lihua Xie. Design and analysis of discrete-time robust kalman filters. Automatica , 38(6):1069–1077, 2002. 9 A Jazwinski. Limited memory optimal filtering. IEEE Transactions on A utomatic Contr ol , 13(5):558–563, 1968. 10 Fred C Schweppe. Uncertain dynamic systems . Prentice Hall, 1973. 11 Christopher V Rao, James B Rawlings, and Jay H Lee. Constrained linear state estimationâ ˘ A ˇ T a moving horizon approach. Automatica , 37(10):1619–1628, 2001. 12 Christopher V Rao, James B Rawlings, and David Q Mayne. Constrained state estimation for nonlinear discrete-time systems: Stability and moving horizon approximations. IEEE transactions on automatic contr ol , 48(2):246–258, 2003. 13 Angelo Alessandri, Marco Baglietto, and Giorgio Battistelli. Robust receding- horizon state estimation for uncertain discrete-time linear systems. Systems & Contr ol Letters , 54(7):627–643, 2005. 14 Angelo Alessandri, Marco Baglietto, and Giorgio Battistelli. Mo ving-horizon state estimation for nonlinear discrete-time systems: New stability results and approximation schemes. Automatica , 44(7):1753–1765, 2008. 15 Angelo Alessandri, Marco Baglietto, and Giorgio Battistelli. Min-max moving- horizon estimation for uncertain discrete-time linear systems. SIAM Journal on Contr ol and Optimization , 50(3):1439–1465, 2012. 16 J Garcia-Tirado, H Botero, and F Angulo. A new approach to state estimation for uncertain linear systems in a moving horizon estimation setting. International Journal of A utomation and Computing , 13(6):653–664, 2016. 17 Hossein Sartipizadeh and T yrone L Vincent. Computationally tractable robust moving horizon estimation using an approximate conve x hull. In Decision and Contr ol (CDC), 2016 IEEE 55th Conference on , pages 3757–3762. IEEE, 2016. 18 G Sánchez, M Murillo, and L Giovanini. Adaptiv e arrival cost update for improving moving horizon estimation performance. ISA transactions , 68:54–62, 2017. 19 James B Rawlings and David Q Mayne. Model predictive control: Theory and design. 2009. 20 James B Ra wlings and Luo Ji. Optimization-based state estimation: Current status and some new results. Journal of Pr ocess Contr ol , 22(8):1439–1444, 2012. 21 Angelo Alessandri, Marco Baglietto, Giorgio Battistelli, and Victor Zavala. Advances in moving horizon estimation for nonlinear systems. In Decision and Contr ol (CDC), 2010 49th IEEE Conference on , pages 5681–5688. IEEE, 2010. 22 Luo Ji, James B Rawlings, Wuhua Hu, Andrew W ynn, and Moritz Diehl. Robust stability of moving horizon estimation under bounded disturbances. IEEE T ransactions on A utomatic Control , 61(11):3509–3514, 2016. 23 Matthias A Müller. Nonlinear moving horizon estimation in the presence of bounded disturbances. Automatica , 79:306–314, 2017. 24 Peter Klaus Findeisen. Moving horizon state estimation of discrete time systems . PhD thesis, Univ ersity of Wisconsin–Madison, 1997. 25 Eduardo D Sontag and Y uan W ang. Output-to-state stability and detectability of nonlinear systems. Systems & Contr ol Letters , 29(5):279–290, 1997. 26 Eduardo D Sontag. Some connections between stabilization and factorization. In Decision and Control, 1989., Proceedings of the 28th IEEE Conference on , pages 990–995. IEEE, 1989. 27 Z-P Jiang, Andrew R T eel, and Laurent Praly . Small-gain theorem for iss sys- tems and applications. Mathematics of Contr ol, Signals and Systems , 7(2):95–120, 1994. 28 Eduardo D Sontag. Input to state stability: Basic concepts and results. In Nonlinear and optimal contr ol theory , pages 163–220. Springer, 2008. 29 Eric L Haseltine and James B Rawlings. Critical evaluation of extended kalman filtering and moving-horizon estimation. Industrial & engineering chemistry r esearc h , 44(8):2451–2460, 2005. IET Research Journals, pp. 1–10 10 c  The Institution of Engineering and T echnology 2015