Multiple Target, Multiple Type Filtering in the RFS Framework

TL;DR

This paper introduces a novel N-type PHD filter within the RFS framework to effectively handle multiple target types and their confusions, improving tracking accuracy in complex scenarios.

Contribution

It extends the standard PHD filter to multiple target types with confusions and develops a Gaussian mixture implementation for practical use.

Findings

Improved tracking performance with target confusions.

Effective discrimination of multiple target types.

Validation through simulations with 16 targets of 4 types.

Abstract

A Multiple Target, Multiple Type Filtering (MTMTF) algorithm is developed using Random Finite Set (RFS) theory. First, we extend the standard Probability Hypothesis Density (PHD) filter for multiple types of targets, each with distinct detection properties, to develop a multiple target, multiple type filtering, N-type PHD filter, where $N\geq2$, for handling confusions among target types. In this approach, we assume that there will be confusions between detections, i.e. clutter arises not just from background false positives, but also from target confusions. Then, under the assumptions of Gaussianity and linearity, we extend the Gaussian mixture (GM) implementation of the standard PHD filter for the proposed N-type PHD filter termed the N-type GM-PHD filter. Furthermore, we analyze the results from simulations to track sixteen targets of four different types using a four-type (quad) GM-PHD filter as a typical example and compare it with four independent GM-PHD filters using the Optimal Subpattern Assignment (OSPA) metric. This shows the improved performance of our strategy that accounts for target confusions by efficiently discriminating them.

Figures

Tables

Table 1. TABLE I: OSPA error at different values of probabilities of confusion p 12 , D subscript 𝑝 12 𝐷 p_{12,D} , p 13 , D subscript 𝑝 13 𝐷 p_{13,D} , p 14 , D subscript 𝑝 14 𝐷 p_{14,D} , p 21 , D subscript 𝑝 21 𝐷 p_{21,D} , p 23 , D subscript 𝑝 23 𝐷 p_{23,D} , p 24 , D subscript 𝑝 24 𝐷 p_{24,D} , p 31 , D subscript 𝑝 31 𝐷 p_{31,D} , p 32 , D subscript 𝑝 32 𝐷 p_{32,D} , p 34 , D subscript 𝑝 34 𝐷 p_{34,D} , p 41 , D subscript 𝑝 41 𝐷 p_{41,D} , p 42 , D subscript 𝑝 42 𝐷 p_{42,D} and p 43 , D subscript 𝑝 43 𝐷 p_{43,D} (0.3, 0.6 and 0.9) for quad GM-PHD filter and 4 independent GM-PHD filters. Time taken is given in brackets.

Method	0.3	0.6	0.9
Quad GM-PHD filter	28.70m (239.82sec)	28.81m (251.27sec)	29.17m (270.13sec)
4 GM-PHD filters	32.18m (147.42sec)	46.47m (152.06sec)	55.86m (160.17sec)

Equations

$$X_{i,k}=\{x_{i,k,1},...x_{i,k,L_i(k)}\}\in \mathcal{F}(\mathcal{X})$$

$$Z_{j,k}=\{z_{j,k,1},...z_{j,k,M_j(k)}\}\in \mathcal{F}(\mathcal{Z})$$

$${\Xi }_{i,k}=S_{i,k}(X_{i,k-1})\cup {\Gamma }_{i,k},$$

$${\Omega }_{ji,k}={\Theta }_{j,k}(X_{i,k})\cup C_{si,k}\cup C_{tiJ,k},$$

$$\begin{array}[]{lll}p_{i,k|k-1}(X_{i,k}|Z_{j,1:k-1})=\\ \int y_{i,k|k-1}(X_{i,k}|X)p_{i,k-1|k-1}(X|Z_{j,1:k-1})\mu(dX)\end{array}$$

$$\begin{array}[]{lll}p_{i,k|k}(X_{i,k}|Z_{j,1:k})=\frac{f_{ji,k}(Z_{j,k}|X_{i,k})p_{i,k|k-1}(X_{i,k}|Z_{j,1:k-1})}{\int f_{ji,k}(Z_{j,k}|X)p_{i,k|k-1}(X|Z_{j,1:k-1})\mu(dX)}\end{array}$$

$$F_i[g,h]=G_{T_i}(h{G}_{L_{i,i}}(g|.))G_{c_i}(g)\prod _{j=1\setminus i}^N{G}_{T_j}(G_{L_{j,i}}(g|.)),$$

$$G_{c_i}(h)=\exp ({\lambda }_i(c_i[g]-1)),$$

$$G_{T_i}(h)=\exp ({\mu }_i(s_i[g]-1)),$$

$$G_{L_{j,i}}(g|x)=1-p_{ji,D}(x)+p_{ji,D}(x)\int g(z)f_{ji}(z|x)dz,$$

$$\begin{array}[]{lll}s_{i}[hG_{L_{i,i}}(g|x)]=&\displaystyle\int s_{i}(x)h(x)\Big{(}1-p_{ii,D}(x)+\\ &p_{ii,D}(x)\int g(z)f_{ii}(z|x)dz\Big{)}dx,\end{array}$$

$$\begin{array}[]{lll}s_{j}[G_{L_{j,i}}(g|x)]=&\displaystyle\int s_{j}(x)\Big{(}1-p_{ji,D}(x)+\\ &p_{ji,D}(x)\int g(z)f_{ji}(z|x)dz\Big{)}dx,\end{array}$$

$$\begin{array}[]{lll}F_{i}[g,h]=&\exp\Bigg{(}\lambda_{i}\Big{(}\int g(z)c_{i}(z)dz-1\Big{)}+\sum_{j=1\setminus i}^{N}\mu_{j}\Big{[}\displaystyle\int s_{j}(x)\Big{(}1-p_{ji,D}(x)+p_{ji,D}(x)\int g(z)f_{ji}(z|x)dz\Big{)}dx-1\Big{]}\\ &+\mu_{i}\Big{[}\displaystyle\int s_{i}(x)h(x)\Big{(}1-p_{ii,D}(x)+p_{ii,D}(x)\int g(z)f_{ii}(z|x)dz\Big{)}dx-1\Big{]}\Bigg{)},\end{array}$$

