Uncertainty Exchange Through Multiple Quadrature Kalman Filtering

One of the major challenges in Bayesian filtering is the curse of dimensionality. The quadrature Kalman filter (QKF) is the method of choice in many real-life Gaussian problems, but its computational complexity increases exponentially with the dimension of the state. As a promising solution to overcome the filter limitations in such scenarios, we further explore the multiple state-partitioning approach, which considers the partition of the original space into several subspaces, with the goal to apply a low-dimensional filter at each partition. In this contribution, the key idea is to take advantage of the estimation uncertainty provided by the QKF to improve the interaction among filters and avoid the point estimate approximation performed in the original Multiple QKF (MQKF). The new filter formulation, named Improved MQKF, considers Gauss-Hermite quadrature rules to propagate the subspaces of interest, together with cubature rules for marginalization purposes. The nested quadrature-cubature approximation provides robustness and improves the filter performance. Simulation results for a multiple target tracking scenario are provided to support the discussion.

and n k are assumed to be independent. The solution to the Bayesian filtering problem is given by the marginal distribution p(x k |y 1:k ), which gathers all the information about the states at time k given by the available measurements, y 1:k = (y 1 , . . . , y k ). A plethora of alternatives are available in the literature, being particle filters (PFs) [2] the most popular in general nonlinear/non-Gaussian systems, and the quadrature Kalman filter (QKF) [3]- [5], under the Gaussian assumption of interest here. One of the major challenges in the development of new filtering methods is the well-known curse of dimensionality [6], that is, the exponential computational complexity increase in high-dimensional systems [7], which implies a loss of accuracy and possible divergence.
To circumvent the curse of dimensionality and associated problems, one may consider several strategies: i) Rao-Blackwellization to reduce the state-space dimension to be explored by the nonlinear filter [8], which can be used in systems with linear substructures; ii) posterior independence assumption PFs [9]- [11]; and iii) the multiple approach based on state partitioning introduced in [12]- [15], leading to the so-called multiple PF (MPF) and multiple QKF (MQKF). In this letter, we propose an improvement of the original MQKF.

