A Novel Structure for a Multi-Bernoulli Filter without a Cardinality Bias

: The original multi-target multi-Bernoulli (MeMBer) filter for multi-target tracking (MTT) is shown analytically to have a significant bias in its cardinality estimation. A novel cardinality balance multi-Bernoulli (CBMeMBer) filter reduces the cardinality bias by calculating the exact cardinality of the posterior probability generating functional (PGFl) without the second assumption of the original MeMBer filter. However, the CBMeMBer filter can only have a good performance under a high detection probability, and retains the first assumption of the MeMBer filter, which requires measurements that are well separated in the surveillance region. An improved MeMBer filter proposed by Baser et al. alleviates the cardinality bias by modifying the legacy tracks. Although the cardinality is balanced, the improved algorithm employs a low clutter density approximation. In this paper, we propose a novel structure for a multi-Bernoulli filter without a cardinality bias, termed as a novel multi-Bernoulli (N-MB) filter. We remove the approximations employed in the original MeMBer filter, and consequently, the N-MB filter performs well in a high clutter intensity and low signal-to-noise environment. Numerical simulations highlight the improved tracking performance of the proposed filter.


Introduction
Considering the multi-target tracking (MTT) environments, the unknown number of targets changes with time because of the presence of target deaths and births. Moreover, the ambiguities in track-to-measurement association, which involves the missed detection and clutters, increase the difficulty of the joint estimation of the target cardinality and target states [1][2][3]. Based on the finite set statistics (FISST) theory, Mahler proposed a rigorous formulation of random finite set (RFS)-type filters. The RFS-based algorithms model the targets and the measurements as RFSs [3]. Compared with the conventional multiple hypothesis tracker (MHT) and joint probability data association (JPDA) algorithms [4][5][6], they avoid data association, which has a high computation burden. In the past few years, several RFS-based methods, such as the well-known probability hypothesis density (PHD) filter [7][8][9], cardinality PHD (CPHD) filter [10][11][12], and MeMBer filter [3] have been created. For example, the PHD filter models both targets and measurements as Poisson RFSs, and sequentially propagates the first moment of the multi-target density (PHD intensity) over time. The CPHD filter improves the cardinality-estimating performance of the PHD filter by propagating the cardinality distribution with the PHD of the multi-target density. The MeMBer filter parameterizes the multitarget distribution that models each potential target as a single Bernoulli RFS [3], which is characterized by the probability of existence and the probability density function (pdf). Accordingly, the multi-Bernoulli (MB) RFS is a union of multiple independent Bernoulli components. Since the Gaussian mixture implementations are proposed for the RFS-based filters, various extensions are where ( ) ( ) Throughout this paper, an MB distribution is abbreviated as: Given a Dirac delta density x δ , we define the functional derivative of F h and given as: The first-order moment, i.e., the PHD function ( ) D ⋅ of a multi-target density ( ) π ⋅ , is given as: The PHD function of an MB is given as [3]: A key result of the functional derivative utilized in our derivation process is the product rule [3]: where  indicates the sum is over all disjoint sets 0 ,..., n W W and 0 , ...
, which permits the calculation of the derivative of a product of an MB PGFl.

Novel Structure of an MB Filter
In this section, the derivation process of the proposed algorithm is presented. Before we commence, the following assumptions are presented, which are used in the derivation of the algorithm.

Assumptions:
a. Birth target distribution is denoted by b. The survival probability conditioned on the target state is

Derivation of the Predictor
At each time 1 k − , the PGFl is also a multi-target multi-Bernoulli process [3]: Here, b k G h     is the PGFl of the birth process, which follows Assumption 1: Also, the predicted parameters of the surviving targets are given as:

Derivation of the Corrector
For notational simplicity, in the following derivation, we omit the time index. As given in the above subsection, the PGFl of the predicted distribution is the MB process: where Denoting the multi-target measurement likelihood function conditioned on the target state set X as The PGFl of the updated density ( ) | k k X π (updated by k Z ) is given as [3]: Using the product rule given in Equation (11) where The set of all assignment vectors θ is denoted as Θ .
The functional derivative of 0 F in 0 g = WRT the set 0 W is given as: where V denotes the "volume" of the surveillance region. For 0 n i ≥ > , the functional derivative of is given as: Furthermore, set 0 g = : Using the association notations and the results of Equations (24) and (26), we rewrite Equation (20) as: Note that the PGFl with the same form as Equation (27) is not an MB PGFl because of the additional sum according to different associations θ . A multi-Bernoulli PGFl is shown in Equation (13), which is the product of single Bernoulli PGFls.
The original MeMBer filter [3] approximates the derivative of , Moreover, in the derivation of the corrector for the original MeMBer filter, the CBMeMBer filter, and the I-MeMBer filter, the approximation that there is no more than one measurement near a true target is indispensable; a performance degradation occurs when targets in the surveillance region are too close, or the clutter intensity is very high. In contrast, the form as given in Equation (27) is exact. However, the combinatorial problem in this formulation is tractable. Thus, we utilize a Gibbs sampler to find the best associations θ .

PHD Intensity of the Posterior
Moreover, we approximate the posterior density with an MB distribution by calculating the PHD intensity (first-order moment). Thus, an MB RFS The PHD of the updated posterior is: We define the function: : The derivative of the functional 0, F h     in 1 h = WRT the target state and measurement set k Z is given as: where ( ) Thus, we rewrite Equation (28) as: Equation (32) shows that the PHD function of the updated posterior has a similar structure to that of an MB PHD function, which is given in Equation (10). Thus, the Bernoulli components of ( ) | k k π ⋅ can be calculated from Equation (32) where ( ) and ( ) From the right-hand side of Equation (33), we can see that there are many choices for the updated Bernoulli RFSs, but we can merge all the choices of associations for the Bernoulli components.
Compared to the LMB filter, our N-MB filter employs a straightforward prediction step. Thus, a different cost matrix needs to be designed to find the best data association. Moreover, the derivation of the LMB filter was based on approximating the GLMB density, which employs the labeled RFS notation. While our N-MB filter models the targets as an unlabeled multi-Bernoulli distribution, the Gibbs sampler is employed to get an approximate multi-Bernoulli distribution for the posterior density. Our N-MB filter was derived by using the theory of probability generating functionals and functional derivatives, which are very important tools for deriving RFS-based filters, for example, the PHD filter, CPHD filter, MeMBer filter, and PMBM filter.

Gaussian Implementation
In the previous section, we derived the recursion for the N-MB filter. Intuitively, the filter is intractable; therefore, approximations are necessary. In this section, we propose a closed-form Gaussian mixture (GM)-based implementation for the derived filter.
We assume that targets follow linear Gaussian dynamics, and measurements follow models with MB births: where F denotes the transition matrix, Q denotes the process noise covariance, H denotes the observation matrix, and R denotes the observation noise covariance. Furthermore, the survival and detection probability were assumed to be constant, which is given as The birth model is also an MB with the parameter set ⋅ is a GM of the form: , where ( )

Update
The prediction step of the proposed filter is the same with that in Baser et al. [25]. Given the predicted MB density and each pdf ( ) then, for each association θ , the updated posterior MB density is: , Then, the updated parameters in Equation (40) can be calculated using: where ( ) The derivation of the above steps involves analytically calculating the products of Gaussians and integrals of Gaussians using the standard results of Gaussian functions [9].

