Optimal State Estimation of Nonlinear Dynamic Systems

Abstract

An optimal estimator for continuous nonlinear systems with nonlinear dynamics, and nonlinear measurement based on the continuous least square error criterion is derived. The solution is exact, explicit, in closed form and gives recursive formulas of the optimal filter. For the derivation of the filter, the following elements are combined: (i) the least squares (LS) criterion based on statistical-deterministic-likelihood approach to estimation; (ii) the state-dependent coefficient (SDC) form representation of the nonlinear system; and (iii) the calculus of variation. The resulting filter is optimal per sample. The filter’s gains need the solution of a nonsymmetric differential matrix Riccati equation. The stability of the estimator is investigated. The performances are demonstrated by simulation of the Van der Pol equation with noisy nonlinear measurement, and system driving noise.

1. Introduction

The Kalman filter and the Kalman-Bucy filter [1, 2] solved the problem of optimal estimation of stochastic and deterministic linear systems. Since then, there is a continuing research on estimation of nonlinear systems.

There are many different approaches for the state reconstruction, estimation, and filtering of nonlinear systems, for a recent review, see [3, 4] and the references within. The space in this chapter is too short to cover them. These approaches can be classified roughly into two types: the stochastic approach and the statistical-deterministic-likelihood approach.

The stochastic approach is based on the Itȏ calculus and computation of the conditional probabilities by the Kolmogorov’s forward/Fokker-Plank equation or Zakai’s equation that are difficult to solve and usually need numerical solution, e.g., see a numerical approach to the filtering problem for a class of nonlinear time-varying systems [5]. The innovations approach to the nonlinear estimation in a white noise is presented in [6]. However, explicit result for a specific nonlinear system is difficult to arrive at. Thus, when a closed-form estimator is sought, the stochastic approach leads, in general, to suboptimal and approximate solutions. The exceptions are [7, 8], where some restricted cases for which closed-form solutions of the optimal filtering equations of continuous systems are presented. Moreover, it was shown that generally the stochastic approach leads to infinite dimensional solution of the optimal estimator [9]. Different classes of nonlinear systems for which there is a closed-form explicit solution are presented in [10, 11] for the nonlinear problem of estimating the parameters of linear system with unknown coefficients. These belong to the specific class of nonlinear systems for which a general solution is presented in [12], Chapter 10.

The Kalman filter [1, 2] was obtained as well by solving the dual of the linear quadratic control problem criterion [13, 14, 15] by calculus of variations within the framework using the statistical-deterministic-likelihood approach. The dual of the LQ criterion is the least squares (LS) criterion also called the mean squares error (MSE) criterion, or joint maximum likelihood (JML) criterion [15, 16, 17], or just maximum likelihood (ML). The statistical-deterministic-likelihood approach has been used to derive filters of linear systems [13, 15]. For linear system, this approach leads to the structure of the Kalman and Kalman-Bucy filters. This shows that the Kalman and Kalman-Bucy filters are not only optimal estimators on the average but also optimal estimators for a single sample. Within the likelihood approach [18], the noises are white and the criterion is the likelihood functional [15]. The deterministic variational approach has been applied in [18] to nonlinear system. Within the statistical-deterministic-likelihood approach [13, 19], the input disturbance and output measurement error are considered as disturbances with unknown statistics ([20], p. 361). This approach is based on the calculus of variations [13] and has been widely used for numerical implicit computations of estimates and smoothers for nonlinear dynamic systems [21].

Thus, the statistical-deterministic-likelihood approach is most tempting for application in developing filters of nonlinear systems [18]. Mortensen [18] derives the general structure of the optimal recursive estimator’s state propagation equation derived from the likelihood approach point of view. This solution has the structure of the state propagation equation of the extended Kalman filter (EKF) thus justifying its usage beyond the heuristic of usage as the first-order Taylor series expansion. However, Mortensen [18] does not derive the respective equation of the gain. Moreover, Mortensen [18] states that the computation of the gain “…suffers from the same kind moment problem or closure problem as does the minimum variance nonlinear filtering.” This means that the derived estimation error gain is not feasible. The solution in this chapter shows that the statistical-deterministic-likelihood approach based on the calculus of variations leads to a solution that is not plagued with the closure problem.

The most popular estimation filter of nonlinear systems is the EKF. The EKF uses the Jacobian f_x of the system’s differential equations function $\dot{x}=f\left(x\right)$ and Jacobian m_x of the measurement’s equations y = m(x) for computation of the estimator’s gain. The stability of the EKF is not guaranteed.

An additional estimation filter of nonlinear systems that has been developed in the recent years with success is the state-dependent differential/difference Riccati equation (SDRE/SDDRE)-based filter of nonlinear system [22, 23, 24, 25]. This has been enabled by the introduction of the state-dependent coefficient (SDC) form [23, 24] approach to filtering. The SDC form represents the nonlinear equation in the quasilinear form $\dot{x}=F\left(x\right)x$ and y = M(x)x. The SDC representation always exists albeit it is not unique. The observability and controllability of the SDC representation are needed; however, for any SDC form, they are not guaranteed. Finding-synthesizing a controllable and observable SDC form representation can be difficult and is not trivial. This problem is dealt with in [26, 27, 28, 29] and some approaches to synthesize feasible SDC forms are proposed. The selection of the “best” SDC is dealt with in [26, 27, 29]. The global uniform stability properties of the SDRE-based filter have been proved only lately in [30, 31, 32, 33].

