Enhancement of the Moving Horizon Estimation Performance Based on an Adaptive Estimation Algorithm

Although moving horizon estimation (MHE) is a very efficient technique for estimating parameters and states of constrained dynamical systems, however, the approximation of the arrival cost remains a major challenge and therefore a popular research topic. +e importance of the arrival cost is such that it allows information from past measurements to be introduced into current estimates. In this paper, using an adaptive estimation algorithm, we approximate and update the parameters of the arrival cost of the moving horizon estimator. +e proposed method is based on the least-squares algorithm but includes a variable forgetting factor which is based on the constant information principle and a dead zone which ensures robustness. We show by this method that a fairly good approximation of the arrival cost guarantees the convergence and stability of estimates. Some simulations are made to show and demonstrate the effectiveness of the proposed method and to compare it with the classical MHE.


Introduction
Since the 1960s, the state estimation of systems by a moving horizon approach has been highlighted. e dazzling success of this estimation technique is based on its ability to explicitly take into account the constraints on the estimate of variables of a dynamic system [1][2][3]. In [4], the problem of estimating the state of a linear system on a peer-to-peer network of linear sensors is discussed. e proposed approach is fully distributed and scalable and allows for taking into account the constraints on noise and state variables by resorting to the moving horizon estimation paradigm. Each node in the network calculates its local state estimate by minimizing a cost function defined over a sliding window of fixed size. e cost function includes a fused arrival cost which is computed in a distributed way by performing a consensus on the local arrival costs. In [5], the authors have proposed a novel type of iterative partition-based moving horizon estimators, the estimates of which approach those of a centralized moving horizon estimator as the number of iterations increases. It uses a deterministic setting and can handle known inputs as well as bounds on the estimated state. ey have derived conditions on the system, its partitions, and the scalar regularization parameter, which guarantee convergence towards the optimal centralized state estimate as well as stability of the estimation error dynamics, even with a finite number of iterations at each time step. e authors, in [6], have investigated on a problem of eventtriggered state estimation for networked linear systems. ey consider that the stochastic system disturbances and noise are bounded and moving horizon estimation (MHE) is used to handle these constraints and establish an event-based state estimation mechanism that aims to provide good state estimates while reducing the frequencies of both the evaluation of the state estimator and networked communication between the plant and the estimator. In [7], the authors focus on distributed moving horizon estimation (DMHE) for a class of two-time-scale nonlinear systems described in the framework of singularly perturbed systems. By taking advantage of the time-scale separation property, a two-timescale system is first decomposed into a reduced-order fast system and a reduced-order slow system. e slow system is further decomposed into several interconnected slow subsystems. In the proposed distributed state estimation scheme, a local estimator is designed for each slow subsystem and for the reduced-order fast system.
In the moving horizon estimation approach, the state estimate is determined by solving an optimization problem online that minimizes the sum of the squared errors. At a sampling time, when a new measurement is obtainable, the older measuring within the estimation window is discarded and the horizon finite optimization problem is solved again to obtain the new state estimate [8][9][10]. To take into account the historical data outside the estimation window, a variable that summarizes the information from this data is included in the objective function of the MHE optimization problem called arrival cost. is arrival cost is crucial for the stability and performance of the moving horizon estimation; therefore, its approximation is an active research topic.
To reduce the size of the estimation window and the size of the optimization problem, a good approximation of the arrival cost is necessary, ensuring good performance and good robustness. e arrival cost is commonly estimated using a weighted norm of states at the beginning of the estimated horizon window [8,11,12]. In [13], the authors proposed to use the Kalman filter approach for updating the term of the arrival cost (in particular the weighting matrix) for the case of linear systems, leading to an inadequate approximation of the weighting matrix leading to relatively poor estimates due to a Gaussian distribution of the Kalman filter. To overcome this problem, the authors in [14] proposed an iterative scheme for updating the arrival cost using the information on the active or inactive constraints from the previous iteration and a quadratic approximation. e idea brings the quadratic approximation as close as possible to the optimal solution of the arrival cost. e fundamental assumption behind the updating scheme is that the constraints should remain unchanged after a certain moment despite the increase of the estimation horizon. For a large estimation window, this assumption is functional giving good results; however, if after smoothing some constraints become inactive, less good results are obtained.
In [15], the author uses nonlinear filters based on sampling, including the unscented Kalman filter based on deterministic sampling, the particle filter based on random sampling, and the cell filter based on the aggregate Markov chain to update the arrival cost of the moving horizon estimation and demonstrate the advantages and performance of those filters over the traditional extended Kalman filter approach. However, the computation time is relatively long and better performance is also obtained due to a fairly large horizon length.
In more recent works, as in [16], the authors proposed a moving horizon estimation algorithm based on multiple estimation windows. e major advantage of this algorithm is that it reduces the size of the optimization problem while maintaining the performance properties of the estimator and its stability by taking advantage of the inactivity of the constraints. In [17], the authors developed a two-step design strategy, namely, decoupling step and convergence step for networked linear systems with unknown inputs under dynamic quantization effects. In the decoupling step, the decoupling parameter of the moving horizon estimator is designed based on certain assumptions on system parameters and quantization parameters. In the convergence step, by employing a special observability decomposition scheme, the convergence parameters of the moving horizon estimator are achieved such that the estimation error dynamics is ultimately bounded. e conventional moving horizon estimation strategy turns out to be unable to guarantee satisfactory performance as the estimation error depends on external disturbances. In [18], a moving horizon estimation strategy is formulated and the corresponding analytical solution of the state estimation is derived by using the completing the square technique. Based on the lifting approach, the dynamics of the estimation error covariance is analyzed and the desired estimator parameters are obtained by solving a set of linear matrix inequalities.
In this work, the weighting matrix is updated using an adaptive estimation algorithm, which is the least-squares algorithm including a dead zone and a variable forgetting factor [19,20]. e dead zone here is used to obtain the robustness; the formula of a forgetting factor based on the principle of constant information can be applied to the process model with unmodeled process uncertainties or external disturbances. Uncertainties can cause the estimates to differ when the standard recursive least-squares algorithm is used, but to prevent this, the inclusion of the dead zone in the estimation algorithm is required. erefore, the weighting matrix will be calculated in closed loop unlike standard estimation techniques. e effectiveness of the proposed approach is evaluated by carrying out simulation studies on benchmark problems found in the literature. e main contributions of this work aim to (i) show that an adaptive estimation method can be implemented to approximate the arrival cost and have much better efficiency than the traditional MHE, (ii) also show that the proposed scheme guarantees robustness to the estimator for sudden disturbances, and (iii) finally show that the bounds of estimates are crucial when updating the arrival cost. e paper is organized as follows. In Section 2, we present the moving horizon estimation problem, followed by the arrival cost update mechanism in Section 3. Section 4 presents a stability analysis; the results are discussed in Section 5, and finally, Section 6 will be devoted to the conclusions.

