Integrating Covariance Intersection Into Bayesian Multitarget Tracking Filters

Multitarget tracking systems typically provide sets of estimated target states as their output. It is challenging to be able to integrate these outputs as inputs to other tracking systems to gain a better picture of the area under surveillance since they do not conform to the standard observation model. Moreover, in cyclic distributed systems, there may be common information between state estimates that would mean that fused estimates may become overconfident and corrupt the system. In this article, we develop a Bayesian multitarget estimator based on the covariance intersection algorithm for multitarget track-to-track data fusion. The approach is integrated into a multitarget tracking algorithm and demonstrated in simulations. The approach is able to account for missed tracks and false tracks produced by another tracking system.


I. INTRODUCTION
Target tracking systems typically use the Kalman filter [1] or some variant of it as the basis for online estimation. Measurements are received sequentially in time and the state estimate is updated according to Bayes' rule. If measurements are not available from a sensor and instead a mean and a covariance from an external tracking system are given, the Kalman filter is no longer applicable. Such systems are known as distributed target tracking systems, and such distributed data fusion algorithms combine the state estimates that are generated by a number of fusion centers or nodes. (Note that this is distinct from distributed systems that are not cyclic, e.g., [2].) A distributed tracking system has numerous advantages over a centralized system [3], [4] including robustness, scalability, and modularity. Such systems often require fusing the outputs from heterogeneous platforms with differing sensing and processing capabilities. Methods for fusing the outputs of different tracking systems have been developed for over three decades [5], [6]. One of the most successful approaches for distributed track-to-track fusion is the covariance intersection method of Julier and Uhlmann [7]. This is highly appealing due to its robustness, simple structure, and applicability to any tracking system that uses Gaussians as the basis for tracking. It has been widely applied to decentralized data fusion (DDF) problems [8] and on the Mars rover [9].
Generalizations of the covariance intersection algorithm to non-Gaussian systems have been proposed based on the exponential mixture density (weighted geometric mean) structure of the algorithm [10], [11]. Extensions for fusing two Gaussian mixtures have also been proposed [12]. Following the proposed multitarget generalization in [10], practical exponential mixture density implementations based on particular multitarget processes [13], [14] have been developed. These have been extended to sensor localization in distributed fusion networks [15] and distributed localization of sensors with partially overlapping field-of-views in fusion networks [16]. Recent applications of distributed multitarget tracking based on exponential mixture densities include consensus exponential mixture density [17], [18], human tracking [19], asynchronous radar system [20], orbit determination and space debris tracking [21], tracking based on multistatic Doppler shifts [22], and localization and tracking of mobile networks [23].
Recently, problems with cardinality estimation using exponential-mixture-density-based fusion of point processes have been highlighted [24]. The approach taken in this article is distinct from these exponential-mixturedensity-based approaches in distributed multitarget tracking in the following respect: It does not attempt to fuse two processes via the generalized version of covariance intersection. Instead, it develops a Bayesian fusion rule that can be integrated at the track level, which enables the exploitation of the robust and effective covariance intersection algorithm for the fusion of two multivariate Gaussian distributions.
This has a number of advantages over the exponential mixture density approaches. First, the approach can take sets of tracks from any type of multitarget tracking algorithm as the input. Second, it can be integrated into any multitarget tracking algorithm that uses Gaussians to represent the track state, either per-track or within a mixture, using the robust and fast covariance intersection algorithm approach. Third, the approach naturally deals with the fact that fields of view may be different via the detection probability inherited from the multitarget tracker. Fourth, it can deal with potentially false information through the inheritance of a clutter process for tracks.
The next section describes Bayesian disintegration and presents the central idea of interpreting the Covariance Intersection algorithm as a form of Bayes rule and calculates the marginal needed to integrate it into Bayesian target tracking filters. Section II-C describes how to use it in specific multitarget tracking filters. The approach is illustrated with a multitarget tracking algorithm. Section IV presents a study with a particular multitarget tracking algorithm.