Since Mortensen’s derivation [18], no progress has been made [4, 34, 35] in explicitly solving the optimal nonlinear filtering problem till [16, 17, 36, 37, 38] for continuous nonlinear systems and [39] for discrete nonlinear systems.

This chapter combines: (i) the LS criterion based on the statistical-deterministic-likelihood approach to estimation; (ii) the SDC form representation of the nonlinear system; and (iii) the calculus of variations; for derivation of a recursive filter in the form of a differential equation as the filter-estimator for nonlinear systems with nonlinear dynamics and nonlinear measurement.

This chapter is based on the preliminary publication [16]. The results for nonlinear time-varying system are presented in [17], for system with input in [37] and for the H_∞ criterion in [38].

The presented approach leads to an optimal, exact, explicit, closed-form, and recursive solution, where state propagation equation is as derived in [18] (and is that as of the EKF). This filter is called here the recursive nonlinear least squares (RNLS) filter. The optimal gain is computed via the solution of a nonsymmetric differential matrix Riccati equation (NDMRE) that uses the respective Jacobians and the SDC form representation.

The importance and novelty of the result in this chapter are:

1. An optimal, exact, explicit, closed-form, and recursive solution to the estimation of nonlinear time-varying systems based on the quadratic least-squares criterion is presented.

2. The fact that the optimal filter of nonlinear systems can be derived by calculus of variations is highlighted.

3. The optimal filter can be taught to students that are familiar with calculus of variations before mastering stochastic calculus.

The RNLS-based filter, the EKF, and the SDDRE-based filter were compared on a common basis in [36, 40].

In the chapter, derivation of the result is presented. The performances of the RNLS-based filter are demonstrated with the Van der Pol differential equation driven by a band-limited noise, and the nonlinear measurement is noise corrupted.

2. Problem statement

A general nonlinear system is dealt with. Let the reality be represented by:

$$\begin{array}{l}\dot{\zeta }\left(t\right)=\phi \left(\zeta \left(t\right),\omega \left(t\right)\right),\zeta \left(t_o\right)={\zeta }_o\\ y\left(t\right)=\eta \left(\zeta \left(t\right),\upsilon \left(t\right)\right)\end{array}$$E1

where ζ(t) is the real state (unknown and of unknown dimension), y(t) is the measured output, υ(t) is the measurement noise, ω(t) is the system driving noise, and the functions φ and η represent the reality. The functions φ and η that describe the real system cannot be either precisely represented or are unknown precisely up to the last detail (e.g., the output measurement function may include some measurement noise or themselves exhibit random uncertain behavior). For the design of the observer, we use the representation model given by:

$$\begin{array}{l}\dot{x}\left(t\right)=f\left(x\left(t\right),w\left(t\right)\right),x\left(t_o\right)=x_o\\ y\left(t\right)=m\left(x\left(t\right),v\left(t\right)\right)\end{array}$$E2

where x(t) ∈ $\mathbb{R}$ⁿ is the state of the model, y(t) ∈ $\mathbb{R}$^p is the model output, w(t) ∈ $\mathbb{R}$^r is the system driving disturbance noise, v(t) ∈ $\mathbb{R}$^p is the measurement noise, f(.):$\mathbb{R}$ⁿ × $\mathbb{R}$^r → $\mathbb{R}$ⁿ and m(.): $\mathbb{R}$ⁿ × $\mathbb{R}$^p → $\mathbb{R}$^p are the representation (model, i.e., exactly known) of the reality and thus approximation of the reality, w(t) and v(t) are the functions of time that represent the difference between the reality and its model. It is assumed that the time functions w(t) and v(t) and the initial conditions, x_o, are of “unknown character” ([15], Section 5.3), i.e., with unknown statistics [18] ([20], p. 361).

The problem: Derive a recursive estimator (in form of a differential equation) for the state of the model, x(t), from the output measurements.

The continuous least square criterion is used [13, 14, 15] in the evaluation of the optimal estimator of linear systems. The covariance constraint and the minimum model error concepts [21] rationalize this approach as well.

The continuous least squares criterion is the dual of the LQ criterion for the control problem.

The objective is ([15], Eq. 5.24)

$$J^{\prime }\left(t\right)=\frac{1}{2}\left\{\begin{array}{l}{\left(x\left(t_o\right)-\overline{x}\left(t_o\right)\right)}^T{P}_{t_o}^{-1}\left(x\left(t_o\right)-\overline{x}\left(t_o\right)\right)\\ +\underset{t_o}{\overset{t}{\int }}\left[\begin{array}{l}{\left[y\left(\tau \right)-m(x\left(\tau \right),v\left(\tau \right))\right]}^T{R}^{-1}\left[y\left(\tau \right)-m(x\left(\tau \right),v\left(\tau \right))\right]\\ +w{\left(\tau \right)}^T{Q}^{-1}w\left(\tau \right)\end{array}\right]\mathit{d\tau }\end{array}\right\}$$E3