Preliminaries and Problem Statement
Notation 1. In the following, we will consider the following notations, ξ ∈ R n is a column vector and ξ T ∈ R n will be its transpose. ξ j|ℓ denotes a finite sequence of vectors over a given index up to ℓ. Let B ∈ R n×n be a n × m matrix and B T is its transpose. If B ∈ R n×n , then B − 1 denotes its inverse. Let P ∈ R n×n be a symmetric real matrix, P is said to be positive definite if Ξ T P Ξ > 0 and ∀Ξ ∈ R n . e Euclidean matrix or 2 Journal of Control Science and Engineering vector norm is denoted by ‖ · ‖, for instance, ξ ∈ R n and P ∈ R n×n are positive definite; it follows that ‖ξ‖ p � ���� ξ T Pξ.
Consider the linear time invariant discrete system: where x ℓ ∈ X θ ⊆ R n is the state vector, w ℓ ∈ W θ ⊆ R p is the process noise vector, y ℓ ∈ Y θ ⊆ R m is the measurement vector, and v ℓ ∈ V θ ⊆ R m is the measurement noise vector. e process noise w ℓ and the measurement noise v ℓ are assumed to be bounded and unknown eventually, i.e., Using all measures for the full information estimator to estimate the system states x ℓ increases the size of the optimization problem over time and leads to heavy calculation.
To overcome this problem, the moving horizon estimation was developed to take into account a certain amount of data and possibly to make dynamic updates of certain parameters in an estimation window moving in time.
A penalty term called arrival cost is essential in the formulation of the MHE because it takes into account past measurements to estimate current states. e problem of the MHE aims to determine, at a time ℓ, an estimate x ℓ|ℓ of the current state and process noises which solve the optimization problem.
For ℓ ≤ N, the full information problem arises: x j+1|ℓ � Ax j|ℓ + w j|ℓ , However, for all ℓ > N, the optimization problem becomes x j+1|ℓ � Ax j|ℓ + w j|ℓ , where the optimal estimated is x ℓ− j|ℓ , w j|ℓ is the process noise estimate based on available measurements at time ℓ (y ℓ− j ), and v j|ℓ � y j − Cx j|ℓ are the estimation residuals. To ensure a fairly good robust stability of the MHE, a judicious and adequate choice of Γ ℓ− N|ℓ (·) and its parameters is necessary [21]. Usually, a relatively long estimation window is used to overcome a poor estimate caused by an inadequate cost function approximation [22]; this case is strongly encountered in the conventional MHE. It can also be noted that, for a short estimation window, an approximation of the optimal cost function is required to properly incorporate the past information into the problem, in case of which a poor estimate will be revealed.
In many dynamical systems, in particular linear systems, a Gaussian distribution is assumed and the arrival cost function is updated using a Kalman filter [23]. Besides, the MHE is equivalent to the Kalman filter for unconstrained linear systems with indeed the same Gaussian distribution. In the following, we will present the proposed mechanism for updating the arrival cost parameters of the MHE.