A. Standard MQKF Background
In the system of interest, we assume that the state in (1) can be partitioned into S subspaces as x k = [x (1) k , . . . , x (S ) k ], such that Q k −1 = diag(Q (1) k −1 , . . . , Q (S ) k −1 ). Under this assumption, we can write the state equation as 1 where functions f (s) k −1 (·) can be different, and the sth process noise is independent and distributed as ν The key idea behind the standard MQKF [15] is to apply an individual filter to each subspace in (3), directly reducing the state dimension each filter must deal with and its computational complexity. Mathematically, this implies that the sth filter is in charge of approximating the subspace marginal predictive and posterior distributions p(x (s) k |y 1:k −1 ) and p(x (s) k |y 1:k ), respectively. To compute such distributions, the sth filter assumes: 1) the joint predictive and posterior distributions, may be 1 x (s ) denotes the sth element (possibly a vector) in a vector x and x (−s ) is the vector of all elements in x except for x (s ) . The dimension of each subspace n (s ) written as, and 2) as done in [12] and [13], the interconnection among subspaces is approximated using point estimates, which are available from the other filters running in parallel, or equivalently, x Under these assumptions, the prediction and update steps of the filter in charge of the sth subspace reduce to approximating Gaussian integrals of the form 2 f (x (s) The MQKF uses Gauss-Hermite quadrature rules (i.e., a set of deterministic sigma-points and weights) [3], [4] to approximate such integrals; therefore, it is a bank of S parallel QKFs, which exchange information at each step.

B. Contribution
The point estimate approximation in (6) and (7) is not a robust solution and may lead to divergence of the filter, because it does not take into account the uncertainty on the predicted/estimated states (i.e., quality of the estimates), which is characterized by the error covariance. Moreover, as the subspace marginal predictive and posterior distributions are Gaussian and fully characterized by their mean and covariance, which are available from the filters running in parallel, we can improve such approximation.
The main goal of this contribution is to further explore and improve the original MQKF in [15]. First of all, we drop the assumption in (6) and (7), and rewrite the general formulation of the algorithm. Then, to obtain a more accurate subspace marginalization, we propose to use a new nested quadraturecubature sigma-point approximation of the integrals in the general solution. The new approach, named improved MQKF (IMQKF), considers Gauss-Hermite quadrature rules [4] to propagate the partition of interest and cubature rules [16] to marginalize the rest of subspaces. The coupled filters only need to exchange the predicted and estimated mean and covariance at each prediction and update steps. We prove that under reasonable conditions, the complexity is lower than that of QKF and, although larger than MQKF, the new method allows inference over large-dimensional systems, where both QKF and MQKF cannot effectively operate. The improved methodology is applied to a representative multiple target tracking scenario.

A. General Filter Formulation
The state-partioning strategy considers a bank of S parallel filters, each one tracking a single subspace. Considering the filter These distributions are recursively computed in two steps: prediction and update. The general formulation of the marginal predictive distribution in this case reads In the prediction step, the filter is interested in the mean and corresponding prediction error covariance of such distribution In the subsequent update step, the filter approximates the mean and covariance of the marginal subspace posterior p(x (s) k |y 1:k ), given by the Kalman update equations [17] x (s) (14) where the Kalman gain is K To compute the updated estimates, the filter needs the predicted measurement and both innovation and cross covariance matrices, which are given by the following integrals:

B. Nested Quadrature-Cubature Sigma-Point Approximation
In general, the problem reduces to approximating the integrals in (11), (12), and (15)- (17). The original MQKF solution in (6) and (7) is the simplest choice, but as already stated, a much better approach is to use the full subspace marginal distributions, i.e., p(x (−s) k −1 |y 1:k −1 ) and p(x (−s) k |y 1:k −1 ). In this contribution, we propose a new formulation to avoid the approximation in [15], exploiting the distributions recursively characterized by each parallel filter. These filtered densities gather knowledge about the uncertainty of the various tracked subspaces, which is now taken into account. Notice that the dimension of the concatenation of the −s subspaces is, potentially, much greater than the dimension of the subspace of interest, n (−s) x n (s) x . Therefore, considering the curse of dimensionality and computational complexity reduction problems, it may not be a good idea to use Gauss-Hermite quadrature rules [4] to marginalize such subspaces. In contrast with the exponential increase in the number of evaluation points in the Gauss-Hermite quadrature rule with the state dimension, the number of sigma points of the cubature rule [16] grows linearly with the state dimension. Therefore, in this contribution, we propose to use a cubature-based marginalization of the subspaces x Notice that at time k, the following distributions and x,k −1|k −1 ) are available from the previous time step.x (−s) k −1|k −1 is build from the concatenation of the different subspace estimates, except for the s-th, and the corresponding covariance matrix Σ (−s) x,k −1|k −1 , which characterizes the filter prediction uncertainty, is build from the individual covariances as a block diagonal matrix. After the prediction step, the filter in charge of the s-th subspace is aware of the information provided by the rest of the filters, thus the distributions p(x x,k |k −1 constructed as before. The characterization of these Gaussian subspace marginal distributions can be used to correctly approximate the integrals in (11), (12), and (15)- (17). The key idea is to use a nested quadrature-cubature sigma-point formulation. In the prediction step, two sets of sigma-points are used for x  (13) and (14) can be approximated using this sigma-points sets aŝ (19) where j implements the marginalization over x (−s) k −1 . Update: As done in the prediction stage, sigma-point rules can be used to approximate the integrals in (15)- (17). First, we compute the transformed sets x,k |k −1 ) . Then, we approximate the integrals of interest using quadrature rules aŝ where j implements the marginalization over x (−s) k . To summarize, considering α points per dimension for the Gauss-Hermite quadrature points, both sets are It is worth saying that it may be useful for numerical stability to implement the square-root version of each QKF [5]. Finally, it is straightforward to use sparse-grid quadrature rules to further reduce the computational complexity with negligible penalty in numerical accuracy [7], [18].