$$G_i(h|z_1,...,z_{M_j})=\frac{\frac{{\delta }^{M_j}}{{\delta }_{{\phi }_{z_1}}...{\delta }_{{\phi }_{z_{M_j}}}}{F}_i[g,h]{|}_{g=0}}{\frac{{\delta }^{M_j}}{{\delta }_{{\phi }_{z_1}}...{\delta }_{{\phi }_{z_{M_j}}}}{F}_i[g,1]{|}_{g=0}},$$

$$\begin{array}[]{lll}\mathcal{D}_{i}(x|z_{1},...,z_{M_{j}})&=\frac{\delta}{\delta_{\varphi_{x}}}G_{i}(h|z_{1},...z_{M_{j}})|_{h=1},\\ &=\mu_{i}s_{i}(x)(1-p_{ii,D}(x))+\sum_{m=1}^{M_{j}}\frac{\mu_{i}s_{i}(x)p_{ii,D}(x)f_{ii}(z_{m}|x)}{\lambda_{i}c_{i}(z_{m})+\sum_{j=1\setminus i}^{N}\mu_{j}\int s_{j}(x)p_{ji,D}(x)f_{ji}(z_{m}|x)dx+\mu_{i}\int s_{i}(x)p_{ii,D}(x)f_{ii}(z_{m}|x)dx},\end{array}$$

$$\begin{array}[]{lll}\frac{\delta}{\delta_{\varphi_{x}}}\lambda_{i}c_{i}(z_{m})&=0,\\ ~{}~{}(\text{No~{}process~{}function~{}h(x)})\\ \frac{\delta}{\delta_{\varphi_{x}}}\mu_{j}\int s_{j}(x)p_{ji,D}(x)f_{ji}(z_{m}|x)dx&=0,\\ ~{}~{}(\text{No~{}process~{}function~{}h(x))}\\ \frac{\delta}{\delta_{\varphi_{x}}}\mu_{i}\int s_{i}(x)h(x)p_{ii,D}(x)f_{ii}(z_{m}|x)dx&=\\ ~{}~{}\mu_{i}s_{i}(x)p_{ii,D}(x)f_{ii}(z_{m}|x).\end{array}$$

$$E[|({\Xi }_1\cup {\Xi }_2...\cup {\Xi }_N)\cap R|]=\sum _{i=1}^N{\int }_R{\mathcal{D}}_{{\Xi }_i}(x)dx$$

$$\begin{array}[]{lll}\mathcal{D}_{i,k|k-1}(x)=&\int p_{i,S,k|k-1}(\zeta)y_{i,k|k-1}(x|\zeta)\mathcal{D}_{i,k-1|k-1}(\zeta)d\zeta\\ &+\gamma_{i,k}(x),\end{array}$$

$$\begin{array}[]{lll}\mathcal{D}_{i,k|k}(x)&=\bigg{[}1-p_{ii,D}(x)+\sum_{z\in Z_{i,k}}\frac{p_{ii,D}(x)f_{ii,k}(z|x)}{c_{s_{i,k}}(z)+c_{t_{i,k}}(z)+\int p_{ii,D}(\zeta)f_{ii,k}(z|\zeta)\mathcal{D}_{i,k|k-1}(\zeta)d\zeta}\bigg{]}\mathcal{D}_{i,k|k-1(x)},\end{array}$$

$$\begin{array}[]{lll}c_{t_{i,k}}(z)=\sum_{j\in\{1,...,N\}\setminus i}\int p_{ji,D}(y)\mathcal{D}_{j,k|k-1}(y)f_{ji,k}(z|y)dy,\end{array}$$

$$c_{s_{i,k}}(z)={\lambda }_i{c}_i(z)={\lambda }_{c_i}A{c}_i(z),$$

$$y_{i,k|k-1}(x|\zeta )=\mathcal{N}(x;F_{i,k-1}\zeta ,Q_{i,k-1})$$

$$f_{ji,k}(z|x)=\mathcal{N}(z;H_{ji,k}x,R_{ji,k})$$

$$\begin{array}{r}{\gamma }_{i,k}(x)=\sum _{v=1}^{V_{{\gamma }_i,k}}{w}_{i,\gamma ,k}^{(v)}\mathcal{N}(x;m_{i,\gamma ,k}^{(v)},P_{i,\gamma ,k}^{(v)})\end{array}$$

$$\begin{array}{r}{\mathcal{D}}_{i,k-1}(x)=\sum _{v=1}^{V_{i,k-1}}{w}_{i,k-1}^{(v)}\mathcal{N}(x;m_{i,k-1}^{(v)},P_{i,k-1}^{(v)}),\end{array}$$

$${\mathcal{D}}_{i,k|k-1}(x)={\mathcal{D}}_{i,S,k|k}(x)+{\gamma }_{i,k}(x),$$

$$\begin{array}[]{lll}\mathcal{D}_{i,S,k|k-1}(x)=&p_{i,S,k}\sum_{v=1}^{V_{i,k-1}}w_{i,k-1}^{(v)}\mathcal{N}(x;\\ &m_{i,S,k|k-1}^{(v)},P_{i,S,k|k-1}^{(v)}),\end{array}$$

$$m_{i,S,k|k-1}^{(v)}=F_{i,k-1}{m}_{i,k-1}^{(v)},$$

$$P_{i,S,k|k-1}^{(v)}=Q_{i,k-1}+F_{i,k-1}{P}_{i,k-1}^{(v)}{F}_{i,k-1}^T,$$

$$\begin{array}{r}{\mathcal{D}}_{i,k|k-1}(x)=\sum _{v=1}^{V_{i,k|k-1}}{w}_{i,k|k-1}^{(v)}\mathcal{N}(x;m_{i,k|k-1}^{(v)},P_{i,k|k-1}^{(v)}),\end{array}$$

Full text

I Introduction