Arrival Cost Update Algorithm
For the constrained MHE problem, it is difficult to find an exact analytical expression of the arrival cost; therefore, an approximation of the arrival cost used for the unconstrained problem will be used; this quadratic approximation is of the following form: where P ℓ− N is the weighting matrix and x ℓ− N is the initial state. Both define the approximation of the arrival cost function. e proposed approach for updating the arrival cost in this work, based on an adaptive algorithm, is formulated as an update of the following elements: x ℓ− N and P ℓ− N . e update of the initial state x ℓ− N will be made by a smooth update, which implies that once the estimate leaves the estimation horizon, the state of the constraints does not change [14,24]: Regarding the update of the weighting matrix P ℓ− N , it is necessary to notice from the outset that, in a stochastic interpretation of the MHE problem, the weighting matrix can be seen as the covariance of the state x ℓ− N [8]. For unconstrained linear systems, analytical solutions exist to compute the weighting matrix; on the contrary, for constrained linear systems, a deduction of the approximate solutions is made, which are fundamentally and based on a recursive process of updating useful information from the system at a given time. In the proposed adaptive estimation method, an update based on a recursive process relying on states (or not) and measurements is used [19,25]: With P 0 > 0, where λ ℓ is the forgetting factor with 0 < λ ℓ ≤ 1, α ℓ is a function which varies between 0 and 1 called the dead zone. Equation (6) produces a recursively Journal of Control Science and Engineering real-time estimate of P ℓ− N , by updating P ℓ− N− 1 (previous estimate) with a dynamic temporal average over time of We can qualify this updating mechanism as a dynamic temporal filter overtime with P 0 , the initial condition and Proof. See Lemma 7 in [26]. Note that the forgetting factor λ ℓ and the function α ℓ each have a different role. When λ ℓ ≈ 1, the algorithm tends to keep the old data in P ℓ− N− 1 ; in this case, all information is included, improving the estimation problem. On the contrary, when λ ℓ ≪ 1, the algorithm tends to discard the old information in P ℓ− N− 1 , considering recent data to estimate P ℓ− N [27].
is case allows avoiding many old data that could negatively affect the estimator performance. e dead zone function α ℓ is used to ensure robustness when estimating.
By the matrix inversion Lemma, we can rewrite (6) as follows: A feedback loop is introduced into the MHE problem by this way to compute P ℓ− N ; this way makes it possible to ensure the convergence of the estimates x ℓ− N|ℓ and w j (j � ℓ − N, . . . , N) and to adjust P − 1 ℓ− N according to the amount of available information. e adaptive algorithm used proposes a variable forgetting factor varying between 0 and 1 and under certain conditions to be respected and a dead zone, both allowing P ℓ− N to be updated as follows [19]: where N 0 is a positive design parameter chosen which is quite large: and d 1 and d 2 are also design parameters, which are chosen positive.
Any sequence λ ℓ is bounded by some σ * > 0. us, the size of dead zone Δ 2 ℓ (1 + μ ℓ /λ ℓ ) is also bounded because μ ℓ is bounded too. e forgetting factor proposed in [19] is based on the constant information principle, which is adequate for external disturbances and unmodeled system uncertainties.