C. Complexity Analysis
In analyzing the complexity of sigma-point methods [7], the number of sigma-points plays the most important role as each sigma-point involves an evaluation of a, potentially, nonlinear function. Let us denote L QKF = α n x and L MQKF = s α n ( s ) x the number of points used by a QKF and a MQKF, respectively, it was proved in [15, Proposition 1] that L MQKF < L QKF when S ≥ 2 and α > 2.
For the IMQKF, one can show that L MQKF < L IMQKF ≤ L QKF which implies that the computational complexity of the IMQKF is between the MQKF and QKF. The first inequality is straightforwardly obtained observing that n x − n (s) x > 0, and thus, L IMQKF = s α n ( s ) The proof for the second inequality, i.e., L IMQKF ≤ L QKF , can be consulted in the supplementary material to this letter.

III. MULTIPLE TARGET TRACKING
Illustrative numerical results of the new method's performance are shown in a multiple target tracking example, where K targets move in a two-dimensional plane and are tracked using a set of N = 100 sensors, uniformly distributed in a deterministic 100 × 100 m 2 grid, receiving signal strength. The states to be tracked for target i are position and velocity, x where the process noise models the possible system errors due to acceleration perturbations, with Q = G ν diag(σ 2 a x , σ 2 a y )G ν T s = 1 s being the sampling period. Notice that MQKF and IMQKF exploit the independence among targets to partition the state space, being each subspace defined as in (30). At time k, the mth sensor signal strength model is given as [14] y m ,k = K i=1 10 log 10 1 x,k , p with P (i) 0 = diag(10, 10, 0.1, 0.1). The performance of the IMQKF, which is compared with the original MQKF [15] and both standard QKF [3], [4] and Cubature KF (CKF) [16], is shown in Fig. 1. Table I summarizes the number of sigma-points at each filter, for α = 3. Additionally, a PF is used to solve the same problem, particularly a sampling importance resampling PF [19] with transitional prior as importance distribution. Two PFs are used for comparison: one with the same number of particles as the QKF and another with much larger number of samples to ensure convergence. We have also implemented the quasi-Monte Carlo Kalman filter (QMCKF) in [20] using a Halton sequence with the same number of points as in the IMQKF. For this number of points, the QMCKF does not finish the simulations due to negative-definite covariance matrices, implying that the IMQKF makes a more efficient use of this number of samples. Fig. 1 (bottom) shows that the QKF-based solutions are equivalent for K = 1 and 2 targets (n x = 4 and 8), and always better than the CKF. The reason why there are no results for the QKF at K ≥ 4 (n x ≥ 16) is clear from Table I, the number of samples is too large. Similarly, PF approaches cannot be efficiently implemented when the state dimension increases. A number of 5000 particles (K = 1) and 50 000 particles (K = 2) are needed for stable operation. For K ≥ 4, the number becomes prohibitive. Even more, when the number of particles is that of QKF, the PF is likely to diverge, as stated in the figure with the percentage of diverged tracks. Fig. 1 (top) and (middle) shows the different KF-like filters' performances in time for K = 4 and 6 targets (n x = 16 and 24). The original MQKF [15] is able to cope with four targets, but fails if the number of targets increases.

IV. CONCLUSION
Multiple QKF is a method that provides accurate tracking solutions in high-dimensional dynamical systems, where the Gaussian assumption holds. The state space is partitioned such that a bank of QKFs is employed, one individual filter per partition, requiring some information exchange. In the original work [15], this interaction was in the form of point estimates. In this letter, a new approach to marginalize out the contribution of the other partitions is provided, accounting for the uncertainty of the individual point estimates and properly fusing them in a Bayesian sense. Particularly, while Gauss-Hermite rules are used at each filter to propagate the partition under interest, cubature rules are used for marginalization purposes due to their lower computational requirements. At the expense of a reasonable complexity increase with respect to [15], the new filter can perform notably in large dimensional systems where QKF cannot due to the curse of dimensionality.