A. Application of the Kalman Filter in Tracking Applications
In this section, we review the Bayesian interpretation of the Kalman filter as used in target tracking applications. The presentation is designed to facilitate the connection with the Bayesian interpretation of the Chernoff fusion rule and the covariance intersection algorithm presented in a novel context in the next section.
Conditional probability and Bayesian estimation relies on the following disintegration of a joint distribution p(x, z) into conditionals and marginals [25], [26] On the left, we consider p(x) to be the prior and p(z | x) to be a likelihood. The essence of Bayesian estimation is to start from the prior and determine a conditional update based on observations of the random variable z, i.e., p(x | z) = p(z|x)p(x) p(z) . In the calculation of the conditional update, we need to compute the marginal p(z), which can be determined by integrating over the joint, so that we can rewrite the decomposition as follows: We now consider the Bayesian disintegration for Gaussians, as used in the Kalman filter, i.e., given matrices H, R, and P and vector m of appropriate dimensions, the following identity holds: (3) where N (x; m, P) is a multivariate Gaussian in vectorvalued variable x with mean m and covariance matrix P, and the Kalman updated mean, state covariance, innovation covariance, and gain terms are given bỹ In the target tracking literature, the marginal p(z) = N (z; Hm, S) corresponds to the probability that the measurement z originated from the target. This is an important quantity in target tracking algorithms since it is used to determine weights of components. For instance, the following statements hold.

1) In Alspach and Sorenson's Gaussian sum filter, it is
used to calculate component weights in the Gaussian mixture posterior density (see [27, eq. (3.3)]). 2) In Reid's algorithm for multitarget tracking [28], this is used to determine the hypothesis weights [cf. (15) and (16)]. 3) In the probabilistic data association (PDA) algorithm of Bar-Shalom, it is used to determine the likelihood ratio that a measurement belongs to a target rather than clutter (cf. [29, eq. (38)]). 4) In the Gaussian mixture Probability Hypothesis Density filter of Vo and Ma [30], it is used to determine the weights of different intensity components [see (41)].
The exponent of the Gaussian is also used to determine a validation region for gating, i.e., [29] where the gate threshold γ corresponds to a gate probability that the validation region contains the true measurement if detected.