Stability Analysis
e stability of the estimator requires that the prediction error converges to zero for the following nominal system: with no measurement and state noise. In [8], the authors require that the system has to evolve according to the constraints because, for the constrained systems' case, a poor range of choice of constraints can prevent convergence to the real values of the system.

Lemma 2.
e estimation algorithm assumes that sup(‖δ ℓ ‖) ≤ d 1 and sup(‖w ℓ ‖) ≤ d 2 , with δ ℓ as the system uncertainty vector; then, the following properties are valid: (i) P ℓ converges and is uniformly bounded Proof. By recursive calculation and the matrix inversion lemma, we get from (7) It follows that which implies the uniform boundedness of P ℓ ; by the property 0 < σ * ≤ σ 0,ℓ , it is easily proven from (11) that P ℓ uniformly converges.
Proof. e update expression of P ℓ− N given by equation (7) can be seen as a particular case of Ricatti equation: It is stable if and only if the eigenvalues of the transfer matrix are strictly inside the unit circle of the complex plane: Hence, All eigenvalues of matrix Λ are greater than 0 and bounded by 1; then, the eigenvalues of (4) are also greater than 0 and bounded by 1. Hence, the matrix P ℓ , therefore, has a steady-state solution P ∞ .

be a Lyapunov function for the moving horizon estimator with the adaptive arrival cost, where x *
ℓ � x ℓ − x ℓ and the matrix P − 1 ℓ is given by equation (6) and A T x ℓ− N|ℓ− 1 x T ℓ− N|ℓ− 1 A are symmetric positive definite and symmetric positive semidefinite, respectively, for all ℓ; then, if (10) is observable, it follows that the MHE with the proposed adaptive arrival cost update is an asymptotically stable observer.

Simulation Results
To demonstrate the performance of the proposed method, we adopted discrete-time models used by [13,29].