Multi-target filtering is a state estimation problem which plays a key role in visual, radar and sonar tracking, robot simultaneous localization and mapping (SLAM), and other signal processing applications. Traditionally, multi-target filters are based on finding associations between targets and measurements using methods including Global Nearest Neighbour (GNN) [1][2], Joint Probabilistic Data Association Filter (JPDAF) [1][3], and Multiple Hypothesis Tracking (MHT) [1][4]. However, these approaches face challenges not only in the uncertainty caused by the data association but also in computational growth exponential to the number of targets and measurements. To address the complexity problem, a unified framework directly extended single- to multi-target tracking by representing multi-target states and observations as a random finite set (RFS) [5]. This estimates both the states and cardinality of an unknown and time varying number of targets in a scene, and includes birth, death, clutter (false alarms), and missed detections. Mahler [5] propagated the first-order moment of the multi-target posterior called Probability Hypothesis Density (PHD) or intensity rather than the full multi-target posterior, thus having much lower computational complexity in the single state space rather than in the joint-state space as in traditional methods. Further developments include the Gaussian mixture (GM-PHD) [6] and the Sequential Monte Carlo (SMC) (particle-PHD filter) implementations [7] which have been applied to visual tracking in [8] and [9], respectively. This approach is flexible, for instance, it has been used to find the detection proposal with the maximum weight as the target position estimate for tracking a target of interest in dense environments by removing the other detection proposals as clutter [10][11]. Furthermore, the PHD filter has also been used for doing visual odometry (VO) [12] and SLAM [13] in robotics. Joint detection, tracking and classification (JDTC) of multiple targets in clutter which jointly estimates the number of targets, their kinematic states, and types of targets (classes) from a sequence of noisy and cluttered observation sets was developed using the PHD filter in [14]. In this approach, the dynamics of each target type (class) is modelled as a class-dependent model set, and the signal amplitude is included within the multi-target likelihood in the PHD-like filter to enhance the discrimination between targets from different classes and false alarms. Similarly, a joint target tracking and classification (JTC) algorithm was developed in [15] using RFS which takes into account extraneous target-originated measurements (of the same type) i.e. multiple measurements originating from a target which can be modeled as a Poisson RFS. In these approaches, the augmented state vector of a target comprises the target kinematic state and class label, i.e. the target type (class) is put into the target state vector. However, all these RFS-based multiple target filters were developed for either a single target type or multiple target type but without taking any account of target confusions between target types at the measurement stage, i.e. measurements originated not only from the same target type but were also confused from the other target types.

Practically, there are many situations where tracking and discrimination of multiple target types is essential, handling confusions between target types. For example, when developing situational awareness for driver assistance and vehicle autonomy [16], a vehicle equipped with a sensor suite must detect and track other road users to select the best sensor focus and course of action, usually concentrating on other vehicles, pedestrians and bicycles. In this particular and many other examples, confusion between target types is common, for example a standard pedestrian detection strategy [17] often provides confused detections between pedestrians and cyclists, and even small cars. Moreover, for sports analysis we often want to track and discriminate sub-groups of the same target type such as players in opposing teams [18]. These types of problems motivate our work to develop a multiple target, multiple type filtering methodology handling target confusions following the RFS based filtering method, particularly the standard PHD filter, without involving any data association.

The main contributions of this paper are as follows.

• •

We model the RFS filtering of $N$ different types of multiple targets with separate but confused detections which we call the N-type PHD filter where $N\geq 2$.

• •

The Gaussian mixture implementation of the standard PHD filter is extended for the proposed N-type PHD filter.

• •

We demonstrate this proposed N-type GM-PHD filter by simulations, specifically for a quad GM-PHD filter ($N=4$) as a typical example under different values of confusion detection probabilities to show that our approach yields improved performance over the standard approach.

We presented preliminary ideas in [19] for three different target types (N=3) for visual tracking applications. In this work, we further develop our approach. We extend from a tri-PHD filter to a N-type PHD filter as well as conducting experiments on more dense simulations. The mathematical proof of the algorithm is presented. We demonstrate the quad GM-PHD filter for four types (N = 4) of sixteen targets with detailed analysis as a typical simulation example under different values of confusion detection probabilities.

The remainder of this paper is organized as follows. Multiple type, multiple target recursive Bayes filtering with RFS is described in section II. A probability generating functional for deriving a N-type PHD filter and the N-type PHD filtering strategy are given in sections III and IV, respectively. In section V, a Gaussian mixture implementation of the N-type PHD filter is described in detail. The experimental results are analyzed and compared in section VI. The main conclusions and suggestions for future work are summarized in section VII.

II Multiple Target, Multiple Type Recursive Bayes Filtering with RFS

A RFS represents a varying number of non-ordered target states and observations which is analogous to a random vector for single target tracking. More precisely, a RFS is a finite-set-valued random variable i.e. a random variable which is random in both the number of elements and the values of the elements themselves. Finite Set Statistics (FISST), the study of statistical properties of RFS, is a systematic treatment of multi-sensor multi-target filtering as a unified Bayesian framework using random set theory [5].

When different detectors run on the same scene to detect different target types, there is no guarantee that these detectors only detect their own type. It is possible to run an independent PHD filter for each target type, but this will not be correct in most cases, as the likelihood of a positive response to a target of the wrong type will in general be different from, usually higher than, the likelihood of a positive response to the scene background. In this paper, we account for this difference between background clutter and target type confusion. This is equivalent to a single sensor (e.g. a smart camera) that has N different detection modes, each with its own probability of detection and a measurement density for N different target types.

So we model a N-type PHD filter to filter N-types of multiple targets in such a way that the first PHD update will filter the first target type treating the others as potential, additional clutter in addition to background clutter, and vice versa. In the joint tracking and classification approaches such as in [14],[15], the target type (class) is put into the target state vector, however, here we follow a different approach which is convenient for handling target confusions from different target types. Accordingly, to derive the N-type PHD filter, it is necessary to first give its RFS representation to extend from a single type single-target Bayes framework to multiple type multi-target Bayes framework. Let the multi-target state space $\mathcal{F}(\mathcal{X})$ and the multi-target observation space $\mathcal{F}(\mathcal{Z})$ be the respective collections of all the finite subsets of the state space $\mathcal{X}$ and observation space $\mathcal{Z}$, respectively. If $L_{i}(k)$ is the number of targets of target type $i$ in the scene at time $k$, then the multiple states for target type $i$, $X_{i,k}$, is the set