B. Chernoff Fusion as Bayesian Estimation
We now apply this classical approach to the Chernoff fusion rule for Gaussians, i.e., the covariance intersection algorithm. To do this, we reinterpret the Chernoff fusion rule as a Bayesian update rule and determine the analogous quantities to use the result as a tracking filter update. The covariance intersection algorithm provides the updated mean and covariance, and thus, we only require to determine the equivalent marginal after the update.
In distributed data fusion, it is common to use a fusion rule [10], [11] based on Chernoff information [31] to determine a density p ω (x) based on two densities describing the same random variable, p(x) and q(x) and a mixing parameter ω, given by Such descriptions are also known as escort distributions in the statistical physics literature [32] and have been used to determine information-related quantities, e.g., [33]. When p(x) and q(x) are multivariate Gaussian distributions, this is equal to the covariance intersection algorithm [7]. The structure of the Chernoff update rule is similar to the Bayes update rule if we permit ourselves to consider the term ( q(x) p(x) ) ω as a kind of unnormalized likelihood. We start with the prior and determine the fused posterior p ω (x) with the relations We note, of course, the nonstandard nature of the proposed (unnormalized) likelihood in that we have not specified a random variable, and that it depends on a whole posterior as input q(x), as well as the prior p(x). This clearly does not factorize into the familiar likelihood-prior and updatemarginal relations due to the dependence on q(x) and p(x). We also note that the likelihood and proposed marginal are both in unnormalized forms, so that we have not explicitly defined the random variable nor, by extension, the joint distribution. However, in Bayesian disintegration [26], p. 325], such convenient decompositions do not always exist. We can see that in this case, if we consider p(x) as the prior and p ω (x) as the update, then we do not have the classical Bayesian decomposition.
The objective of the current work is to determine a Bayesian disintegration that would permit us to use the expression q(x) p(x) ω as a likelihood and find the related marginal. If this is possible, then we can integrate the approach into mixture-based filters. Due to nice properties of Gaussians, and their utility for target tracking applications, we shall concentrate on Gaussian descriptions of q(x) and p(x).
In this case, we have the nice property that all of the terms are in the form of multivariate Gaussians, which enables us to determine analytic expressions. Let us now consider a disintegration for the Chernoff update rule with two Gaussians N (x; a, A), N (x; b, B) with the same dimensions for the means a, b and covariances A, B, i.e., (8) where d and D are the updated mean and covariance determined by the covariance intersection algorithm [7] If we can find a suitable normalization for the marginal on the right, then we have specified a joint distribution and a prospective Bayesian disintegration. Note that the unnormalized marginal can be determined from the Chernoff information [34], so that we can write it explicitly as which is equal to where we define C ω (A, B) to be If we consider the marginal to be a probability density in a, then this is a scaled Gaussian distribution in a, and we can write where This description of the joint distribution is the key result in the article that enables us to integrate covariance intersection into target tracking filters. More specifically, the specification of the marginal N (a; b, V ) permits us to integrate the covariance intersection algorithm into Bayesian tracking filters since it enables us to determine weight calculations determined via Bayes' rule. See Appendix A for an example of how to employ the result in the Gaussian sum filter of Alspach and Sorenson [27]. While the update is the same as the covariance intersection algorithm for determining the updated mean and covariance, the calculation of the marginal enables us to calculate weights in target tracking filters in a Bayesian manner in the same way as described in the previous section. We note here that here we have taken the marginal to be the probability density function in vector a and consider the covariance A to be a parameter for the likelihood, i.e., we do not have a probability distribution over mean and covariance pairs (a, A). The appropriateness of this choice of likelihood function will be assessed through statistical simulations in the following sections.
This line of reasoning has enabled us to recast the Chernoff update rule as a Bayesian update rule. With this new interpretation, we can take the mean a and covariance A of another target track estimate and use them to determine an updated distribution. In the Gaussian case, this is exactly the covariance intersection rule. The exponent of the Gaussian is also used to determine a validation region for gating, i.e., The covariance intersection algorithm provides the updated mean d and covariance D in the updated Gaussian N (x; d, D).

C. Application to Target Tracking Algorithms
We now describe how to apply the approach in different tracking filters by replacing the Kalman fusion rule with the new rule formed by the new joint distribution with reference to particular algorithms. In a multitarget context, we assume that we are given a set of thresholded tracks from another tracking system. By analogy with classical multitarget tracking models, we assume that we have a detection probability that gives the probability that a target under surveillance gives rise to a track and that a clutter model describes the distribution of false tracks. We note that the key result itself is straightforward and can be applied in multitarget tracking systems of a different nature since they use the same basic Gaussian components. Although this relies on known results for the multiplication of Gaussians, the application of the covariance intersection method has not been applied in this way before. In the next section, we illustrate the use of the method with a multitarget tracking algorithm in a simulated study.