Example 1.
where w ℓ � |s ℓ | with s ℓ a sequence of independent, normally distributed random values with zero mean and covariance equal to 1, and v ℓ is also a sequence of independent, normally distributed random values with zero mean and covariance equal to 0.1. e initial state x 0 with zero mean is also normally distributed with a unit covariance. e constrained estimation problem here is formulated with Q � I, R � 0.01, P 0 � 0.5I, and x 0 � [0.5, 0.5] T . We choose N � 10 and d 1 � 0.01 and d 2 � 0.005 and N 0 � 50 as design parameters. e proposed algorithm is compared with the classical MHE (MHE KF ) which propagates the initial state and the arrival cost estimate using Kalman filter. e sum square estimation error is used as a benchmark: where x (i) denotes the i th component of the vector x. Based on 70 trials, the average sum square estimation error was computed. e results are shown in Table 1, where it can be seen that the performance of MHE AD is superior to the classical MHE. Figures 1(a) and 1(b) show the states and their estimates x (1) and x (2) with the horizon length N � 10.
We can see that despite the relatively short horizon length chosen, the MHE AD estimate is closer than MHE KF . e estimates of MHE AD tend to follow MHE KF for x (1) . On the contrary, for the state x (2) , we observe a divergent behavior. At every sample time, the square error is shown in Figures 1(c) and 1(d); we can observe that the error of the proposed method is lower than MHE KF , and we can also notice that the performance of MHE KF is quite degraded, probably, because of the relatively short horizon length.
In Figures 2(a) and 2(b), we can observe the evolution of the forgetting factor and the dead zone function, respectively. e forgetting factor is varying between 0 and 1, providing a good adaptation capability and allowing the proposed algorithm to insert pertinent information into the data and follow sudden variations. e dead zone function makes some variations between 0 and 1 to allow a rapid adaptation of the weighting matrix when a sudden change occurs. We can note that when the estimation error is large, the dead zone function is equal to 1, and when the error is relatively low, the dead zone function is equal to 0. Figures 3(a) and 3(b) allow appreciating the evolution of the trace of the weighting matrix P ℓ for both estimators.
e MHE KF weighting matrix converges to the steadystate value and does not change anymore, while MHE AD adapts its weighting matrix according to variations of the process noise.
A particular case is studied and presented in Figure 4, where the process noise covariance σ w varies with time. e covariance will take values 0.9 when 0 ≤ ℓ < 5, 0.55 when 5 ≤ ℓ < 20, 0.8 when 20 ≤ ℓ < 40, 1.3 when 40 ≤ ℓ < 60, and 0.4 when 60 ≤ ℓ < 80. We can see in Figures 4(a) and 4(b) that the proposed approach still performs efficiently better than the conventional MHE despite the variation with time of the covariance. We can also note that compared to the previous scenario, the estimate is less good for both estimators, but this can be explained by the fact that the estimators have no knowledge of this happening. Figures 4(c) and 4(d) show the trace of the matrix P ℓ used by MHE KF and MHE AD , respectively. We can observe that the MHE KF weighting matrix converges to a steady value at 20s and does not change anymore, while the proposed method tries to adapt its weighting matrix according to the variations of the process noise.
It can be seen that the trace of P ℓ tends to stabilize when the error is small. On the contrary, the trace of the weighting matrix becomes smaller when the difference between the estimate and the real state grows, ensuring good adaptation despite the influence of varying disturbances.
In Figures 4(e) and 4(f ), we can see the evolution of the forgetting factor and the dead zone, respectively.
We can see that the evolution of the forgetting factor is independent of the dead zone, but it evolves according to the variation of the process noise. We also notice that the forgetting factor tries to adapt to the variations while keeping its value close to 1. On the contrary, the behavior of the dead zone is such that when the estimation error is small, the dead zone function tends to 0 while when the error grows, the dead zone is equal to 1, so to allow a good adaptation and provide an acceptable estimation. Figure 5 shows the performance of MHE AD with different values of the horizon length. We can notice that, for small values of N, the estimation error is not negligible, while for large values of N, the estimation error is very appreciable. Table 2 allows us to assess the results quantitatively.

Example 2.
Still in the same logic of showing the performance of the proposed approach, the example of the following system is provided: Table 1: Average sum square estimation error (example 1).