[TABLE]

where $i\in\{1,...,N\}$. Similarly, if $M_{j}(k)$ is the number of received observations from detector $j$, then the corresponding multiple target measurements is the set

[TABLE]

where $j\in\{1,...,N\}$. Some of these observations may be false i.e. due to clutter (background) or confusion (response due to another target type).

The uncertainty in the state and measurement is introduced by modeling the multi-target state and the multi-target measurement using RFS. Let $\Xi_{i,k}$ be the RFS associated with the multi-target state of target type $i$, then

[TABLE]

where $S_{i,k}(X_{i,k-1})$ denotes the RFS of surviving targets of target type $i$, and $\Gamma_{i,k}$ is the RFS of the new-born targets of target type $i$. We do not consider target spawning in this paper.

Further, the RFS $\Omega_{ji,k}$ associated with the multi-target measurements of target type $i$ from detector $j$ is

[TABLE]

where $J=\{1,...,N\}\setminus i$ and $\Theta_{j,k}(X_{i,k})$ is the RFS modeling the measurements generated by the target $X_{i,k}$ from detector $j$, and $C_{si,k}$ models the RFS associated with the clutter (false alarms) for target type $i$ which comes from the scene background (using detector $j$). However, we now must also include $C_{tiJ,k}$ which is the RFS associated with measurements of all target types $J=\{1,...,N\}\setminus i$ being treated as confusion while filtering target type $i$ i.e. measurements of all target types are included into clutter (confusion) except measurement of target type $i$ while filtering target type $i$.

Analogously to the single-target case, the dynamics of $\Xi_{i,k}$ are described by the multi-target transition density $y_{i,k|k-1}(X_{i,k}|X_{i,k-1})$, while $\Omega_{ji,k}$ is described by the multi-target likelihood $f_{ji,k}(Z_{j,k}|X_{i,k})$ for multiple target type $i\in\{1,...,N\}$ from detector $j\in\{1,...,N\}$. The recursive equations are

[TABLE]

[TABLE]

where $\mu$ is an appropriate dominating measure on $\mathcal{F}(\mathcal{X})$ [5]. Though a Monte Carlo approximation of this optimal multi-target types Bayes recursion is possible considering multiple targets of a single type [7], the number of particles required is exponentially related to the number of targets and their types in the scene. To make it computationally tractable, we extended Mahler’s method of propagating the first-order moment of the multi-target posterior instead of the full multi-target posterior as its approximation called the Probability Hypothesis Density (PHD) [5], $\mathcal{D}_{i,k|k}(x|Z_{j,1:k})=\mathcal{D}_{i,k|k}(x)=\int\delta_{x}(x)p_{i,k|k}(X|Z_{j,1:k})\delta X$ where $\delta_{x}(x)=\sum_{w\in x}\delta_{w}(x)$, for $N\geq 2$ types of multiple targets by deriving the updated PHDs from Probability Generating Functionals (PGFLs) starting from the standard predicted PHDs for each target type.

III Probability generating functional (PGFL)

The probability generating functional is a convenient representation for stochastic modelling with a point process [5], a type of random process for which any one realisation consists of a set of isolated points either in time or space. Now, we model joint (probability generating) functionals which take into account the clutter due to the other target types in addition to the background for deriving the updated PHDs. Starting from the standard proved predicted PHDs in [5] (refer to the appendix VIII for the proof) for each multi-target type, we will derive novel extensions for the updated PHDs of N-type PHD filter from PGFLs for each target type for handling confusions between target types.

The joint functional for target type $i$ treating all other target types as clutter is given by

[TABLE]

where $i\in\{1,...,N\}$ denotes target type, $g$ is related to the target measurement process and $h$ is related to the target state process.

[TABLE]

where $G_{c_{i}}(h)$ is the Poisson PGFL [5] for false alarms where $\lambda_{i}$ is the average number of false alarms for target type $i$ and the functional $c_{i}[g]=\int g(z)c_{i}(z)dz$ where $c_{i}(.)$ is the uniform density over the surveillance region;

[TABLE]

where $G_{T_{i}}(h)$ is the prior PGFL and $\mu_{i}$ is the average number of targets, each of which is distributed according to $s_{i}(x)$ for target type $i$; and

[TABLE]

where $G_{L_{j,i}}(g|x)$ is the Bernoulli detection process for each target of target type $i$ using detector $j$ with probability of detection for target type $i$ by detector $j$, $p_{ji,D}$, and $f_{ji}(z|x)$ is a likelihood defining the probability that $z$ is generated by the target type $i$ conditioned on state $x$ from detector $j$ [5]. Expanding $s_{i}[hG_{L_{i,i}}(g|x)]$ and $s_{j}[G_{L_{j,i}}(g|x)]$ as

[TABLE]

and

[TABLE]

Accordingly, $F_{i}[g,h]$ is expanded as

[TABLE]

The updated PGFL $G_{i}(h|z_{1},...z_{M_{j}})$ for target type $i$ is obtained by finding the $M_{j}^{th}$ functional derivative of $F_{i}[g,h]$ [5] and is given by

[TABLE]

The updated PHD for target type $i$ treating all other target types as clutter can be obtained by taking the first-order moment (mean) [5] of Eq. (14) and setting $h=1$,

[TABLE]

Thus, $\mathcal{D}_{i}(x|z_{1},...,z_{M_{j}})$ in Eq. (15), is the updated PHD for target type $i$ treating all other target types as clutter. The term $\mu_{i}s_{i}(x)$ in Eq. (15) is the predicted PHD for target type $i$ (refer to the appendix VIII for the proof). The following equational facts are important for the derivation of Eq. (14) to get Eq. (15) (refer to [20], Chapter 4, for the full derivation):

[TABLE]

IV N-type PHD Filtering Strategy