where Q is an a priori estimate of the driving force errors, w(t), Q ∈ $\mathbb{R}$^r × r, Q > 0, R is an a priori estimate of the measurement noise errors, v(t), R ∈ $\mathbb{R}$^p × p, R > 0, $P_{t_o}$ is an a priori covariance estimate of the initial conditions errors, $P_{t_o}$ ∈ $\mathbb{R}$^n × n, $P_{t_o}>0$, $\overline{x}\left(t_o\right)$ is an a priori estimate of the initial conditions.

We wish to minimize the continuous least squares error objective (3) with respect to w(τ), t_o ≤ τ ≤ t subject to the model (2) in order to find an estimate of x(τ). That is, we are looking for the representation-realization of the difference between the reality and the model, w(t), that best fits, the observations. In other words and roughly speaking, “we want to pass the solution to Eq. (3) as closely as possible, through the observations.” The presented approach also constitutes the statistical methods approach to filtering ([15], Section 5.3).

The problem above is solvable by a batch solution [21] that will minimize the objective (3). Here, we look for a recursive solution in the form of differential equations.

Throughout the chapter, it is assumed that all functions satisfy the necessary boundedness, smoothness, and differentiability conditions for existence of solution.

3. The state-dependent coefficient—SDC form

In this chapter, we deal with a specific structure of the model of the nonlinear system (2). It is assumed that:

1. I. Eq. (2) is partitioned as (with a slight abuse of notation):

$$\begin{array}{l}f\left(x,w\right)≕f\left(x\left(t\right)\right)+\mathit{Gw}\left(t\right);x\left(t_o\right)=x_o,\\ m\left(x,v\right)≕m\left(x\left(t\right)\right)+v\left(t\right)\end{array}$$E4

1. II. At the origin, we have

$$\begin{array}{l}f\left(0\right)=0\\ m\left(0\right)=0\end{array}$$E5

Then, by defining the state-dependent coefficient form (SDC) [23] as:

$$\begin{array}{l}f\left(x\left(t\right)\right)≕F\left(x\left(t\right)\right)x\left(t\right)\\ m\left(x\left(t\right)\right)≕M\left(x\left(t\right)\right)x\left(t\right)\end{array}$$E6

The dynamic equations of the system (4) are written as

$$\begin{array}{l}\dot{x}\left(t\right)=F\left(x\left(t\right)\right)x\left(t\right)+\mathit{Gw}\left(t\right),x\left(t_o\right)=x_o,\\ y\left(t\right)=M\left(x\left(t\right)\right)x\left(t\right)+v\left(t\right)\end{array}$$E7

where F ∈ $\mathbb{R}$^n × n, G ∈ $\mathbb{R}$^n × r, M ∈ $\mathbb{R}$^p × n. The SDC form (6) always exists albeit is not unique. It is assumed that all matrices F(ξ), M(ξ), are piecewise continuous and uniformly bounded with respect to all variables.¹ An important property of the SDC representation, that is needed, is its observability and controllability as a time-varying system along all trajectories that the RNLS filter can attain. The observability and/or controllability of a specific SDC form are not guaranteed. Finding-synthesizing a controllable and an observable SDC form representation can be difficult and is not trivial. This problem is dealt with in [26, 27, 28, 29] where some approaches to synthesize feasible SDC forms are proposed.

4. Derivation of the main result

In this section, the main result is derived for the specific structure of the nonlinear system (7), i.e., nonlinear dynamics, f(x(t)), nonlinear measurement, m(x(t)), that are represented in the SDC form given in Eq. (6), and the quadratic criterion

$$J\left(t\right)=\frac{1}{2}\left\{\begin{array}{l}{\left(x\left(t_o\right)-\overline{x}\left(t_o\right)\right)}^T{P}_{t_o}^{-1}\left(x\left(t_o\right)-\overline{x}\left(t_o\right)\right)\\ +\underset{t_o}{\overset{t}{\int }}\left[{\left[y\left(\tau \right)-m\left(x\left(\tau \right)\right)\right]}^T{R}^{-1}\left[y\left(\tau \right)-m\left(x\left(\tau \right)\right)\right]+w{\left(\tau \right)}^T{Q}^{-1}w\left(\tau \right)\right]\mathit{d\tau }\end{array}\right\}$$E8

that is minimized with respect to, w(t), subject to Eq. (7). Calculus of variations is applied in derivation of the main result for nonlinear systems (7). The Hamiltonian is

$$\begin{array}{l}H\left(x,\lambda ,t\right)=\frac{1}{2}{\left[y\left(t\right)-m\left(x\left(t\right)\right)\right]}^T{R}^{-1}\left[y\left(t\right)-m\left(x\left(t\right)\right)\right]\\ +\frac{1}{2}w{\left(t\right)}^T{Q}^{-1}w\left(t\right)-\lambda {\left(t\right)}^T\left[f\left(x\left(t\right)\right)+\mathit{Gw}\left(t\right)\right]\hfill \end{array}$$E9

