- Email: [email protected]
Optimization Based Constrained Unscented Gaussian Sum Filter
Krist V. Gernaey, Jakob K. Huusom and Rafiqul Gani (Eds.), 12th International Symposium on Process Systems Engineering and 25th European Symposium on Computer Aided Process Engineering. c 2015 Elsevier B.V. All rights reserved. 31 May - 4 June 2015, Copenhagen, Denmark.
Optimization Based Constrained Unscented Gaussian Sum Filter Krishna Kumar Kottakki, Sharad Bhartiya, Mani Bhushan Department of Chemical Engineering; Indian Institute of Technology Bombay; Mumbai, India [email protected], [email protected], [email protected]
Abstract Bayesian state estimation for constrained nonlinear systems involves two main issues: (i) appropriate representation of the underlying non-Gaussian densities, and (ii) use of appropriate techniques to bound the estimates so that they satisfy the give constraints. In this work we present a new constrained state estimation approach, named as Optimization based Constrained Unscented Gaussian Sum Filter. The proposed approach uses Unscented Gaussian Sum Filter to represent the non-Gaussian densities and an optimization approach to incorporate the state constraints. Benefits of the proposed approach are illustrated by applying it to a benchmark case study. Keywords: Unscented transformation, Unscented Gaussian sum filter, Sum of Gaussians
1. Introduction Recursive Bayesian estimation for nonlinear dynamical systems uses process and measurement models along with available process measurements to obtain conditional density of states. These process models are often derived from first principles, that is material and energy balances, and are typically nonlinear. Further, the process states inherently satisfy thermodynamic properties/constraints, for e.g. mole fraction of individual species/components in a multi component mixture are positive and their sum is unity. Hence, a practical nonlinear constrained state estimation approach should have the following two features: (i) an ability to approximate non-Gaussian densities arising out of nonlinear transformation or non-Gaussian noise, and (ii) an ability to incorporate constraints on the state estimates. Several nonlinear state estimation approaches such as Extended Kalman Filter (EKF) (Anderson and Moore, 1979), Unscented Kalman Filter (UKF) (Julier and Uhlmann, 2004), Ensemble Kalman Filter (EnKF) (Gillijns et al., 2006), Gaussian Sum Filters (Sorenson and Alspach, 1971) and Particle Filters (Arulampalam et al., 2002) have been reported in literature that account for the nonlinear transformations and/or the non-Gaussian nature of the densities. While EKF involves linearization of the nonlinear process and measurement models, other approaches utilize deterministically or randomly chosen samples thereby avoiding the explicit linearization. Among these sampling based approaches, UKF is based on unscented transformation (UT) that requires only (2n + 1) deterministically chosen sigma points, with n being the number of states. This use of limited number of samples opens up possibilities of applying UKF for solving practical (industrial) state estimation problems (Kolas et al., 2009; Vachhani et al., 2006). However, similar to Kalman Filter (Kalman, 1960), UKF also assumes that the underlying probability densities are Gaussian, an assumption commonly violated for nonlinear systems and/or non-Gaussian noise. In literature, a few modifications of UKF have been proposed that avoid the Gaussianity assumpˇ tion. In particular, Gaussian-sum UKF (GS-UKF) (Simandl and Dun´ık, 2006) is based on the
1716
Kottakki et al.
result that a sum of Gaussians can approximate any density to an arbitrary degree of accuracy (Sorenson and Alspach, 1971). GS-UKF uses N(2n+1) sigma points to obtain a sum of Gaussians based representation of non-Gaussian prior, where N refers to the number of Gaussians used and is a user defined parameter. Recently Unscented Gaussian Sum Filter (UGSF) (Kottakki et al., 2014) has been proposed that enables a sum of Gaussians representation of the non-Gaussian prior using only 2n+1 sigma points. In particular, the computational effort in UGSF is similar to that in UKF and it can thus be potentially applied to practical systems (Kottakki et al., 2014). The sum of Gaussians prior in UGSF is subsequently updated using Bayes’ rule. The updated Gaussian sum posterior in UGSF is reapproximated as a single Gaussian at the next time instant using the first two moments of the updated Gaussian sum density. This reapproximation avoids degeneracy of the individual Gaussians in the sum of Gaussians representation. Implementation on several case studies demonstrated superior performance of UGSF when compared to UKF. While UGSF addresses the issue of non-Gaussianity of the prior density, it does not incorporate constraints on the state estimates making them potentially infeasible and thus unusable. In this work, we propose an extension of the UGSF approach so that it can effectively incorporate state constraints in the estimation procedure. This is achieved by modifying the sigma point generation step to incorporate constraints as well as using explicit constrained optimization formulation to obtain the posterior density. The resulting approach is labeled as Optimization Based Constrained Unscented Gaussian Sum Filter (OCUGSF). In particular, the proposed OCUGSF solves 2n+1 optimization problems to obtain a constrained posterior estimate. Organization of this paper is as follows: Problem statement with a brief summary of elements used in the proposed method is presented in Section 2 and Section 3, respectively. Benefits of the proposed OCUGSF are demonstrated by comparing its performance with other state estimation techniques on a benchmark three state isothermal batch reactor case study in Section 4. The paper is concluded in Section 5.
2. Problem Statement Consider a sampled data system consisting of nonlinear process dynamics, a linear measurement function and interval constraints on the states as, tk
x(tk ) = x(tk−1 ) +
f (x(t), u(t))dt + wk , wk ∼ N (0, Q)
(1)
tk−1
yk = Hx(tk ) + vk , vk ∼ N (0, R) dk ≤ xk ≤ ek
(2) (3)
where, x(t) ∈ Rn , u(t) ∈ R p , represent the state and input vectors at time t while yk ∈ Rm , wk ∈ Rn , vk ∈ Rm represent observation, state noise and measurement noise, respectively at time tk . Function f : Rn × R p → Rn represents the nonlinear state dynamics and H ∈ Rmxn represents the linear observation model. Measurements yk are assumed to be available at regularly spaced sampling instants tk , at k = 0, 1, 2, 3, . . . with Ts = tk −tk−1 being the sampling interval. For ease of notation, we define xk x(tk ). Eq.(3) specifies the interval constraints on each component of state vector x. The problem is to find a point estimate for xk , whose dynamics are described in Eq.(1), using available measurements Yk = y1 , . . . , yk , in Eq.(2) subject to the interval constraints given in Eq.(3). We limit our attention to the linear measurement model Eq.(2) for ease of presentation. For nonlinear measurement models, the UGSF approach can be extended (Kottakki et al., 2014) using ideas presented in the current work. Now we briefly review the proposed OCUGSF for the given system (Eqs.(1), (2)) and the state constraints (Eq.(3)).
Optimization Based Constrained Unscented Gaussian Sum Filter
1717
3. Optimization Based Constrained Unscented Gaussian Sum Filter (OCUGSF) The proposed OCUGSF approach uses a three step framework to obtain constrained state estimate for the problem (Eqs.(1) to (3)). These steps are: (i) use of Interval Constrained Unscented Transformation (ICUT) approach (Vachhani et al., 2006) to generate the constrained sigma points from the given posterior density, (ii) propagation of sigma points through the process model, and (iii) a constrained optimization based update step. These are discussed next. 3.1. Interval Constrained Unscented Transformation (ICUT) Based Sigma Point Generation: Given that the random variable xk−1 has a Gaussian distribution with mean xˆk−1|k−1 and covariance Pk−1|k−1 at time tk−1 , the Interval Constrained Unscented Transformation approach deterministically selects 2n+1 constrained sigma points and corresponding weights as (Vachhani et al., 2006), ⎧ ⎪ ⎨ xˆk−1|k−1 , (i),c xˆk−1|k−1 + θi,k−1 [ Pk−1|k−1 ]i , χk−1|k−1 = ⎪ ⎩ xˆ k−1|k−1 − θi,k−1 [ Pk−1|k−1 ]i−n ,
w(i),c = bk−1 , i=0 w(i),c = ak−1 + θi,k−1 bk−1 , i = 1, . . . , n w(i),c = ak−1 + θi,k−1 bk−1 , i = n + 1, . . . , 2n (4)
ak−1 =
(2κ − 1)
2n √ 2(n + κ) ∑ θi,k−1 − (2n + 1) n + κ
√ 1 − ak−1 n + κ , bk−1 = 2(n + κ)
(5)
i=1
θi,k−1 = min (Θ j,i ), i = 1, 2, . . . , 2n 1≤ j≤n
where, for j = 1, 2, . . . , n; i = 1, 2, . . . , 2n ⎧ √ n+ if S j,i = 0 ⎪ ⎪ κ ⎪ ⎪ √ e j,k−1 −xˆ j,k−1|k−1 ⎨ n + κ, min if S j,i > 0 S( j,i) Θ j,i ⎪ ⎪ √ d −xˆ ⎪ ⎪ n + κ, j,k−1 S( j,i)j,k−1|k−1 if S j,i < 0 ⎩ min with S
Pk−1|k−1
− Pk−1|k−1
(6)
(7)
(8)
(i),c
In Eq.(4), χk−1|k−1 , w(i),c represent the ith constrained sigma point and its corresponding weight, respectively. It can be noted that these sigma points have been projected away from the mean using a factor θi,k−1 . Since this factor has been chosen according to the state constraints (Eqs.(6) to (8)), the sigma points will be feasible with respect to the interval constraints Eq.(3). Further, if the constraints are not active (dk = −∞ and ek = +∞), the expressions (Eqs.(4) to (8)) result in unconstrained sigma points obtained using UT as have been reported in the literature (Julier and Uhlmann, 2004; Vachhani et al., 2006). 3.2. Prediction Step: The constrained sigma points obtained in Eq.(4) are then transformed through the process model (Eq.(1)) to obtain the predicted sigma points as follows: (i),c
(i),c
χk|k−1 = χk−1|k−1 +
tk
f (x(t), u(t))dt, i = 0, . . . , 2n
(9)
tk−1
Without loss of generality, it is assumed that since the sigma points generated by Eq.(4) satisfy the (i),c interval constraints, the predicted sigma points χk|k−1 (Eq.(9)) will also satisfy these constraints. The proposed OCUGSF uses the UGSF design choices (Kottakki et al., 2014) to obtain nonGaussian prior as a sum of Gaussians. These design choices are: (i) Number of Gaussians in prior
1718
Kottakki et al. (i),c
is the same as the number of sigma points i.e 2n + 1, (ii) χk|k−1 is the mean of ith Gaussian in the sum of Gaussians prior, (iii) the covariance of each ith Gaussian, i = 0, . . . , 2n is identical and equal to the process noise covariance Q and (iv) weight w(i),c of the ith sigma point is the weight of ith Gaussian, i = 0, . . . , 2n in the sum of Gaussians. These design choices result in a constrained non-Gaussian prior: pxk |Yk−1 (ξk |Yk−1 ) =
2n
1
∑ w(i),c (2π)n/2 |Q|1/2 exp
#
i=0
$ −1 (i),c (i),c [ξk − χk|k−1 ]T Q−1 [ξk − χk|k−1 ] 2
(10)
Now we illustrate a procedure to obtain a constrained posterior corresponding to the prior (Eq.(10)), measurement model (Eq.(2)) and the state constraints (Eq.(3)). 3.3. Main Contribution: Constrained Update Using Optimization For the linear measurement model with Gaussian uncertainties (Eq.(2)), the likelihood density is given as, pyk |xk ,Yk−1 (ζk |ξk ,Yk−1 ) =
#
1 (2π)m/2 |R|1/2
exp
$ −1 [ζk − Hξk ]T R−1 [ζk − Hξk ] 2
(11)
where ζk ∈ Rm is the random variable corresponding to the measurement yk . Since the above likelihood density is an unconstrained Gaussian, direct application of Bayes’ rule with the nonGaussian prior (Eq.(10)) and likelihood densities (Eq.(11)), can yield an unconstrained sum of Gaussians posterior density (Vachhani et al., 2006). To avoid this situation, we propose to update the means of the individual Gaussians in the sum of Gaussians prior by using a constrained optimization formulation. This formulation replaces the conventional Bayes’ rule based update step (i),c as presented by Kottakki et al. (2014). Thus the constrained updated mean (χk|k ) and updated covariances (Ξ ) of ith individual Gaussians are given by, (i),c
χk|k = arg min
χk :dk ≤χk ≤ek
(i),c (i),c (yk − Hχk )T R−1 (yk − Hχk ) + (χk − χk|k−1 )T Q−1 (χk − χk|k−1 )
Ξ = Q − QH T [HQH T + R]−1 HQ
(12) (13)
(i)
while the updated weights δk are given by, (i),c
w˜ k
(i),c
δk
# $ 1 (i,c) (i,c) = w(i),c exp − (yk − Hχk|k−1 )T [HQH T + R]−1 (yk − Hχk|k−1 ) , i = 0, . . . , 2n 2 (i),c
= w˜ k
2n
/ ∑ w˜ k
( j),c
, i = 0, . . . , 2n
(14) (15)
j=0
Conventional Gaussian sum filters often encounter degeneracy in weights in subsequent iterations (Sorenson and Alspach, 1971). The proposed OCUGSF overcomes this problem by approximating the sum of Gaussians posterior density as a single Gaussian for subsequent step of sigma point generation at the next time instant (Kottakki et al., 2014). The mean and covariance of this single Gaussian are chosen to be the mean and covariance of the sum of Gaussians posterior density which are given as: c xˆk|k =
2n
∑
i=0
%
& (i),c (i),c χk|k
δk
c = , Pk|k
2n
%
∑ δk
i=0
(i),c
& (i),c
(i),c
c c T ][χk|k − xˆk|k ] Ξ + [χk|k − xˆk|k
(16)
To summarize, Eqs.(4) to (16) illustrate the steps involved in the proposed OCUGSF approach for obtaining the constrained posterior moments at kth time instant for the given posterior moments at (k − 1)th time instant.
1719
Optimization Based Constrained Unscented Gaussian Sum Filter
It can be shown that if the state constraints are not active, then the proposed OCUGSF approach reduces to UGSF approach. Further, if the non Gaussian prior density (Eq.(10)) is approximated by a single Gaussian instead of sum of Gaussians, then for the unconstrained and the constrained scenarios, OCUGSF reduces to the UKF (Julier and Uhlmann, 2004) and the Unscented Recursive Nonlinear Dynamic Data Reconciliation (URNDDR) (Vachhani et al., 2006; Narasimhan and Rengaswamy, 2009) approaches, respectively.
4. Results and Discussions A three state isothermal batch process is chosen to validate the performance of proposed OCUGSF. The process involves gas-phase reversible reactions under isothermal conditions, A B+C, 2B C. The dynamical process and measurement models are (Kolas et al., 2009): ⎡ ⎤ ⎡ ⎤ −kc1 x1 + kc2 x2 x3 x˙1 2 ⎣x˙2 ⎦ = ⎣kc1 x1 − kc2 x2 x3 − 2kc3 x + 2kc4 x3 ⎦ 2 x˙3 kc1 x1 − kc2 x2 x3 + kc3 x22 − kc4 x3 T yk = 32.84 32.84 32.84 x1 x2 x3
(17) (18)
Eqs.(17) and (18) represent the nonlinear process model and linear measurement model, respectively. In Eq.(17), x1 , x2 and x3 represent the concentrations of the species A, B and C, respectively, and kci (i = 1, . . . , 4) represents the rate constant for ith reaction. True states were generated for 100 different noise realizations consisting of 720 sampling time instants, i.e. for 180 seconds. The parameters used for true state generation and for state estimation are reported in Table 1. Table 1: Process and simulation parameters for isothermal batch reactor (Kolas et al., 2009) x0 = 0.5 0.05 0 Δt = 0.25seconds 1 =x xˆ0|0 0
kc = 0.5 0.05 0.2 0.01 Q = diag(0.0012 0.0012 0.0012 ) T dk = 0.0129 0.0743 0.0243
2 xˆ0|0
P = 4I3×3 0 T = 0.6 0.1 0.05
ek = 0.4758
0.3575
T
0.6601
The performance of estimation algorithms, namely, OCUGSF, URNDDR, UGSF, and UKF are 1 and investigated for two different scenarios of initial posterior densities, i.e Case I: xˆ0|0 = xˆ0|0 2 P0|0 = P0 , and Case II: xˆ0|0 = xˆ0|0 and P0|0 = P0 . Average of sum of the square roots of the estimation errors (ASSREE) was chosen %' to compare the performance& of estimation algorithms and is defined as: ASSREE =
1 100
100
∑
l=1
1 720
720
(l) (l) (l) (l) (l) (l) ∑ [xˆk|k − xk ]T [xˆk|k − xk ] . Here xk and xˆk|k represent k=1 l th simulation run, respectively. The ASSREE values for
the true and the estimated states for the four estimation algorithms along with their computational times are reported in Table 2 for both Cases I and II. The average number of constraint violations is also reported under the heading #violations in this Table. From the ASSREE values reported in Table 2, it can be observed that the proposed OCUGSF results in superior state estimates compared to the other three algorithms with higher, but acceptable computational cost. It can also be noted that OCUGSF estimates have satisfied the state constraints for all the time instants across the simulation runs. Even though the URNDDR approach has reported no constraint violation, the Gaussian approximation of nonGaussian prior (Narasimhan and Rengaswamy, 2009) has resulted in inferior performance compared to OCUGSF. For UGSF, while the prior approximation is non-Gaussian, it has violated the state constraints at several time instances as it is an unconstrained estimator. Gaussian approximation of non-Gaussian prior and unconstrained state estimation has resulted in inferior performance of UKF. State x2 tracking for one particular simulation run has been shown in Fig. 1 and illustrates the superior performance of OCUGSF over UKF, UGSF and URNDDR approaches.
1720
Kottakki et al.
Table 2: Performance comparison for three state isothermal batch reactor
Algorithm UKF UGSF URNDDR OCUGSF
Case-I (xˆ0|0 = x0 ) ASSREE #violations CPU Seconds 159.39 705 0.0271 5.9164 645 0.0273 3.5503 0 0.1131 2.1247 0 0.1400
ASSREE 155.29 6.0094 3.4396 2.7890
2 ) Case-II (xˆ0|0 = xˆ0|0 #violations CPU Seconds 705 0.0320 643 0.0322 0 0.1781 0 0.2172
Figure 1: Tracking of state x2 by UKF, UGSF, URNDDR and OCUGSF for Case-I (xˆ0|0 = x0 )
5. Conclusions In this work, a new approach, labeled OCUGSF is proposed for constrained state estimation of nonlinear dynamical systems. The approach uses the attractive features of UGSF, in better representation of non-Gaussian densities using a sum of Gaussians, and explicitly incorporates the state constraints. Comparison with UKF, UGSF and URNDDR on a literature case study demonstrated the superior performance of the proposed OCUGSF approach.
References Anderson, B. D. O., Moore, J. B., 1979. Optimal filtering. Prentice-Hall Englewood Cliffs, N.J. Arulampalam, M. S., Maskell, S., Gordon, N., Clapp, T., 2002. A tutorial on particle filters for online nonlinear/nongaussian bayesian tracking. IEEE Transactions on Signal Processing 50 (2), 174–188. Gillijns, S., Mendoza, O., Chandrasekar, J., De Moor, B. L. R., Bernstein, D., Ridley, A., 2006. What is the ensemble kalman filter and how well does it work? American Control Conference (2006), 4448–4453. Julier, S., Uhlmann, J., March 2004. Unscented filtering and nonlinear estimation. Proceedings of the IEEE 92 (3), 401 – 422. Kalman, R., 1960. A new approach to linear filtering and prediction problems. Transactions of the ASME–Journal of Basic Engineering 82, 35–45. Kolas, S., Foss, B. A., Schei, T., 2009. Constrained nonlinear state estimation based on the ukf approach. Computers and Chemical Engineering 33 (8), 1386 – 1401. Kottakki, K. K., Bhartiya, S., Bhushan, M., 2014. State estimation of nonlinear dynamical systems using nonlinear update based unscented gaussian sum filter. Journal of Process Control 24 (9), 1425 – 1443. Narasimhan, S., Rengaswamy, R., 2009. Reply to comments on robust and reliable estimation via unscented recursive nonlinear dynamic data reconciliation (urnddr). Journal of Process Control 19 (4), 719 – 721. Sorenson, H., Alspach, D., 1971. Recursive bayesian estimation using gaussian sums. Automatica 7 (4), 465 – 479. Vachhani, P., Narasimhan, S., Rengaswamy, R., 2006. Robust and reliable estimation via unscented recursive nonlinear dynamic data reconciliation. Journal of Process Control 16 (10), 1075 – 1086. ˇ Simandl, M., Dun´ık, J., 2006. Sigma point gaussian sum filter design using square root unscented filters. Proceedings of the 16th IFAC World Congress 16 (1), 1000–1005.
Optimization Based Constrained Unscented Gaussian Sum Filter Krishna Kumar Kottakki, Sharad Bhartiya, Mani Bhushan Department of Chemical Engineering; Indian Institute of Technology Bombay; Mumbai, India [email protected], [email protected], [email protected]
Abstract Bayesian state estimation for constrained nonlinear systems involves two main issues: (i) appropriate representation of the underlying non-Gaussian densities, and (ii) use of appropriate techniques to bound the estimates so that they satisfy the give constraints. In this work we present a new constrained state estimation approach, named as Optimization based Constrained Unscented Gaussian Sum Filter. The proposed approach uses Unscented Gaussian Sum Filter to represent the non-Gaussian densities and an optimization approach to incorporate the state constraints. Benefits of the proposed approach are illustrated by applying it to a benchmark case study. Keywords: Unscented transformation, Unscented Gaussian sum filter, Sum of Gaussians
1. Introduction Recursive Bayesian estimation for nonlinear dynamical systems uses process and measurement models along with available process measurements to obtain conditional density of states. These process models are often derived from first principles, that is material and energy balances, and are typically nonlinear. Further, the process states inherently satisfy thermodynamic properties/constraints, for e.g. mole fraction of individual species/components in a multi component mixture are positive and their sum is unity. Hence, a practical nonlinear constrained state estimation approach should have the following two features: (i) an ability to approximate non-Gaussian densities arising out of nonlinear transformation or non-Gaussian noise, and (ii) an ability to incorporate constraints on the state estimates. Several nonlinear state estimation approaches such as Extended Kalman Filter (EKF) (Anderson and Moore, 1979), Unscented Kalman Filter (UKF) (Julier and Uhlmann, 2004), Ensemble Kalman Filter (EnKF) (Gillijns et al., 2006), Gaussian Sum Filters (Sorenson and Alspach, 1971) and Particle Filters (Arulampalam et al., 2002) have been reported in literature that account for the nonlinear transformations and/or the non-Gaussian nature of the densities. While EKF involves linearization of the nonlinear process and measurement models, other approaches utilize deterministically or randomly chosen samples thereby avoiding the explicit linearization. Among these sampling based approaches, UKF is based on unscented transformation (UT) that requires only (2n + 1) deterministically chosen sigma points, with n being the number of states. This use of limited number of samples opens up possibilities of applying UKF for solving practical (industrial) state estimation problems (Kolas et al., 2009; Vachhani et al., 2006). However, similar to Kalman Filter (Kalman, 1960), UKF also assumes that the underlying probability densities are Gaussian, an assumption commonly violated for nonlinear systems and/or non-Gaussian noise. In literature, a few modifications of UKF have been proposed that avoid the Gaussianity assumpˇ tion. In particular, Gaussian-sum UKF (GS-UKF) (Simandl and Dun´ık, 2006) is based on the
1716
Kottakki et al.
result that a sum of Gaussians can approximate any density to an arbitrary degree of accuracy (Sorenson and Alspach, 1971). GS-UKF uses N(2n+1) sigma points to obtain a sum of Gaussians based representation of non-Gaussian prior, where N refers to the number of Gaussians used and is a user defined parameter. Recently Unscented Gaussian Sum Filter (UGSF) (Kottakki et al., 2014) has been proposed that enables a sum of Gaussians representation of the non-Gaussian prior using only 2n+1 sigma points. In particular, the computational effort in UGSF is similar to that in UKF and it can thus be potentially applied to practical systems (Kottakki et al., 2014). The sum of Gaussians prior in UGSF is subsequently updated using Bayes’ rule. The updated Gaussian sum posterior in UGSF is reapproximated as a single Gaussian at the next time instant using the first two moments of the updated Gaussian sum density. This reapproximation avoids degeneracy of the individual Gaussians in the sum of Gaussians representation. Implementation on several case studies demonstrated superior performance of UGSF when compared to UKF. While UGSF addresses the issue of non-Gaussianity of the prior density, it does not incorporate constraints on the state estimates making them potentially infeasible and thus unusable. In this work, we propose an extension of the UGSF approach so that it can effectively incorporate state constraints in the estimation procedure. This is achieved by modifying the sigma point generation step to incorporate constraints as well as using explicit constrained optimization formulation to obtain the posterior density. The resulting approach is labeled as Optimization Based Constrained Unscented Gaussian Sum Filter (OCUGSF). In particular, the proposed OCUGSF solves 2n+1 optimization problems to obtain a constrained posterior estimate. Organization of this paper is as follows: Problem statement with a brief summary of elements used in the proposed method is presented in Section 2 and Section 3, respectively. Benefits of the proposed OCUGSF are demonstrated by comparing its performance with other state estimation techniques on a benchmark three state isothermal batch reactor case study in Section 4. The paper is concluded in Section 5.
2. Problem Statement Consider a sampled data system consisting of nonlinear process dynamics, a linear measurement function and interval constraints on the states as, tk
x(tk ) = x(tk−1 ) +
f (x(t), u(t))dt + wk , wk ∼ N (0, Q)
(1)
tk−1
yk = Hx(tk ) + vk , vk ∼ N (0, R) dk ≤ xk ≤ ek
(2) (3)
where, x(t) ∈ Rn , u(t) ∈ R p , represent the state and input vectors at time t while yk ∈ Rm , wk ∈ Rn , vk ∈ Rm represent observation, state noise and measurement noise, respectively at time tk . Function f : Rn × R p → Rn represents the nonlinear state dynamics and H ∈ Rmxn represents the linear observation model. Measurements yk are assumed to be available at regularly spaced sampling instants tk , at k = 0, 1, 2, 3, . . . with Ts = tk −tk−1 being the sampling interval. For ease of notation, we define xk x(tk ). Eq.(3) specifies the interval constraints on each component of state vector x. The problem is to find a point estimate for xk , whose dynamics are described in Eq.(1), using available measurements Yk = y1 , . . . , yk , in Eq.(2) subject to the interval constraints given in Eq.(3). We limit our attention to the linear measurement model Eq.(2) for ease of presentation. For nonlinear measurement models, the UGSF approach can be extended (Kottakki et al., 2014) using ideas presented in the current work. Now we briefly review the proposed OCUGSF for the given system (Eqs.(1), (2)) and the state constraints (Eq.(3)).
Optimization Based Constrained Unscented Gaussian Sum Filter
1717
3. Optimization Based Constrained Unscented Gaussian Sum Filter (OCUGSF) The proposed OCUGSF approach uses a three step framework to obtain constrained state estimate for the problem (Eqs.(1) to (3)). These steps are: (i) use of Interval Constrained Unscented Transformation (ICUT) approach (Vachhani et al., 2006) to generate the constrained sigma points from the given posterior density, (ii) propagation of sigma points through the process model, and (iii) a constrained optimization based update step. These are discussed next. 3.1. Interval Constrained Unscented Transformation (ICUT) Based Sigma Point Generation: Given that the random variable xk−1 has a Gaussian distribution with mean xˆk−1|k−1 and covariance Pk−1|k−1 at time tk−1 , the Interval Constrained Unscented Transformation approach deterministically selects 2n+1 constrained sigma points and corresponding weights as (Vachhani et al., 2006), ⎧ ⎪ ⎨ xˆk−1|k−1 , (i),c xˆk−1|k−1 + θi,k−1 [ Pk−1|k−1 ]i , χk−1|k−1 = ⎪ ⎩ xˆ k−1|k−1 − θi,k−1 [ Pk−1|k−1 ]i−n ,
w(i),c = bk−1 , i=0 w(i),c = ak−1 + θi,k−1 bk−1 , i = 1, . . . , n w(i),c = ak−1 + θi,k−1 bk−1 , i = n + 1, . . . , 2n (4)
ak−1 =
(2κ − 1)
2n √ 2(n + κ) ∑ θi,k−1 − (2n + 1) n + κ
√ 1 − ak−1 n + κ , bk−1 = 2(n + κ)
(5)
i=1
θi,k−1 = min (Θ j,i ), i = 1, 2, . . . , 2n 1≤ j≤n
where, for j = 1, 2, . . . , n; i = 1, 2, . . . , 2n ⎧ √ n+ if S j,i = 0 ⎪ ⎪ κ ⎪ ⎪ √ e j,k−1 −xˆ j,k−1|k−1 ⎨ n + κ, min if S j,i > 0 S( j,i) Θ j,i ⎪ ⎪ √ d −xˆ ⎪ ⎪ n + κ, j,k−1 S( j,i)j,k−1|k−1 if S j,i < 0 ⎩ min with S
Pk−1|k−1
− Pk−1|k−1
(6)
(7)
(8)
(i),c
In Eq.(4), χk−1|k−1 , w(i),c represent the ith constrained sigma point and its corresponding weight, respectively. It can be noted that these sigma points have been projected away from the mean using a factor θi,k−1 . Since this factor has been chosen according to the state constraints (Eqs.(6) to (8)), the sigma points will be feasible with respect to the interval constraints Eq.(3). Further, if the constraints are not active (dk = −∞ and ek = +∞), the expressions (Eqs.(4) to (8)) result in unconstrained sigma points obtained using UT as have been reported in the literature (Julier and Uhlmann, 2004; Vachhani et al., 2006). 3.2. Prediction Step: The constrained sigma points obtained in Eq.(4) are then transformed through the process model (Eq.(1)) to obtain the predicted sigma points as follows: (i),c
(i),c
χk|k−1 = χk−1|k−1 +
tk
f (x(t), u(t))dt, i = 0, . . . , 2n
(9)
tk−1
Without loss of generality, it is assumed that since the sigma points generated by Eq.(4) satisfy the (i),c interval constraints, the predicted sigma points χk|k−1 (Eq.(9)) will also satisfy these constraints. The proposed OCUGSF uses the UGSF design choices (Kottakki et al., 2014) to obtain nonGaussian prior as a sum of Gaussians. These design choices are: (i) Number of Gaussians in prior
1718
Kottakki et al. (i),c
is the same as the number of sigma points i.e 2n + 1, (ii) χk|k−1 is the mean of ith Gaussian in the sum of Gaussians prior, (iii) the covariance of each ith Gaussian, i = 0, . . . , 2n is identical and equal to the process noise covariance Q and (iv) weight w(i),c of the ith sigma point is the weight of ith Gaussian, i = 0, . . . , 2n in the sum of Gaussians. These design choices result in a constrained non-Gaussian prior: pxk |Yk−1 (ξk |Yk−1 ) =
2n
1
∑ w(i),c (2π)n/2 |Q|1/2 exp
#
i=0
$ −1 (i),c (i),c [ξk − χk|k−1 ]T Q−1 [ξk − χk|k−1 ] 2
(10)
Now we illustrate a procedure to obtain a constrained posterior corresponding to the prior (Eq.(10)), measurement model (Eq.(2)) and the state constraints (Eq.(3)). 3.3. Main Contribution: Constrained Update Using Optimization For the linear measurement model with Gaussian uncertainties (Eq.(2)), the likelihood density is given as, pyk |xk ,Yk−1 (ζk |ξk ,Yk−1 ) =
#
1 (2π)m/2 |R|1/2
exp
$ −1 [ζk − Hξk ]T R−1 [ζk − Hξk ] 2
(11)
where ζk ∈ Rm is the random variable corresponding to the measurement yk . Since the above likelihood density is an unconstrained Gaussian, direct application of Bayes’ rule with the nonGaussian prior (Eq.(10)) and likelihood densities (Eq.(11)), can yield an unconstrained sum of Gaussians posterior density (Vachhani et al., 2006). To avoid this situation, we propose to update the means of the individual Gaussians in the sum of Gaussians prior by using a constrained optimization formulation. This formulation replaces the conventional Bayes’ rule based update step (i),c as presented by Kottakki et al. (2014). Thus the constrained updated mean (χk|k ) and updated covariances (Ξ ) of ith individual Gaussians are given by, (i),c
χk|k = arg min
χk :dk ≤χk ≤ek
(i),c (i),c (yk − Hχk )T R−1 (yk − Hχk ) + (χk − χk|k−1 )T Q−1 (χk − χk|k−1 )
Ξ = Q − QH T [HQH T + R]−1 HQ
(12) (13)
(i)
while the updated weights δk are given by, (i),c
w˜ k
(i),c
δk
# $ 1 (i,c) (i,c) = w(i),c exp − (yk − Hχk|k−1 )T [HQH T + R]−1 (yk − Hχk|k−1 ) , i = 0, . . . , 2n 2 (i),c
= w˜ k
2n
/ ∑ w˜ k
( j),c
, i = 0, . . . , 2n
(14) (15)
j=0
Conventional Gaussian sum filters often encounter degeneracy in weights in subsequent iterations (Sorenson and Alspach, 1971). The proposed OCUGSF overcomes this problem by approximating the sum of Gaussians posterior density as a single Gaussian for subsequent step of sigma point generation at the next time instant (Kottakki et al., 2014). The mean and covariance of this single Gaussian are chosen to be the mean and covariance of the sum of Gaussians posterior density which are given as: c xˆk|k =
2n
∑
i=0
%
& (i),c (i),c χk|k
δk
c = , Pk|k
2n
%
∑ δk
i=0
(i),c
& (i),c
(i),c
c c T ][χk|k − xˆk|k ] Ξ + [χk|k − xˆk|k
(16)
To summarize, Eqs.(4) to (16) illustrate the steps involved in the proposed OCUGSF approach for obtaining the constrained posterior moments at kth time instant for the given posterior moments at (k − 1)th time instant.
1719
Optimization Based Constrained Unscented Gaussian Sum Filter
It can be shown that if the state constraints are not active, then the proposed OCUGSF approach reduces to UGSF approach. Further, if the non Gaussian prior density (Eq.(10)) is approximated by a single Gaussian instead of sum of Gaussians, then for the unconstrained and the constrained scenarios, OCUGSF reduces to the UKF (Julier and Uhlmann, 2004) and the Unscented Recursive Nonlinear Dynamic Data Reconciliation (URNDDR) (Vachhani et al., 2006; Narasimhan and Rengaswamy, 2009) approaches, respectively.
4. Results and Discussions A three state isothermal batch process is chosen to validate the performance of proposed OCUGSF. The process involves gas-phase reversible reactions under isothermal conditions, A B+C, 2B C. The dynamical process and measurement models are (Kolas et al., 2009): ⎡ ⎤ ⎡ ⎤ −kc1 x1 + kc2 x2 x3 x˙1 2 ⎣x˙2 ⎦ = ⎣kc1 x1 − kc2 x2 x3 − 2kc3 x + 2kc4 x3 ⎦ 2 x˙3 kc1 x1 − kc2 x2 x3 + kc3 x22 − kc4 x3 T yk = 32.84 32.84 32.84 x1 x2 x3
(17) (18)
Eqs.(17) and (18) represent the nonlinear process model and linear measurement model, respectively. In Eq.(17), x1 , x2 and x3 represent the concentrations of the species A, B and C, respectively, and kci (i = 1, . . . , 4) represents the rate constant for ith reaction. True states were generated for 100 different noise realizations consisting of 720 sampling time instants, i.e. for 180 seconds. The parameters used for true state generation and for state estimation are reported in Table 1. Table 1: Process and simulation parameters for isothermal batch reactor (Kolas et al., 2009) x0 = 0.5 0.05 0 Δt = 0.25seconds 1 =x xˆ0|0 0
kc = 0.5 0.05 0.2 0.01 Q = diag(0.0012 0.0012 0.0012 ) T dk = 0.0129 0.0743 0.0243
2 xˆ0|0
P = 4I3×3 0 T = 0.6 0.1 0.05
ek = 0.4758
0.3575
T
0.6601
The performance of estimation algorithms, namely, OCUGSF, URNDDR, UGSF, and UKF are 1 and investigated for two different scenarios of initial posterior densities, i.e Case I: xˆ0|0 = xˆ0|0 2 P0|0 = P0 , and Case II: xˆ0|0 = xˆ0|0 and P0|0 = P0 . Average of sum of the square roots of the estimation errors (ASSREE) was chosen %' to compare the performance& of estimation algorithms and is defined as: ASSREE =
1 100
100
∑
l=1
1 720
720
(l) (l) (l) (l) (l) (l) ∑ [xˆk|k − xk ]T [xˆk|k − xk ] . Here xk and xˆk|k represent k=1 l th simulation run, respectively. The ASSREE values for
the true and the estimated states for the four estimation algorithms along with their computational times are reported in Table 2 for both Cases I and II. The average number of constraint violations is also reported under the heading #violations in this Table. From the ASSREE values reported in Table 2, it can be observed that the proposed OCUGSF results in superior state estimates compared to the other three algorithms with higher, but acceptable computational cost. It can also be noted that OCUGSF estimates have satisfied the state constraints for all the time instants across the simulation runs. Even though the URNDDR approach has reported no constraint violation, the Gaussian approximation of nonGaussian prior (Narasimhan and Rengaswamy, 2009) has resulted in inferior performance compared to OCUGSF. For UGSF, while the prior approximation is non-Gaussian, it has violated the state constraints at several time instances as it is an unconstrained estimator. Gaussian approximation of non-Gaussian prior and unconstrained state estimation has resulted in inferior performance of UKF. State x2 tracking for one particular simulation run has been shown in Fig. 1 and illustrates the superior performance of OCUGSF over UKF, UGSF and URNDDR approaches.
1720
Kottakki et al.
Table 2: Performance comparison for three state isothermal batch reactor
Algorithm UKF UGSF URNDDR OCUGSF
Case-I (xˆ0|0 = x0 ) ASSREE #violations CPU Seconds 159.39 705 0.0271 5.9164 645 0.0273 3.5503 0 0.1131 2.1247 0 0.1400
ASSREE 155.29 6.0094 3.4396 2.7890
2 ) Case-II (xˆ0|0 = xˆ0|0 #violations CPU Seconds 705 0.0320 643 0.0322 0 0.1781 0 0.2172
Figure 1: Tracking of state x2 by UKF, UGSF, URNDDR and OCUGSF for Case-I (xˆ0|0 = x0 )
5. Conclusions In this work, a new approach, labeled OCUGSF is proposed for constrained state estimation of nonlinear dynamical systems. The approach uses the attractive features of UGSF, in better representation of non-Gaussian densities using a sum of Gaussians, and explicitly incorporates the state constraints. Comparison with UKF, UGSF and URNDDR on a literature case study demonstrated the superior performance of the proposed OCUGSF approach.
References Anderson, B. D. O., Moore, J. B., 1979. Optimal filtering. Prentice-Hall Englewood Cliffs, N.J. Arulampalam, M. S., Maskell, S., Gordon, N., Clapp, T., 2002. A tutorial on particle filters for online nonlinear/nongaussian bayesian tracking. IEEE Transactions on Signal Processing 50 (2), 174–188. Gillijns, S., Mendoza, O., Chandrasekar, J., De Moor, B. L. R., Bernstein, D., Ridley, A., 2006. What is the ensemble kalman filter and how well does it work? American Control Conference (2006), 4448–4453. Julier, S., Uhlmann, J., March 2004. Unscented filtering and nonlinear estimation. Proceedings of the IEEE 92 (3), 401 – 422. Kalman, R., 1960. A new approach to linear filtering and prediction problems. Transactions of the ASME–Journal of Basic Engineering 82, 35–45. Kolas, S., Foss, B. A., Schei, T., 2009. Constrained nonlinear state estimation based on the ukf approach. Computers and Chemical Engineering 33 (8), 1386 – 1401. Kottakki, K. K., Bhartiya, S., Bhushan, M., 2014. State estimation of nonlinear dynamical systems using nonlinear update based unscented gaussian sum filter. Journal of Process Control 24 (9), 1425 – 1443. Narasimhan, S., Rengaswamy, R., 2009. Reply to comments on robust and reliable estimation via unscented recursive nonlinear dynamic data reconciliation (urnddr). Journal of Process Control 19 (4), 719 – 721. Sorenson, H., Alspach, D., 1971. Recursive bayesian estimation using gaussian sums. Automatica 7 (4), 465 – 479. Vachhani, P., Narasimhan, S., Rengaswamy, R., 2006. Robust and reliable estimation via unscented recursive nonlinear dynamic data reconciliation. Journal of Process Control 16 (10), 1075 – 1086. ˇ Simandl, M., Dun´ık, J., 2006. Sigma point gaussian sum filter design using square root unscented filters. Proceedings of the 16th IFAC World Congress 16 (1), 1000–1005.