Here we state PHD recursions in a generic form for multiple target, multiple type filtering with $Z_{1,k},...,Z_{N,k}$ separate but confused multi-target measurements between different target types, i.e. the N-type PHD filter, where $N\geq 2$. For N types of multiple targets, the PHDs, $\mathcal{D}_{\Xi_{1}}(x)$, $\mathcal{D}_{\Xi_{2}}(x)$, …, $\mathcal{D}_{\Xi_{N}}(x)$, are the first-order moments of RFSs, $\Xi_{1}$, $\Xi_{2}$, … $\Xi_{N}$, and they are intensity functions on a single state space $\mathcal{X}$ whose peaks identify the likely positions of the targets. For any region $R\subseteq\mathcal{X}$

[TABLE]

where$|.|$ is used to denote the cardinality of a set. In practice, Eq.(17) means that by integrating the PHDs on any region $R$ of the state space, we obtain the expected number of targets (cardinality) in $R$.

Accordingly, the Bayesian iterative prediction and update of the N-type PHD filtering strategy is given as follows.

The PHD prediction for target type $i$ is defined as

[TABLE]

where $\gamma_{i,k}(.)$ is the intensity function of a new target birth RFS $\Gamma_{i,k}$, $p_{i,S,k|k-1}(\zeta)$ is the probability that the target still exists at time $k$, $y_{i,k|k-1}(.|\zeta)$ is the single target state transition density at time $k$ given the previous state $\zeta$ for target type $i$.

Thus, following Eq. (15), the final updated PHD for target type $i$ is obtained by setting $\mu_{i}s_{i}(x)=\mathcal{D}_{i,k|k-1}(x)$,

[TABLE]

The clutter intensity $c_{t_{i,k}}(z)$ due to all types of targets $j\in\{1,...,N\}$ except target type $i$ in Eq.(19) is given by

[TABLE]

This means that when we are filtering target type $i$, all the other target types will be included as clutter. Eq.(20) converts state space to observation space by integrating the PHD estimator $\mathcal{D}_{j,k|k-1}(y)$ and likelihood $f_{ji,k}(z|y)$ which defines the probability that $z$ is generated by detector $j$ conditioned on state $x$ of the target type $i$ taking into account the confusion probability $p_{ji,D}(y)$, the detection probability for target type $i$ by detector $j$. Hence, it maps the state space of wrongly detected targets to the measurement space and treats them as clutter.

The clutter intensity due to the scene $i$, $c_{s_{i,k}}(z)$, in Eq. (19) is given by

[TABLE]

where $c_{i}(.)$ is the uniform density over the surveillance region $A$, and $\lambda_{c_{i}}$ is the average number of clutter returns per unit volume for target type $i$ i.e. $\lambda_{i}=\lambda_{ci}A$. While the PHD filter has linear complexity with the current number of measurements ($m$) and with the current number of targets ($n$) i.e. computational order of $O(mn)$, the N-type PHD filter has linear complexity with the current number of measurements ($m$), with the current number of targets ($n$) and with the total number of target types ($N$) i.e. computational order of $O(mnN)$ as can be seen from Eq. (19).

In general, the clutter intensity due to the background for target type $i$, $c_{s_{i,k}}(z)$, can be different for each target type as they depend on the receiver operating characteristic (ROC) curves of the detection processes. Moreover, the probabilities of detection $p_{ii,D}(x)$ and $p_{ji,D}(x)$ may all be different although assumed constant across both the time and space continua.

The pictorial summary of recursive N-type PHD filter structure is given in Fig. 1. The intersection (common meeting point) of the N-type PHD filtering process is at Eq. (20) which involves the predicted PHDs and likelihoods of all other target types along with their corresponding probabilities of confusion when filtering target type $i$. This equation converts state space to observation space by integrating the PHD estimators (predicted PHDs) and likelihoods along with their corresponding probabilities of confusion to treat as clutter as clutter operates on the measurement space. Actually, Eq. (20) is called in Eq. (19) as a clutter along with a clutter due to the scene for updating PHD of target type $i$. At this point, the confusion of target types will be solved. Except this, every updated PHD retains its own type i.e. the updated PHDs are not combined. Initialization in Fig. 1 shows the measurement-driven birth of targets. Thus, given the measurements of all target types at time $k$ and the updated PHDs at time $k-1$ ($\mathcal{D}_{i,k-1}=\mathcal{D}_{i,k-1|k-1}$ where $i\in\{1,...,N\}$), the filtering process occurs recursively to get the updated PHDs at time $k$ ($\mathcal{D}_{i,k}=\mathcal{D}_{i,k|k}$ where $i\in\{1,...,N\}$). Obviously, the updated PHDs at time $k-1$ are the PHDs after merging and pruning of the updated PHDs at the previous iteration (refer to section V) which will be predicted as existing targets in the next iteration (at time $k|k-1$). The estimated states are extracted from the updated PHDs of each target type (with weights greater than 0.5, refer to section V). In fact, after states extraction from the updated PHDs, the extracted states can be combined which do not affect the filtering process afterwards. However, for labelling of each target as in our preliminary work [19], keeping extracted states of each target type separately is important i.e. we apply the Hungarian algorithm to the extracted states of each target type separately for frame to frame labelling in [19], of course keeping the unique identities of all targets of the target types.

V Gaussian Mixture-Based N-type PHD Filter Implementation

The Gaussian mixture implementation of the standard PHD (GM-PHD) filter [6] is a closed-form solution of the PHD filter with the assumptions of a linear Gaussian system. In this section, this standard implementation is extended for the N-type PHD filter, more importantly solving Eq. (20). We assume each target follows a linear Gaussian model.

[TABLE]

[TABLE]

where $\mathcal{N}(.;m,P)$ denotes a Gaussian density with mean $m$ and covariance $P$; $F_{i,k-1}$ and $H_{ji,k}$ are the state transition and measurement matrices, respectively. $Q_{i,k-1}$ and $R_{ji,k}$ are the covariance matrices of the process and the measurement noises, respectively, where $i\in\{1,...,N\}$ and $j\in\{1,...,N\}$. Besides, a current measurement driven birth intensity inspired by but not identical to [21] is introduced at each time step, removing the need for the prior knowledge (specification of birth intensities) or a random model, with a non-informative zero initial velocity. The intensity of the spontaneous birth RFS is $\gamma_{i,k}(x)$ for target type $i$