A. Gaussian Sum Filter
In this section, we look at different update mechanisms for Gaussian mixture priors and representations of the joint distributions. We show that it is possible to determine different conditional distributions depending on the choice of joint distribution. Suppose that w (i) for i = 1, . . . , N are nonnegative and sum to 1. Suppose further that p (i) (x) for i = 1, . . . , N are probability density functions. Then, the usual approach for updating mixture distributions with a likelihood can be described with the following Bayesian disintegration: For instance, the Gaussian sum filter of Alspach and Sorenson [27] is determined through the relation Each integral on the right can be calculated with the Kalman filter relation.
As we have seen in the previous section, the likelihood that we are interested in is inherently coupled to the prior. Hence, we need to carefully consider how to apply the new Bayesian approach for mixture distributions. However, we note that though mixture distributions are conventionally used in the form above, there is no special reasoning needed to define joint mixture distributions in the following way: In this formulation, we consider N hypotheses, where each hypothesis has its own likelihood. For instance, in the Gaussian case, we can consider For intuition purposes, we could think of this as a scenario where there is uncertainty about which sensor provided observation z, represented with the mixture weights w (i) , and where each sensor has its own observation matrix H (i) and noise covariance R (i) , described through likelihood N (z; H (i) x, R (i) ). We will use this description as a motivation for finding an analogous result for the Chernoff rule in the following section.

B. Chernoff Rule With Mixture Priors
We now consider relations for mixtures based on different applications of the Chernoff fusion rule. We shall be particularly interested, in the following section, in fusing a Gaussian distribution with a Gaussian mixture. The usual approach for updating mixture distributions with a the Chernoff rule can be described with the following relation: (20) and in the case where q(x) is a Gaussian (input track) and p(x) is the prior, we have This is problematic since raising the Gaussian mixture to a power does not have a nice analytic form. Julier proposed the following approximation [37] for this scenario: and This has been used in some applications [17], [18]. Based on the approach described in the previous section, we propose a different approximation based on the joint distribution for Gaussians given in the previous section, i.e., we consider the relation (25) or, using the relations in the previous section (26) where This leads to the following description for p ω (x): where the updated weights are determined witĥ This does not take the original form of the Chernoff fusion rule. However, it has the advantage that it does not have the complication of raising a Gaussian mixture to a power. As we have seen, this takes the form of a conditional probability update. However, since the different modeling approach does not impose the same likelihood for each term in the mixture, we are able to achieve a closed-form solution in terms of Gaussian mixtures. In the next section, we use these results in the context of Bayesian estimation for point processes.

IV. EXPERIMENTS
In this section, we shall analyze the performance of the new covariance intersection update method through two scenarios using simulated data. To avoid the computational complexities when calculating the ω weighting fusion parameter, a fast intersection method is used to compute the exponential weighting parameter [38]. Each of the multitarget filters used in the experiments below (with the exception of Scenario A) are LCC [39] filters, following an optimized Gaussian mixture multitarget filter implementation as described in [36]. The Mahalanobis distance is used to evaluate potential similar components in the merging step with threshold τ merge .
Each scenario shares the same underlying target dynamic and sensor model with shared parameters as seen in Table I. It should be noted that some of these parameters are only applicable for the multitarget scenarios. The target space is over a 2-D region of dimensions X × Y (m). M targets are initially generated uniformly across the region at time k = 1 and are detected with a probability P d with survival probability P s . New targets are birthed every timestep with a Poisson process with mean μ γ . These targets follow a near-constant velocity motion model with state dimensions [x,ẋ, y,ẏ], transition matrix F k , process noise Q k , and process noise standard deviation σ q . Targets are observed by the sensor with an observation matrix H k , noise R k and measurement noise standard deviation σ r . Each simulation is run for T timesteps with uniformly spatially distributed false alarms being generated at each timestep according to a Poisson process with mean μ λ . All results shown are averaged over 100 Monte Carlo (MC) runs of each scenario. The Optimal SubPattern Assignment (OSPA) metric is used to assess the performance [40]. In the following scenarios, an OSPA cutoff c parameter of 100 and an order p-value of 2 are used. Only targets that possess a weight above an estimation threshold τ track are used in calculating the OSPA. The lines in the graphs represent the averaged results over the MC runs. For the execution time per timestep and OSPA results, the shaded areas show the 2σ standard deviation of the runs. For the cardinality results, the shaded areas show the mean of the 2σ standard deviation of the multiobject filter estimate. Each of the simulations has been running on a PC with a 12-core AMD Ryzen 3900X CPU with 64 GB of RAM.