where λ(t) is the costate.

The necessary conditions for optimality ([15], Example 7.11) are

$$\begin{array}{l}{H}_w=0\hfill \\ \dot{\lambda }\left(t\right)=H_x^T;\\ \lambda \left(t_o\right)=\frac{1}{2}\frac{\partial }{\partial \widehat{x}\left(t_o\right)}{\left(\widehat{x}\left(t_o\right)-\overline{x}\left(t_o\right)\right)}^T{P}_{t_o}^{-1}\left(\widehat{x}\left(t_o\right)-\overline{x}\left(t_o\right)\right)\\ \lambda \left(t\right)=0\text{since}x\left(t\right)\mathit{\text{is\hspace{0.17em}free}}\\ Q^{-1}>0,P_{t_o}^{-1}>0,R^{-1}>0\end{array}$$E10

This gives

$$\begin{array}{c}{H}_w=w{\left(t\right)}^T{Q}^{-1}-\lambda {\left(t\right)}^{T}G=0\\ \dot{\lambda }\left(t\right)={\left[{\left[y\left(t\right)-m\left(\widehat{x}\left(t\right)\right)\right]}^T{R}^{-1}\left[-m_{\widehat{x}}\left(\widehat{x}\left(t\right)\right)\right]-\lambda {\left(t\right)}^T{f}_{\widehat{x}}\left(\widehat{x}\left(t\right)\right)\right]}^T\\ \lambda \left(t_o\right)={\left(\widehat{x}\left(t_o\right)-\overline{x}\left(t_o\right)\right)}^T{P}_{t_o}^{-1}\end{array}$$E11

This leads to the nonlinear two-point boundary value problem (TPBVP) for $t_o\le \tau \le t,$

$$\widehat{w}\left(\tau \right)=Q{G}^T\lambda \left(\tau \right)$$

$$\frac{d}{\mathit{d\tau }}\widehat{x}\left(\tau \right)=f\left(\widehat{x}\left(\tau \right)\right)+{\mathit{GQG}}^T\lambda \left(\tau \right);\widehat{x}\left(t_o\right)=\overline{x}\left(t_o\right)+P_{t_o}\lambda \left(t_o\right)$$E12