[TABLE]

where $V_{\gamma_{i},k}$ is the number of birth Gaussian components for target type $i$ where $i\in\{1,...,N\}$, $m_{i,\gamma,k}^{(v)}$ is the current measurement (noisy version of position) and zero initial velocity used as mean and $P_{i,\gamma,k}^{(v)}$ is the birth covariance for a Gaussian component $v$ of target type $i$.

It is assumed that the posterior intensities for target type $i$ at time $k-1$ are Gaussian mixture of the form

[TABLE]

where $i\in\{1,...,N\}$ and $V_{i,k-1}$ is the number of Gaussian components of $\mathcal{D}_{i,k-1}(x)$. Under these assumptions, the predicted intensities at time $k$ for target type $i$ are given following Eq. (18) by

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

where $\gamma_{i,k}(x)$ is given by (24).

Since $\mathcal{D}_{i,S,k|k-1}(x)$ and $\gamma_{i,k}(x)$ are Gaussian mixtures, $\mathcal{D}_{i,k|k-1}(x)$ can be expressed as Gaussian mixture of the form

[TABLE]

where $w_{i,k|k-1}^{(v)}$ is the weight accompanying the predicted Gaussian component $v$ for target type $i$ and $V_{i,k|k-1}$ is the number of predicted Gaussian components for target type $i$ where $i\in\{1,...,N\}$.

Now, assuming the probabilities of detection to be constant i.e. $p_{ii,D}(x)=p_{ii,D}$, the final updated PHD for target type $i$ is given as follows. Accordingly, the posterior intensity for target type $i$ at time $k$ (updated PHD) treating all other target types as clutter is also a Gaussian mixture which corresponds to Eq. (19), and is given by

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

$c_{s_{i,k}}(z)$ is given in Eq. (21). Therefore, all that left is to formulate the implementation scheme for $c_{t_{i,k}}(z)$ which is given in (20) and is given again as

[TABLE]

where $\mathcal{D}_{j,k|k-1}(y)$ is given in Eq. (27), $f_{ji,k}(z|y)$ is given in Eq. (23) and $p_{ji,D}(y)$ is assumed constant i.e. $p_{ji,D}(y)=p_{ji,D}$. Since $w_{j,k|k-1}^{(v)}$ is independent of the integrable variable $y$, Eq. (29) becomes

[TABLE]

This can be simplified further using the following equality given that $P_{1}$ and $P_{2}$ are positive definite

[TABLE]

Therefore, Eq. (30) becomes,

[TABLE]

where $i\in\{1,...,N\}$.

The key steps of the N-type GM-PHD filter are summarised in Algorithms 1 and 2. The number of Gaussian components in the posterior intensities may increase without bound as time progresses. To keep the number of N-type GM-PHD components to a reasonable level after the measurement update, it is necessary to prune weak and duplicated components. First, weak components with weight $w_{k}^{v}<10^{-5}$ are pruned. Further, Gaussian components with Mahalanobis distance less than $U=4m$ from each other are merged. These pruned and merged Gaussian components, the output of Algorithm 2, will be predicted as existing targets in the next iteration. Finally, Gaussian components of the posterior intensity, the output of Algorithm 1, with means corresponding to weights greater than 0.5 as a threshold are selected as multi-target state estimates.

VI Experimental Results

In this section, a simulation filtering example using a quad GM-PHD filter for four different types of multiple targets is analyzed. We have also made experimental simulation analyses of a dual GM-PHD filter ($N=2$) and a tri-GM-PHD filter ($N=3$) and then applied the tri-GM-PHD filter for visual tracking applications in [19] for three types of targets. In this experiment, we demonstrate the quad GM-PHD filter (N = 4) with detailed analysis as a typical simulation example under different values of confusion detection probabilities. Accordingly, we define a sequence of 120 frames with sixteen trajectories that emanate from four types of targets that appear in the scene at different positions and time steps (frames), as shown in Fig. 3 and Fig. 4. This is a typical example of not only a higher number of target types (four) but also an example of a dense scene i.e. it consists of trajectories of 16 targets in the same scene with many crossings. Our results are obtained after running 50 simulations i.e. the number of Monte Carlo (MC) is 50. Obviously, the goal of a N-type PHD filter is to handle confusions among $N\geq 2$ different target types; not to deal with sparse or dense targets in the scene. With regards to sparse or dense targets in the scene, it has the same characteristics as the standard PHD filter.

The initial locations and covariances for all target types are given by Eq. (33) as follows.

[TABLE]

The state vector $x_{k}=[p_{x,xk},p_{y,xk},\dot{p}_{x,xk},\dot{p}_{y,xk}]^{T}$ consists of position $(p_{x,xk},p_{y,xk})$ and velocity $(\dot{p}_{x,xk},\dot{p}_{y,xk})$, and the measurement is a noisy version of the position, $z_{k}=[p_{x,zk},p_{y,zk}]^{T}$. Each of the target trajectories follows a linear Gaussian dynamic model of Eq. (22) with matrices

[TABLE]

[TABLE]

where $I_{n}$ and $0_{n}$ denote the $n\times n$ identity and zero matrices, respectively. $\Delta=1s$ is the sampling period. $\sigma_{v_{i}}=5m/s^{2}$ where $i\in\{1,2,3,4\}$ is the standard deviation of the process noise for target type $i$.

For the algorithm, we assume each target has a survival probability $p_{1,S}=p_{2,S}=p_{3,S}=p_{4,S}=0.99$. The probabilities of detection are $p_{11,D}=0.90$ , $p_{22,D}=p_{33,D}=0.92$, $p_{44,D}=0.91$, and different values of confusion detection probabilities (0.0, 0.3, 0.6 and 0.9) are analyzed for $p_{12,D}$, $p_{13,D}$, $p_{14,D}$, $p_{21,D}$, $p_{23,D}$, $p_{24,D}$, $p_{31,D}$, $p_{32,D}$, $p_{34,D}$, $p_{41,D}$, $p_{42,D}$ and $p_{43,D}$.