MHE AD MHE KF
x (1) 44.01 83.16 x (2) 70.97 125. 10 6 Journal of Control Science and Engineering Here, disturbances are modeled with the same distribution as example 1. e constrained estimation problem is formulated with Q � I, R � 20 0 0 2 , P 0 � 10 I, and are the control input and the disturbance, respectively. H represents the Heaviside function. We also choose the horizon N � 10 in this example and d 1 � 0.00005, d 2 � 0.001, and N 0 � 5000 as design parameters. Figures 6(a)-6(c) show the states and their estimates x (1) , x (2) , and x (3) . Figures 6(d)-6(f ) show the square error of states and their estimates. From these figures, we can see that the proposed method (MHE AD ) gives a good estimate compared to the conventional MHE (MHE KF ); this is because, in this example, the inputs of the system are bounded and persistent. We can simply say the MHE AD estimates can track properly the states. e divergence of the conventional MHE is due to the inclusion of a large value of the disturbance (sudden variation) at t � 4.5 s; on the contrary, the  proposed algorithm adapts to slow or/and sudden changes.
is can easily be explained by the fact that the crucial parameters of the proposed algorithm can follow slow or/ and sudden changes or disturbances and provide good adaptability. Figures 7(a) and 7(b) show the evolution of the forgetting factor and the dead zone function, respectively. e same remarks as in the previous example are also considered in this case. It can be seen that the best performance of the proposed method occurred when the forgetting factor was equal to 1, and the dead zone function       (1) and x (2) .

MHE AD
x (1) x (2)      ) 2 (f ) Figure 6: Comparison of estimators for model (22) with horizon length N � 10. (a) Estimation of x (2) . (b) Estimation of x (2) . (c) Estimation of x (3)  Journal of Control Science and Engineering equal to 0, when the disturbances were less important. When the disturbances were considerable, the forgetting factor tried to adjust to them. Figures 8(a) and 8(b) allow appreciating the evolution of the trace of the weighting matrix P ℓ for both estimators. We can see that the evolution of the trace of the classical MHE weighting matrix has difficulty reaching the steady state unlike the previous example, which justifies the weak reconstruction; on the contrary, the evolution of the weighting matrix of MHE AD is smoother; that is, when the disturbances are less consequent, the evolution of the trace of the weighting matrix of MHE AD behaves as a classical constrained optimization problem. When the noise increases, the proposed estimator does not follow the system and the forgetting factor decreases and Table 3: Average sum square estimation error (example 2).

MHE AD
x (1) x (2) x ( moves away from 1 to allow a fairly fast adaptation of P ℓ . Based on 70 trials, the average sum square estimation error was computed, and the results are shown in Table 3. Likewise, in Figure 9, we show the performance of MHE AD with different values of the horizon length. It can be seen that MHE AD errors do not change significantly with the different horizon sizes, showing that this algorithm seems to be a good method for the arrival cost approximation, even for short horizons. Table 4 allows us to assess the results quantitatively.

Conclusion
In this paper, an alternative method to approximate the arrival cost for the moving horizon estimation has been studied on linear systems. e proposed approach is based on the least-squares algorithm but includes a variable forgetting factor which is based on the constant information principle and a dead zone which ensures robustness. e use of CasADi toolbox [30] was used for the implementation of the optimization algorithm for real-time operation. e results obtained show that the adaptive estimation algorithm allows a good approximation of the arrival cost for the moving horizon estimation. is can easily be explained by the fact that the crucial parameters of the proposed method can follow slow or/and sudden changes or disturbances and provide good adaptability. Since the update mechanism does not depend on the model, the proposed method can be extended to nonlinear systems.

Notations
FIE: Full information estimator x(t): State vector y ℓ : Measurement outputs Ψ ℓ : Objective function G: Noise weighting matrix w: Process noise vector v: Measurement noise vector f: Process vector h:

Measurement vector A:
Jacobian matrix with respect to x ℓ C: Jacobian matrix of h with respect to x ℓ N: Horizon length Γ(.): Arrival cost Q: Weighting matrix representing the confidence in the dynamic model R: Weighting matrix representing the confidence in the measurements x: Estimate vector x 0|ℓ : Estimate vector in the FIE problem x 0 : Initial state x ℓ− N : Initial state P ℓ− N : Weighting matrix λ ℓ : Forgetting factor α ℓ : Dead zone function N 0 : Design parameter e ℓ : Estimation error d 1,2 : Design parameters.

Data Availability
No data were used to support the findings of this study.