$$\frac{d}{\mathit{d\tau }}\lambda \left(\tau \right)=-{\left[f_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T\lambda \left(\tau \right)-{\left[m_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T{R}^{-1}\left[y\left(\tau \right)-m\left(\widehat{x}\left(\tau \right)\right)\right];\lambda \left(t\right)=0$$

4.1. Explicit solution of the TPBVP

The system’s dynamic equation (Eq. (4)) in the SDC form is Eq. (7). Thus, the optimal solution is given by the TPBVP.

$$\widehat{w}\left(\tau \right)=Q{G}^T\lambda \left(\tau \right)$$

$$\frac{d}{\mathit{d\tau }}\widehat{x}\left(\tau \right)=F\left(\widehat{x}\left(\tau \right)\right)\widehat{x}\left(\tau \right)+{\mathit{GQG}}^T\lambda \left(\tau \right);\widehat{x}\left(t_o\right)=\overline{x}\left(t_o\right)+P_{t_o}\lambda \left(t_o\right)$$E13

$$\frac{d}{\mathit{d\tau }}\lambda \left(\tau \right)=-{\left[f_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T\lambda \left(\tau \right)-{\left[m_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T{R}^{-1}\left[y\left(\tau \right)-M\left(\widehat{x}\left(\tau \right)\right)\widehat{x}\left(\tau \right)\right];\lambda \left(t\right)=0$$E14

The usage of the SDC form converts the nonlinear TPBVP (Eq. (12)) to a time-varying TPBVP (Eq. (13)) thus enables a causal solution. This is as up to the current time, as the solution propagates forward in time, $\widehat{x}\left(t\right)$ is a known function of time and the integration goes forward in time. The solution follows [16]. For illustration, the “homogeneous” case is presented here. In this case, the TPBVP is

$$\begin{array}{l}\frac{d}{\mathit{d\tau }}\widehat{x}\left(\tau \right)=F\left(\widehat{x}\left(\tau \right)\right)\widehat{x}\left(\tau \right)+{\mathit{GQG}}^T\lambda \left(\tau \right);\widehat{x}\left(t_o\right)=\overline{x}\left(t_o\right)+P_{t_o}\lambda \left(t_o\right)\\ \frac{d}{\mathit{d\tau }}\lambda \left(\tau \right)=-{\left[f_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T\lambda \left(\tau \right)+{\left[m_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T{R}^{-1}M\left(\widehat{x}\left(\tau \right)\right)\widehat{x}\left(\tau \right);\lambda \left(t\right)=0\end{array}$$E15

By setting $\widehat{x}\left(\tau \right)=P\left(\tau \right)\lambda \left(\tau \right)$ in Eq. (15), the nonsymmetric differential matrix Riccati equation is given by:

$$\dot{P}=F\left(\widehat{x}\left(\tau \right)\right)P+P{\left[f_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T+{\mathit{GQG}}^T-P{\left[m_{\widehat{x}}\left(\widehat{x}\left(\tau \right)\right)\right]}^T{R}^{-1}M\left(\widehat{x}\left(\tau \right)\right)P,P\left(t_o\right)=P_{t_o}$$E16

The solution of the nonhomogeneous time-varying TPBVP (Eqs. (13) and (14)) is hinted by the necessary condition $\widehat{x}\left(t_o\right)=\overline{x}\left(t_o\right)+P_{t_o}\lambda \left(t_o\right)$. The derivation then follows closely [16].

4.2. The main result

The solution in the form of differential equations, the continuous recursive nonlinear least squares (RNLS) filter, is given by:

$$\begin{array}{l}\dot{\widehat{x}}\left(t\right)=f\left(\widehat{x}\left(t\right)\right)+K\left(\widehat{x}\left(t\right),t\right)\left[y\left(t\right)-m\left(\widehat{x}\left(t\right)\right)\right],\widehat{x}\left(t_o\right)={\widehat{x}}_o\\ \mathit{\text{or}}\\ \dot{\widehat{x}}\left(t\right)=F\left(\widehat{x}\left(t\right)\right)\widehat{x}\left(t\right)+K\left(\widehat{x}\left(t\right),t\right)\left[y\left(t\right)-M\left(\widehat{x}\left(t\right)\right)\widehat{x}\left(t\right)\right],\widehat{x}\left(t_o\right)={\widehat{x}}_o\end{array}$$E17

where the filter’s gain is

$$K\left(\widehat{x}\left(t\right),t\right)=P\left(\widehat{x}\left(t\right),t\right){\left[m_{\widehat{x}}\left(\widehat{x}\left(t\right)\right)\right]}^T{R}^{-1}$$E18

and $P\left(\widehat{x}\left(t\right),t\right)$ is given by the nonsymmetric differential matrix Riccati equation

$$\begin{array}{l}\dot{P}\left(\widehat{x}\left(t\right),t\right)=F\left(\widehat{x}\left(t\right)\right)P\left(\widehat{x}\left(t\right),t\right)+P\left(\widehat{x}\left(t\right),t\right){\left[f_{\widehat{x}}\left(\widehat{x}\left(t\right)\right)\right]}^T\hfill \\ +{\mathit{GQG}}^T-P\left(\widehat{x}\left(t\right),t\right){\left[m_{\widehat{x}}\left(\widehat{x}\left(t\right)\right)\right]}^T{R}^{-1}M\left(\widehat{x}\left(t\right)\right)P\left(\widehat{x}\left(t\right),t\right);P\left(t_o\right)=P_{t_o}\hfill \end{array}$$E19

where $\widehat{x}\left(t\right)$ is the estimated state and $f\left(x\right)≕F\left(x\right)x$, $f_x\left(x\right)=\frac{\partial f\left(x\right)}{\partial x}$, $m\left(x\right)≕M\left(x\right)x$, $m_x\left(x\right)=\frac{\partial m\left(x\right)}{\partial x}$.

Notice:

1. The first term of the right-hand side of Eq. (19) includes the SDC form and the second term includes the Jacobian and same is in the last term. The SDC and the Jacobian are equal for linear systems only.

2. $\widehat{x}\left(t\right)$ is known up to the current time t. Thus, Eq. (17) can be propagated forward in time.

3. The solution requires the solution of the nonsymmetric differential matrix Riccati equation (Eq. (19)) and the solution, P, is nonsymmetric.

4. The solution of the nonsymmetric Riccati matrix equation depends on the estimated state $\widehat{x}\left(t\right),$ and is formally denoted $P\left(\widehat{x}\left(t\right),t\right)$.

5. Notice that the state propagation Eq. (19) has exactly the same structure as derived by Mortensen [18] and used by the EKF. The solution of Eqs. (18) and (19) gives explicitly the filter’s gain.

6. In [18], it is claimed that computation of the filter’s optimal gain, P, (Eqs. (18) and (19)) suffers from “…the moment or closure problem…”. In this chapter, it is shown that the filter’s optimal gain is solved completely and explicitly by the NDMRE (Eq. (19)).

4.3. A compact form of the optimal solution

In order to enable better understanding of Eq. (17–19), the following presents Eq. (17–19) by suppressing the explicit and implicit dependence on time.² The optimal filter is

$$\dot{\widehat{x}}=F\widehat{x}+K\left[y-M\widehat{x}\right],\widehat{x}\left(t_o\right)={\widehat{x}}_o$$E20

$$K={\mathit{Pm}}_{\widehat{x}}^T{R}^{-1}$$E21

$$\dot{P}=\mathit{FP}+{\mathit{Pf}}_{\widehat{x}}^T+{\mathit{GQG}}^T-{\mathit{Pm}}_{\widehat{x}}^T{R}^{-1}\mathit{MP};P\left(t_o\right)=P_{t_o}$$E22

or

$$\dot{P}=\mathit{FP}+P{\left[f_{\widehat{x}}\right]}^T+{\mathit{GQG}}^T-\mathit{KMP};P\left(t_o\right)=P_{t_o}$$E23

One can clearly see that for linear system, Eq. (20–22) gets the structure of the Kalman filter as then $F=f_{\widehat{x}}$ and $M=m_{\widehat{x}}$.

5. Stability analysis of the RNLS estimator

The deterministic stability of the RNLS estimator-filter along the filter’s trajectories Eq. (20–22) is considered. Recall that optimality does not guarantee stability. The stability of the RNLS filter is connected to the stability of the NDMRE equation. Let us consider the system/observer:

$$\dot{\widehat{x}}=F\widehat{x}+K\left[y-M\widehat{x}\right]=\left[F-\mathit{KM}\right]\widehat{x}+\mathit{Ky},\widehat{x}\left(t_o\right)={\widehat{x}}_o$$E24

$$K={\mathit{PM}}^T{R}^{-1}$$E25

$$\dot{P}=\mathit{FP}+{\mathit{PF}}^T+{\mathit{GQG}}^T-{\mathit{PM}}^T{R}^{-1}\mathit{MP};P\left(t_o\right)=P_{t_o}$$E26

where the explicit and implicit time dependency is suppressed as in the previous section. Eqs. (24–26) are actually the deterministic SDDRE-based observer of Eq. (7) whose stability is treated in [31, 32, 33]. The matrix Riccati equation (Eq. (26)) is symmetric.

First, existing result on the stability of optimal estimators of system Eqs. (24–26) as a linear time-varying system is cited. The following result is valid for linear time-invariant and time-variant systems.Theorem 1.

[31, 32, 41] Consider the symmetric Riccati equation (Eq. (26)) where $Q\ge 0,R>0$ and P_o ≥ 0 are symmetric, $\left(F,M\right)$ is detectable, and $\left(F,{\mathit{GQ}}^{1/2}\right)$ is stabilizable. Then, there exists $K={\mathit{PM}}^T{R}^{-1}$ such that $F-\mathit{KM}$ is asymptotically stable.

A Lyapunov function for the autonomous system Eqs. (24–26) (i.e. y = 0) is

$$V\left(t\right)=\frac{1}{2}\widehat{x}{\left(t\right)}^T{P}^{-1}\left(t\right)\widehat{x}\left(t\right)$$E27

For which

$$\dot{V}\left(t\right)=-\widehat{x}{\left(t\right)}^T{P}^{-T}\left[{\mathit{GQG}}^T+{\mathit{PM}}^T{R}^{-1}\mathit{MP}\right]P^{-1}\widehat{x}\left(t\right)$$E28

where $\left[{\mathit{GQG}}^T+{\mathit{PM}}^T{R}^{-1}\mathit{MP}\right]$ is positive definite.

Next, the NDMRE equation is considered. It is dealt with in [42, 43, 44, 45, 46]. The only reference that is directly addressing the stability issue of an NDMRE is [42] (Chapter 9). The Riccati equation related to the time-invariant control problem is dealt with in [42] (Theorem 9.1.23 and Remark 9.1.24). Although not explicitly stated, these results apply as well to time-varying systems. Motivated by this theorem and remark, translated by duality to the estimation problem, the following conjecture is formulated.Conjecture 1.

Consider the nonsymmetric differential Riccati matrix equation.

$$\dot{P}=\mathit{FP}+{\mathit{Pf}}_{\widehat{x}}^T+{\mathit{GQG}}^T-{\mathit{Pm}}_{\widehat{x}}^T{R}^{-1}\mathit{MP},P\left(t_o\right)=P_o,$$E29

where $Q\ge 0,R>0$ are symmetric, $\left(F,M\right)$ and $\left(f_{\widehat{x}},m_{\widehat{x}}\right)$ are detectable, and $\left(F,{\mathit{GQ}}^{1/2}\right)$ and $\left(f_{\widehat{x}},{\mathit{GQ}}^{1/2}\right)$ are stabilizable. Then, there exist $K_1={\mathit{PM}}^T{R}^{-1}$ and $K_2={\mathit{Pm}}_{\widehat{x}}^T{R}^{-1}$ such that $F-K_{1}M$ and $f_{\widehat{x}}-K_2{m}_{\widehat{x}}$ are stable.

This conjecture is supported by [42] (Theorem 9.1.23 and Remark 9.1.24). The requirement of detectability (observability) and stabilizability (controllability) is not explicitly required in [42] (supposedly they appear implicitly). This conjecture means that the filter given by Eqs. (20–22) is stable. An issue under research is (loosely): in addition to the conditions in Conjecture 1, the boundedness conditions of all matrices and variables (the output and system driving noise and measurement noise) are sufficient conditions for this stability, as for the SDDRE-based filter [31, 32, 33]?

Notice that for the symmetric case, this well-known result for linear system results in Theorem 1.

The stability of the RNLS filter is investigated via Lyapunov analysis. As the solution of the nonsymmetric Riccati equation in Eq. (19) is eventually not symmetric, the following symmetric Lyapunov function is dealt with here:

$$V=\frac{1}{2}{x}^T\left(P^{-1}+P^{-T}\right)x$$E30

The derivative of the Lyapunov function is [47]

$$\dot{V}=-\frac{1}{2}{x}^T\left[\begin{array}{l}{P}^{-T}P\left(M^T{R}^{-1}M+m_x^T{R}^{-1}{m}_x\right)\\ -P^{-T}P{\left(M-m_x\right)}^T{R}^{-1}\left(M-m_x\right)\\ +P^{-1}{\mathit{GQG}}^T{P}^{-1}+P^{-T}{\mathit{GQG}}^T{P}^{-T}\\ +{\left(f_x-F\right)}^T{P}^{-1}+P^{-T}\left(f_x-F\right)\end{array}\right]x$$E31

For linear system, $F=f_{\widehat{x}},M=m_{\widehat{x}}$, we have Eq. (28).

The first terms in Eq. (31) are potentially nonnegative definite

$$P^{-T}P\left(M^T{R}^{-1}M+m_x^T{R}^{-1}{m}_x\right)\ge 0$$E32

The second term in Eq. (31) is negative (nonpositive) definite

$$-P^{-T}P{\left(M-m_x\right)}^T{R}^{-1}\left(M-m_x\right)\ge 0$$E33

The next two terms in Eq. (31) are indefinite and can be negative

$$P^{-1}{\mathit{GQG}}^T{P}^{-1}+P^{-T}{\mathit{GQG}}^T{P}^{-T}$$E34

The last two terms in Eq. (31)

$${\left(f_x-F\right)}^T{P_{\alpha }}^{-1}+{P_{\alpha }}^{-T}\left(f_x-F\right)$$E35

are indefinite.

The discussion above hints that for small nonsymmetry, for sure, the NDMRE stabilizes the RNLS filter. The stability of the RNLS filter is summarized in the following conjecture. Further results are beyond the scope of this chapter.Conjecture 2:

If

1. The nonlinearities are such that $‖f_x-F‖$ and $‖m_x-M‖$ are bounded/uniformly bounded and sufficiently small,

2. The observability and controllability conditions are satisfied along the filter trajectories, then the RNLS filter is asymptotically stable.

Remark: Simulation results show/hint that as long as the incremental observability matrices

$$\mathit{Ob}\left(F\left(x\right),M\left(x\right)\right)=\left[\begin{array}{l}M\left(x\right)\\ M\left(x\right)F\left(x\right)\\ ⋮\\ M\left(x\right)F{\left(x\right)}^{n-1}\end{array}\right],\mathit{Ob}\left(f_x\left(x\right),m_x\left(x\right)\right)=\left[\begin{array}{l}{m}_x\left(x\right)\\ m_x\left(x\right)f_x\left(x\right)\\ ⋮\\ m_x\left(x\right)f_x{\left(x\right)}^{n-1}\end{array}\right]$$

and the incremental controllability matrices

$$\begin{array}{l}\mathit{Co}\left(F\left(x\right),{\mathit{GQ}}^{1/2}\right)=\left[{\mathit{GQ}}^{1/2}F\left(x\right){\mathit{GQ}}^{1/2}\cdots F{\left(x\right)}^{n-1}{\mathit{GQ}}^{1/2}\right],\\ \mathit{Co}\left(f_x\left(x\right),{\mathit{GQ}}^{1/2}\right)=\left[\begin{array}{llll}{\mathit{GQ}}^{1/2}& f_x\left(x\right){\mathit{GQ}}^{1/2}& \cdots & f_x{\left(x\right)}^{n-1}{\mathit{GQ}}^{1/2}\end{array}\right]\end{array}$$

along the estimator’s trajectory of the RNLS filter are nonsingular, i.e.,

$$\begin{array}{l}\mathit{rank}\left[\mathit{Ob}\left(F\left(x\right),M\left(x\right)\right)\right]=n,\mathit{rank}\left[\mathit{Ob}\left(f_x\left(x\right),m_x\left(x\right)\right)\right]=n,\mathit{rank}\left[\mathit{Co}\left(F\left(x\right),{\mathit{GQ}}^{1/2}\right)\right]=n,\text{and}\\ \mathit{rank}\left[\mathit{Co}\left(f_x\left(x\right),{\mathit{GQ}}^{1/2}\right)\right]=n\end{array},$$

then: (i) the estimation errors of the filter for the deterministic case, i.e., w(t) = 0 and v(t) = 0, converge to zero; and (ii) for the case with bounded disturbance and bounded measurement noise, the estimation errors are bounded, i.e., do not diverge.

6. Example

This section demonstrates the performance of the RNLS-based estimator on a generalized nonlinear time-varying Van der Pol differential equation driven by band-limited noise and noise-corrupted nonlinear measurement. The state is $x={\left[\begin{array}{ll}\xi & \dot{\xi }\end{array}\right]}^T$ interpreted as position and velocity. The Van der Pol equation is

$$\mu \ddot{\xi }+2c\left({\xi }^2-1\right)\dot{\xi }+\mathit{k\xi }=w$$

That can be put in matrix form as:

$$\frac{d}{\mathit{dt}}\left[\begin{array}{l}\xi \\ \dot{\xi }\end{array}\right]=\left[\begin{array}{ll}0& 1\\ -\frac{k}{\mu }& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]\left[\begin{array}{l}\xi \\ \dot{\xi }\end{array}\right]+\left[\begin{array}{l}0\\ 1\end{array}\right]w$$

The noisy measurement is

$$y=\frac{\xi }{\sqrt{1+{\xi }^2}}+v$$

Then, we have

$$f\left(x\right)=\left[\begin{array}{ll}0& 1\\ -\frac{k}{\mu }& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]\left[\begin{array}{l}\xi \\ \dot{\xi }\end{array}\right]=\left[\begin{array}{l}\dot{\xi }\\ -\frac{\mathit{k\xi }}{\mu }-\frac{2c}{\mu }\left({\xi }^2-1\right)\dot{\xi }\end{array}\right]$$

The SDC form system matrix is selected as:

$$F\left(x\right)=\left[\begin{array}{ll}0& 1\\ -\frac{k}{\mu }& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]$$

and the Jacobian is

$$f_x\left(x\right)=\left[\begin{array}{ll}0& 1\\ -\left(\frac{k}{\mu }+\frac{4c}{\mu }\xi \dot{\xi }\right)& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]$$

$$\begin{array}{l}m\left(x\right)=\frac{\xi }{\sqrt{1+{\xi }^2}}\\ M\left(x\right)=\left[\begin{array}{ll}\frac{1}{\sqrt{1+{\xi }^2}}& 0\end{array}\right]\\ m_x\left(x\right)=\left[\begin{array}{ll}\frac{1}{{\left(1+{\xi }^2\right)}^{3/2}}& 0\end{array}\right]\end{array}$$

The observability matrices are

$$\mathit{Ob}\left(F\left(x\right),M\left(x\right)\right)=\left[\begin{array}{ll}\frac{1}{\sqrt{1+{\xi }^2}}& 0\\ 0& \frac{1}{\sqrt{1+{\xi }^2}}\end{array}\right],$$

$$\mathit{Ob}\left(f_x\left(x\right),m_x\left(x\right)\right)=\left[\begin{array}{ll}\frac{1}{{\left(1+{\xi }^2\right)}^{3/2}}& 0\\ 0& \frac{1}{{\left(1+{\xi }^2\right)}^{3/2}}\end{array}\right]$$

and controllability matrices are

$$\mathit{Co}\left(F\left(x\right),{\mathit{GQ}}^{1/2}\right)=\left[\begin{array}{ll}0& 1\\ 1& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]Q^{1/2}$$

$$\mathit{Co}\left(f_x\left(x\right),{\mathit{GQ}}^{1/2}\right)=\left[\begin{array}{ll}0& 1\\ 1& -\frac{2c}{\mu }\left({\xi }^2-1\right)\end{array}\right]Q^{1/2}$$

The observability and controllability matrices have full rank for all bounded trajectories.

The system and the RNLS estimator were implemented in SIIMULINK^® with the following parameters:

Ts = 0.1 msec	Sampling interval
μ = 1	Mass
c = 0.01	Damping coefficient
k = 0.1	Spring stiffness
R = 1e-5 [1/Hz]	Spectral density of the measurement noise—v
Q = 1e0 [(1/sec^2)^2/Hz]	Spectral density of the system driving noise—w
Po = [0.001 0; 0 0.001]	Initial condition of the P matrix
x(to) = [2 0]^T	Initial conditions of the state

The measurement noise and system driving noises are white in 100 [rad/sec] bandwidth.

The following figures present the performances of the RNLS filter. Figure 1 presents the measured output—y and the estimated output versus time. Figure 2 presents the real (true) position—ξ and the estimated position—$\widehat{\xi }$ versus time. Figure 3 presents the real (true) velocity—$\dot{\xi }$ and the estimated velocity—$\dot{\widehat{\xi }}$ versus time. The transient performance is demonstrated. Figure 4 presents the filter’s gains: gain of the position state, K1, and the gain of the velocity state, K2, versus time. Figure 5 shows the solution of the Riccati equation matrix, P, versus time. One can clearly see that the P matrix is nonsymmetric P12 ≠ P21.

Figure 1.

The measured output—y and the estimated output versus time.

Figure 2.

The real position—x and the estimated position state—$\widehat{x}$ versus time.

Figure 3.

The real velocity—$\dot{x}$ and the estimated velocity—$\widehat{\dot{x}}$, versus time.

Figure 4.

Filter’s gains, K1 gain of the position state, K2 gain of the velocity state, versus time.

Figure 5.

The terms of the solution of the Riccati equation—P matrix, versus time.

Figure 6 presents the phase plane plot of the velocity estimation error versus the position estimation errors. One can see that following the initial transient, the estimation errors concentrate around the origin.

Figure 6.

Phase plane plot of velocity versus position estimation errors.

7. Conclusions

The mean least square error criterion has been used to derive the optimal estimator for continuous nonlinear systems with nonlinear dynamics and nonlinear measurement. The solution is exact, explicit, in closed form, and in recursive form. Simulation example demonstrates the performance.

Acknowledgments

This research was partially supported by the Prof. Pazy Research Foundation.