The measurement follows the observation model of Eq. (23) with matrices

[TABLE]

where $\sigma_{r_{ii}}=\sigma_{r_{ji}}=6m$ ($i\in\{1,2,3,4\}$ and $j\in\{1,2,3,4\}$) is the standard deviation of the measurement noise.

Since there are many targets in the scene, we use about 40 clutter returns (10 for each target type) over the surveillance region. A current measurement driven birth intensity inspired by but not identical to [21] is introduced at each time step, removing the need for the prior knowledge (specification of birth intensities) or a random model, with a non-informative zero initial velocity. At birth, Gaussian components of each target type have a corresponding initial weight $w_{1,\gamma,k}^{(i)}=w_{2,\gamma,k}^{(i)}=w_{3,\gamma,k}^{(i)}=w_{4,\gamma,k}^{(i)}=3\times 10^{-6}$. This very small initial weight is assigned to the Gaussian components for new births as this is effective for high clutter rates. This is basically equivalent to the average number of appearing (birth) targets per time step ($n_{b}$) divided uniformly across the surveillance region ($A$).

The configuration of the detectors is shown Fig. 2. As shown in this figure, detector 1 detects target type 1 (targets 1, 2, 3 and 4) with probability of detection $p_{11,D}$, target 5 which is of target type 2 with probability of detection $p_{12,D}$, target 9 which is of target type 3 with probability of detection $p_{13,D}$ and target 13 which is of target type 4 with probability of detection $p_{14,D}$. Detector 2 detects target type 2 (targets 5, 6, 7 and 8) with probability of detection $p_{22,D}$, target 1 which is of target type 1 with probability of detection $p_{21,D}$, target 10 which is of target type 3 with probability of detection $p_{23,D}$ and target 14 which is of target type 4 with probability of detection $p_{24,D}$. Similarly, detector 3 detects target type 3 (targets 9, 10, 11 and 12) with probability of detection $p_{33,D}$, target 2 which is of target type 1 with probability of detection $p_{31,D}$, target 6 which is of target type 2 with probability of detection $p_{32,D}$ and target 15 which is of target type 4 with probability of detection $p_{34,D}$. Moreover, detector 4 detects target type 4 (targets 13, 14, 15 and 16) with probability of detection $p_{44,D}$, target 3 which is of target type 1 with probability of detection $p_{41,D}$, target 7 which is of target type 2 with probability of detection $p_{42,D}$ and target 11 which is of target type 3 with probability of detection $p_{43,D}$. This means that targets 1, 2 and 3 from target type 1, targets 5, 6 and 7 from target type 2, targets 8, 9 and 10 from target type 3, and targets 13, 14 and 15 from target type 4 are detected two times. Our main goal is to filter out confused measurements which correspond to a specific target i.e. doubly detected targets are estimated once, not twice. Therefore, the number of targets in the scene is 16 maximum, not 28.

The figures shown in Fig. 3 and Fig. 4 show the comparisons of the outputs of both the quad GM-PHD filter (Fig. 4) and four independent GM-PHD filters (Fig. 3) for detection probability of confusion of 0.6. For both approaches, the simulated ground truths are shown in red for target type 1, black for target type 2, yellow for target type 3 and magenta for target type 4 while the estimates are shown in blue circles for target type 1, green triangles for target type 2, cyan asterisks for target type 3 and black circles for target type 4. Accordingly, for simulated measurements, the quad GM-PHD filter outputs estimates of target type 1 (targets 1, 2, 3 and 4), target type 2 (targets 5, 6, 7 and 8), target type 3 (targets 9, 10, 11 and 12) and target type 4 (targets 13, 14, 15 and 16) being well differentiated as shown in Fig. 4. However, using four independent GM-PHD filters, estimates of targets 1, 2, 3, 4, 5, 9 and 13 are obtained from GM-PHD filter 1, estimates of targets 1, 5, 6, 7, 8, 10 and 14 from GM-PHD filter 2, estimates of targets 2, 9, 10, 11, 12 and 15 from GM-PHD filter 3, and estimates of targets 3, 7, 11, 13, 14, 15 and 16 from GM-PHD filter 4 i.e. targets 1, 2, 3, 5, 6, 7, 9, 10, 11, 13, 14 and 15 are estimated twice though intermittently, as shown in Fig. 3 overlayed. Even if the confusion rates increase to 0.6, the proposed method filters out the target confusions effectively, discriminating the target types. This confusion problem is solved by using our proposed approach with a computation time of 251.27 seconds for 120 iterations when compared to 152.06 seconds using four independent GM-PHD filters implemented on a Core i7 2.30 GHz processor and a 8 GB RAM laptop using MATLAB when setting detection probabilities of confusion to 0.6, for example, as given in Table I. In Fig. 3 and Fig. 4, E stands for the end point whereas the other side of each target is the starting point of the simulation.

When we set the probabilities of confusion to 0.0 i.e. no target confusions, the quad GM-PHD filter performs similar to four GM-PHD filters i.e. it degrades to four GM-PHD filters. However, if the values of confusion detection probabilities are very close to the values of the true detection probabilities, the quad GM-PHD filter still effectively filters the confusion in the detection of targets though it sometimes fails to discriminate the target types. If the probability of confusion is the same as of true detection, then the result is random on first guess (sometimes fails to discriminate the target types) though it still filters out the confusions effectively. Therefore, the values of the confusion probabilities ($p_{12,D}$, $p_{13,D}$, $p_{14,D}$, $p_{21,D}$, $p_{23,D}$, $p_{24,D}$, $p_{31,D}$, $p_{32,D}$, $p_{34,D}$, $p_{41,D}$, $p_{42,D}$ and $p_{43,D}$) should be less than the values of the true detection probabilities ($p_{11,D}$, $p_{22,D}$, $p_{33,D}$ and $p_{44,D}$) to discriminate the target types. However, in real applications (e.g. visual tracking), this does not happen, i.e. the confusion detection probabilities can never become equal in values to the true detection probabilities as object detectors are at least becoming more accurate than random guessing; nobody would employ a random detector.