Gibbs Sampling
We employ the Gibbs sampling algorithm to solve the data association problem, which is an efficient algorithm for finding the best choices for the assignment of targets to measurements [28]. This subsection presents the implementation of the integrated Gibbs sampler in detail.
After the prediction step, there are N predicted Bernoulli components (representing N targets), and the number of measurements is z N . We define: , Consider a realization of a random variable θ that distributes according to a probability distribution π on { } 0 : . A Gibbs sampler is used to generate positive one-to-one vectors by independently sampling from π . The distribution π is: where Γ is the set of one-to-one mapping vectors in { } 0, It is difficult to sample from the distribution in Equation (49) directly. The Markov chain Monte Carlo (MCMC) method is the well-known algorithm for sampling from a complex distribution [29]. The Gibbs sampling algorithm is an efficient case of the so-called Metropolis-Hasting MCMC algorithm [30]. The proposed samples in a Gibbs sampler are always adopted; therefore, the Gibbs sampler is considered to be an efficient tool for generating the best associations [22,28,30,31]. Given results obtained using Equation (48), the assignment problem can be immediately solved by utilizing a Gibbs sampler.
Finally, the pseudocode for the filter is summarized in Algorithm 1.

Simulation
This section demonstrates the simulation results of the proposed N-MB filter conducted on linear Gaussian models. We set the simulation environments similar those in Baser et al. [25].  Figure 1 shows the true tracks of targets.

Targets Model
The kinematic state is given as ζ  x y y p p p p . The Gaussian target dynamics were modeled using the constant velocity (CV) model:  ; , Measurements followed the observation function in Equation (37) with parameters: where 10 m ε σ = was the std of the measurement noise. Clutters were modeled as a Poisson process with a clutter rate C λ and uniformly distributed according to the spatial probability density 1 /V over the surveillance region, where = × 6 2 4 10 m V was the "volume" of the surveillance region. We used a threshold = TH 0.7 r for the target state extraction and we removed those Gaussian components whose weights were lower than −5 10 . Those Gaussian components in a Bernoulli component whose distances were lower than 4 m were merged.

Performance Evaluation
Simulations were performed using MTT filters implemented in MATLAB (9.1.0.441655 R2016b, MathWorks, Beijing, China) on computers with an Intel Core i7-7700K CPU @ 4.20GHz and 16 GB RAM. We analyzed the performance in terms of the position estimates obtained by the proposed N-MB filer and CBMeMBer filter proposed in Baser et al. [25], the CPHD filter proposed in Vo et al. [13], and the track-oriented marginal MeMBer-Poisson (TOMB/P) filter proposed in García-Fernández et al. [17] using Monte Carlo simulations. We employed the optimal sub-pattern assignment (OSPA) distance as the error metric [32] to compare the performance of these filters. The OSPA is a distance metric for denoting the difference between two sets of points. The generalized OSPA (GOSPA), which was proposed recently, was also employed to verify the simulation's performance. GOSPA penalizes localization errors for detected targets and the errors due to missed and false targets [33], which is not considered in OSPA. The GOSPA metric with parameter 2 α = is given by which denotes the false-detection error.