A. Scenario A
Scenario A integrates the new rule into a Gaussian sum filter [27] and compares with the centralized Gaussian sum filter and the previously proposed approach by Julier [12]. Since we use linear-Gaussian dynamic and observation models, the filters typically converge very quickly. Hence, rather than compare their performance over time, we analyze their performance at the first time step, i.e., we focus on the accuracy of the Bayesian estimation step. In this scenario, we consider the performance of the update for varying measurement noise of the filter feeding into the update. Ten values of σ r , the standard deviation of the measurement noise, were logarithmically spaced between 0.01 and 1.
At time k = 1, a single target is generated and for each Gaussian sum filter, N = 1000 Gaussian components are uniformly distributed across the state space. No Gaussian mixture reduction techniques are performed and every component is updated with each obtained measurement.
The following four scenarios are considered with the results presented in Fig. 1. 1) A Gaussian sum filter with a measurement update from one sensor (solid blue line). 2) A "centralized" Gaussian sum filter with a measurement update from two different sensors (dotted blue line). 3) A Gaussian sum filter with one measurement update that is then updated with Julier's method with a Kalman filtered mean and covariance from another sensor (red line). 4) A Gaussian sum filter with one measurement update that is then updated with the new method with a Kalman filtered mean and covariance from another sensor (green line).
The weighted average estimate is then extracted from each filter and compared in terms of root-mean-square error (RMSE) from the ground truth target in Fig. 1. The results indicate that the new rule provides more accurate estimates than the previously proposed approach, i.e., that the calculation of new mixture weights enables more accurate estimates of the filter, as in the previous approach distant components would be given the same weight are closer components. It should also be noted that while the two covariance intersection methods perform significantly better than the single Gaussian filter, they do perform worse than the centralized filter for each measurement noise level.

B. Scenario B
Scenario B simulates a challenging scenario in order to test the new fusion method. Here, the simulated scenario has been extended to a multitarget case following the parameters and structure detailed above, and where the probability of detection P d = 0.7, the probability of survival P s = 0.99, the initial number of targets is set at 100 with clutter rate μ λ = 200.
The results indicate that the cardinality estimates [see Fig. 2(a)] of the new approach enable more accurate calculations than a single Gaussian mixture filter, though it is not as accurate as the centralized approach, as indicated in the OSPA calculations in Fig. 2(b). The results for the cardinality estimate indicate that in the scenario considered, there is not a large difference between the new approach and the centralized approach, though the centralized approach does outperform the new approach.
C. Scenario C Scenario C analyzes the effect of varying the probability of detection on the fusion methods. Initially M = 10 targets are generated with a high probability of survival P s = 0.999  and zero clutter rate μ λ = 0. Starting from a low probability of detection P d = 0.5, 25 values are linearly spaced up to P d = 0.9. For each value of P d , the scenario is run for T timesteps. The OSPA metric for the multitarget filter is calculated at each timestep and then is averaged over the T timesteps for each P d , as displayed in Fig. 3. The following three multitarget filter variations are considered.
1) An LCC filter is updated with a single set of measurements (solid red line). 2) A "centralized" LCC filter, which is updated twice at each timestep with two separate sets of measurements (dashed blue line). 3) An LCC filter is updated with a single set of measurements and then updated using the new rule and The results in Fig. 3 show that as to be expected, both the centralized and new fusion approach multitarget filters perform noticeably better than the single multitarget filter and that there is no significant difference in performance between the centralized filter and the new approach.