On the other hand, if each target is regarded as a type (e.g. each of the four target types in this example has only one target), the N-type GM-PHD filter is used as a labeler of each target i.e. it discriminates those targets from frame to frame (from one time step to another) whether or not confusions between targets exist rather than simply degrading to the standard GM-PHD filter(s). Obviously, the number of target types ($N\geq 2$) needs to be known in advance so that measurements from each target type will be collected and then provided to the algorithm.

Furthermore, we assess tracking accuracy using the cardinality (number of targets) and Optimal Subpattern Assignment (OSPA) metric [22] (using order $p=1$ and cutoff $c=100$). From Fig. 5a (when setting the probabilities of confusion to 0.6), we observe that the cardinality of targets estimated using four independent GM-PHD filters (blue) has much more deviation from the ground truth (in red) when compared to the one obtained using our proposed quad GM-PHD filter (green). Similarly, the OSPA error of using four independent GM-PHD filters (blue) is much greater than that of using the quad GM-PHD filter (green) as shown in Fig. 5b. The overall average value of the OSPA error for four independent GM-PHD filters is 46.47m compared to 28.81m when using our proposed quad GM-PHD filter as given in Table I. The OSPA error and time taken (in brackets) when using probabilities of confusion of 0.3 and 0.9 are also given in Table I. As can be observed from Fig. 6 and Table I, as we increase the probabilities of confusions from 0.0 to 0.9, the OSPA error for quad GM-PHD filter is almost constant which shows how efficient the quad GM-PHD filter is in handling target confusions. However, for the case of using four independent GM-PHD filters, the OSPA error increases significantly as we increase the probabilities of confusions from 0.0 to 0.9 which is due to the increase in target confusions. The time taken (given in Table I) also increases slightly for both methods with the increment of target confusions.

VII Conclusion

In this work, we propose a novel, N-type PHD filter where $N\geq 2$, which is an extension of the standard PHD filter in the RFS framework to account for many different types of targets with separate but confused observations of the same scene. In this approach, we assume that there are confusions between detections, i.e. clutter arises not just from background false positives, but also from target confusion. Under the Gaussianity and linearity assumptions, the Gaussian mixture (GM) implementation is proposed for this N-type PHD filter, N-type GM-PHD filter. We evaluate the quad GM-PHD filter as a typical example and compare it to four independent GM-PHD filters, indicating that our approach shows better performance determined using cardinality, OSPA metric and discrimination rate among the different target types. Even though, we show the simulation analysis for $N=4$, in principle the methodology can be applied to $N$ types of targets where $N$ is a variable in which the number of possible confusions may rise as $N(N-1)$. For instance, after experimenting the dual GM-PHD filter ($N=2$) and the tri-GM-PHD filter ($N=3$) by simulation and making sure that they show similar behaviour as for $N=4$, we applied the tri-GM-PHD filter for visual tracking applications in [19]. In case there is no target confusion, the N-type GM-PHD filter performs similar to $N$ independent GM-PHD filters. On the other hand, if each target is regarded as a type, for instance each of the four target types in this given example has only one target, the N-type GM-PHD filter is used as a labeler of each target. This means it discriminates each target from frame to frame (from one time step to another) rather than simply degrading to the standard GM-PHD filter(s).

VIII Appendix A

PHD filter for multi-target tracking of a single target type was proposed by Mahler in [5] where its derivation has been given. We have derived novel extensions for the updated PHDs of a N-type PHD filter from PGFLs of each target type handling confusions among target types starting from the standard proved predicted PHDs in section III. In this case we assume the standard predicted PHD can be applied to each target type as influence of the measurement confusions among the target types we are solving occurs at the update stage of the PHD recursion. To make this paper self-sufficient, we derive the predicted PHDs for each target type $i\in\{1,...,N\}$ starting from the PGFLs of the predicted processes. We use the product rule and chain rule for functional derivatives to find PHD prediction equation for target type $i\in\{1,...,N\}$ (differentiate and set $h=1$) as follows:

[TABLE]

where (taking $\varphi_{x}=\delta_{x}$)

[TABLE]

[TABLE]

and

[TABLE]

the aim is to compute the PHD of the predicted process $D_{i,k|k-1}(x)$ for target type $i\in\{1,...,N\}$.

Deriving Eq. 36 in $\varphi_{x}=\delta_{x}$ and using the product rules yields:

[TABLE]

Now, using the chain rule $A$ reads:

[TABLE]

where $\frac{\delta G_{i,S}(h|.)}{\delta\varphi_{x}}$ is found by deriving Eq. 39 and using the linearity of $h\rightarrow\int h(y)y_{i,k|k-1}(y|.)dy$ as follows:

[TABLE]

Let us write $g(.)=p_{i,S,k}(.)y_{i,k|k-1}(x|.)$ for simplicity’s sake. Now all we need to do is to compute \frac{\delta G_{i,k-1}}{\delta g(.)}(G_{i,S}(h|.))\bigg{|}_{h=1}. Using the Janossy densities we can write:

[TABLE]

Deriving $G_{i,k-1}$ in $g(.)$ gives:

[TABLE]

Using the definition of functional derivatives with $F[h]=h(x_{l})$ gives:

[TABLE]

That is, since the Janossy densities are symmetric functions:

[TABLE]

Or, using the characterization of the PHD $\mathcal{D}_{i,k-1}(x)$ as the first-moment density of the process:

[TABLE]

Combining Eq. 40, Eq. 41, Eq. 42 and Eq. 47 yields the result:

[TABLE]

Acknowledgment

We would like to acknowledge the support of the Engineering and Physical Sciences Research Council (EPSRC), grant references EP/K009931 and a James Watt Scholarship. We would also like to thank Dr. Daniel Clark for sharing his expertise and understanding of RFS methodology.