Simulation Results
In Figure 2, the decomposed x -and y -components of the true tracks, state estimates generated by the N-MB filter, and measurements across all the time steps are shown versus time. The results showed that the proposed filter could correctly track the targets with individual motions throughout the targets' births and deaths. Three targets crossed with each other at time step 40 k = , and two targets crossed with each other at time step 60 k = , and the proposed filter had no difficulty handling these situation.  Figure 1) but randomly generated measurement data (the generated measurements over all time steps for one MC run is shown in Figure 2). Considering an environment with a high clutter intensity (the clutter rate was set to 80 C λ = and the detection probability was set as 0.98 D p = ), the mean and std of the estimated cardinality versus time for the CBMeMBer, CPHD, TOMB/P, and N-MB filters are shown in Figure 3. The CBMeMBer filter suffers from false estimations greatly. This phenomenon is easy to explain because the CBMeMBer employed the assumption that measurements were well separated [3,25]. In contrast, both TOMB/P and the proposed N-MB filter performed well in the cardinality estimation, and had a lower variance than the CPHD filter. Moreover, the MC average of the RMS OSPA distances are given in Figure 4, where the clutter rate was set to 20 C λ = and the detection probability was set to 0.7 D p = (i.e., a low signal-to-noise environment). The CBMeMBer filter suffered a great performance degradation when the detection probability was low. This was because the filter required the approximation of a high detection probability to correct for the negative term involved in the existence probability of the resulting MB components [25,34]. In contrast, the proposed N-MB filter performed well in the environment with a low detection probability. Compared to the TOMB/P filter, the N-MB filter had a better performance in most of time steps, but the TOMB/P filter outperformed the N-MB filter in terms of removing death targets at 71 k = . The CPHD filter was also not suitable for application in the low-signal-to-noise environment. In Figure 5, we give the averages of 200 MC trials for the RMS OSPA distance for the CBMeMBer, CPHD, TOMB/P, and N-MB filters for various clutter rates from 10 C λ = to 80 C λ = . In Figure 6, we give the results for various detection probabilities from 0.5 As expected, the OSPA error increased with the increased clutter rates and a lower detection probability. Overall, these results showed that the N-MB filter had a better performance compared to the other filters. The performance of the CBMeMBer filter was greatly influenced by the increased clutter rate  Although OSPA is a good method for measuring the performance of RFS-based MTT filters, we do not know whether the difference in performance of OSPA error comes from a better detection ability, a lower number of false targets, or a lower number of missed targets, which is key information regarding the analysis of filters. As mentioned in Section 5.2, we employed decomposed GOSPA components to verify the performance of filters in terms of location estimations, false estimations, and missed estimations. The RMS OSPA/GOSPA error averaged over all time steps of the filters in different situations are given in Figure 7 and Figure 8. These results prove that the N-MB filter had a similar performance to the TOMB/P filter regarding location accuracy, but outperformed the latter in evaluation of missed targets. The N-MB filter was more sensitive regarding initializing birth tracks; however, this feature meant the N-MB filter produced more false estimates. Overall, the N-MB filter outperformed the TOMB/P filter, especially in the high signal-to-noise environment. Again, the CBMeMBer filter showed its efficient performance in environment with a high signal-to-noise ratio and low clutter intensity, but the performance degraded dramatically in challenging environments.

Comparison with the LMB Filter
As summarized in Section 3.3, the proposed N-MB filter has a similar structure with the LMB filter, while it has a straightforward prediction step because the N-MB filter was derived based on the theory of PGFls and functional derivatives, while the LMB filter is based on approximating the GLMB density. An efficient implementation method that combines the prediction and update steps of the LMB filter is proposed in Mahler [23], which avoids the complex calculation burden of the prediction step. We compared the time-averaged GOSPA error and computation times (seconds) of the proposed N-MB filter and the LMB filter [23] in different situations, and the results are given in Table 1. Both filters employ the gating technique, which is widely used in many papers [15,16,22].
Intuitively, the proposed N-MB filter had a similar performance to the LMB filter in high signalto-noise environments. However, the N-MB filter had a more stable performance in more challenging environments (low signal-to-noise environments with a high clutter intensity). In contrast, a serious leak-tracking problem occurred in the tracking results of the LMB filter when applied in challenging environments.

Conclusions
Although the CBMeMBer filter was cardinality balanced, the derivation of which removes the second approximation employed in the original MeMBer filter [3]. It employs the first approximation in Mahler [3], which assumes that measurements are well-separated in the surveillance region. A filter based on this approximation cannot be applied in an environment with high clutter intensity. Moreover, the approximation of a high detection probability is required in the derivation of the CBMeMBer filter; therefore, the filter cannot be applied in the low signal-to-noise environment.
In this paper, we derived a novel-structured MB filter based on PGFl theory and unlabeled RFSs. The proposed filter removed the assumptions utilized in deriving the MeMBer filter. Thus, the filter performed well in the low signal-to-noise and high clutter density environments. Compared to filters based on a GLMB distribution, the proposed filter did not extract explicit tracks for targets. The theory of tracking based on sets of trajectories [35] was utilized to solve this problem. Various simulation results showed that the proposed filter outperformed the LMB filter in low signal-to-noise environments with a little higher computational cost.