D. Scenario D
Scenario D follows the same simulation structure as the previous scenario, except now the clutter rate μ λ is varied instead of the probability of detection. Still for M = 10 targets, μ λ is gradually increased from μ λ = 1 to μ λ = 20 M = 200 over 50 linearly spaced values. Similarly to Scenario B, the OSPA metrics shown in Fig. 4 show the averaged OSPA metric over the T timesteps for each value of μ λ .
The results in Fig. 4 indicate that for low clutter, the centralized approach performs best. However, as the clutter level increases, the fusion approach can give more accurate estimates than the centralized approach. This is likely due to the fact that the new approach is using data that have already filtered out the clutter and has a better estimate of the number of targets.
We would typically expect that a centralized tracker should outperform a distributed tracker. However, we need to consider what is being fused in order to draw such a conclusion.
1) In the centralized tracker, we take the output of two sensors and estimate the target population based on applying Bayes' rule twice, where the observations are sets of measurements from the two different sensors (including false alarms and missed detections). 2) In the distributed scenario, the targets are updated with the measurements from one sensor with Bayes' rule, and then again with the new fusion rule (also based on Bayes' rule at the target number level) with tracks that have been estimated over time with another tracker.
We can reasonably expect that this other tracker has a better estimate the number of targets better than a full measurement set from another sensor, which is true. Hence, it is actually not surprising that the estimate in the number of targets can get better when we apply the new fusion approach, which is also based on the Bayes rule. If we feed better data, we can expect a better outcome: The data from the other tracker have already filtered out the false alarms. The authors believe that this is why there is a better estimate of the number of targets in this scenario as indicated by the experiments.

E. Scenario E
Finally, the effect of the data incest problem is investigated. In a DDF context, data incest is when the same information is used multiple times along a data processing chain. A covariance intersection approach [7] has inherent protection against this through the nature of its fusion; however, the Gaussian mixture approximations of the multitarget filters used above have no such protection since the approach is Bayesian. So in order to test the robustness of the proposed strategy, the following scenario is used. Two LCC filters are run simultaneously and are colocated both spatially and temporally. At each timestep, a full tracking cycle is executed: a local prediction and update followed by broadcasting its local estimate and then an update with neighbor's estimates. The second half of this cycle, broadcast and update, is then repeated an R CI number of times. In the following results (see Fig. 5), 1, 2, and 5 repeats are used.
This failure of the method in this scenario is confirmed by the cardinality and OSPA results, Fig. 5(a) and (b), where the 5 repeats run, shows a severe degradation in performance when compared to the other runs. This performance degradation is also shown in the 2 repeat run where it consistently overestimates the number of targets present. This overconfidence can be explained by the results shown in Fig. 5(c), whereas more updates are performed, the underlying multitarget filters increase their confidence in the existing targets by increasing the weights and decreasing the Kalman covariances of the corresponding Gaussian components. However, due to the repeated nature of the updates and the lack of new information these Kalman covariances decrease to a point, where the gating operations fail and the input tracks are no longer associated correctly. This leads to the R CI = 5 results, where by this time, the weights of the components have fallen so significantly that the components have pruned from the mixture. These results clearly show that the Gaussian mixture implementations cannot adequately deal with the data incest problem alone and, therefore, should not be used in cyclical networks where there is a large sharing of data between nodes. Instead, if the network topology is that of a directed acyclic graph, this approach should prove useful.

V. CONCLUSION
A new fusion rule for distributed multitarget tracking is proposed based on Chernoff information and covariance intersection. The method has been developed to exploit sets of tracks produced by other tracking systems. The approach has the advantage that it can be integrated directly into a range of existing multitarget tracking systems through the replacement of Kalman filtering updates with covariance intersection updates with a fast fusion weighting calculation. The approach takes advantage of the multitarget tracking system model to account for missed tracks and false tracks produced by the other tracking system, hence providing a robust approach for distributed sensor fusion. It also has the advantage of being able to integrate tracks from multitarget tracking systems of a different nature. Since the approach is Bayesian, there still remain potential problems of data incest in cyclic networks. However, when it can be reasonably assumed that there is low communication, then the approach offers a rapid way of integrating tracks produced by another system.