# Multitarget Tracking with Spatial Nonmaximum Suppressed Sensor Selection

Authors: Liang Ma; Kai Xue; Ping Wang
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

Multitarget tracking is one of the most important applications of sensor networks, yet it is an extremely challenging problem since multisensor multitarget tracking itself is nontrivial and the difficulty is further compounded by sensor management. Recently, random finite set based Bayesian framework has opened doors for multitarget tracking with sensor management, which is modelled in the framework of partially observed Markov decision process (POMDP). However, sensor management posed as a POMDP is in essence a combinatorial optimization problem which is NP-hard and computationally unacceptable. In this paper, we propose a novel sensor selection method for multitarget tracking. We first present the sequential multi-Bernoulli filter as a centralized multisensor fusion scheme for multitarget tracking. In order to perform sensor selection, we define the hypothesis information gain (HIG) of a sensor to measure its information quantity when the sensor is selected alone. Then, we propose spatial nonmaximum suppression approach to select sensors with respect to their locations and HIGs. Two distinguished implementations have been provided using the greedy spatial nonmaximum suppression. Simulation results verify the effectiveness of proposed sensor selection approach for multitarget tracking.


---

## Body

## 1. Introduction

With recent advances in microelectromechanical systems, various kinds of sensors with strong communication ability and accurate data allocation have been manufactured at surprisingly low cost. Different applications of the sensor network have received increasing research interest, such as environmental monitoring, target tracking, and event detection [1, 2]. A static sensor network is a special type of sensor networks, which is composed of densely distributed sensors with fixed and known locations. In the static sensor network, sensor selection is of crucial importance for applications to fulfill a specific task in the optimal way, especially for multitarget tracking [3]. However, very little progress has been made in this area since multisensor multitarget tracking itself is nontrivial and the difficulty is further compounded by sensor selection [4].In the literature, there are extensive studies on the sensor selection (also named as sensor management) problem for target tracking. For the single target tracking case, [5] adopted the decentralized posterior Cramr-Rao lower bound to measure tracking accuracy and used an iterative local search technique for sensor selection. [3] used convex relaxation for sensor selection and then adopted the sequential Kalman filter to track each target in a distributed manner. Recently, random finite set (RFS) based Bayesian framework has provided an elegant solution for multitarget tracking with sensor management using mathematical tools from the finite set statistics (FISST) [6]. In [4, 7], sensor management has been modelled by a partially observed Markov decision process (POMDP) given the multitarget state, which has been shown effective in single sensor control and single sensor selection [8, 9]. However, multisensor management modelled as a POMDP is in essence a combinatorial optimization problem which is NP-hard and computationally unacceptable for a sensor network.This paper considers sensor selection problem for multitarget tracking from a new perspective:which sensors should not be selected. It is intuitive that we should only select informative sensors in use and ignore noninformative sensors in order to alleviate data transmission bandwidth in the sensor network. In this paper, we show that this perspective offers a better alternative to exploit the sensor selection problem which is effective and easy to implement. We first present a centralized multisensor fusion scheme in RFS based Bayesian framework for multitarget tracking and then propose the sequential multi-Bernoulli filter as a feasible approximation of multisensor multitarget state estimation. With respect to sensor selection, we define thehypothesis information gain (HIG) of a sensor to measure its information quantity when it is selected alone. Two types of HIG are presented using the Rényi divergence between the prior and posterior target density, respectively, for individual target state and multitarget state. Then, we propose thespatial nonmaximum suppression (SNMS) approach to suppress sensors with low HIG and subsequently select sensors with high HIG. Two sensor selection approaches via SNMS have been provided in the sequential Monte Carlo (SMC) implementation corresponding to two defined HIGs. Numerical simulations have verified the capability of the proposed sensor selection approach for multitarget tracking.The reminder of this paper is organized as follows. Section2 presents a general description of the sensor selection problem for multitarget tracking to lay a solid foundation. In Section 3, we illustrate the RFS based Bayesian framework, the Cardinality-Balanced multi-Bernoulli filter and then provide the sequential multi-Bernoulli filter as a centralized multisensor fusion scheme. The sensor selection approach via SNMS is discussed in Section 4. Section 5 provides numerical results, followed by the conclusions in Section 6.


## 2. Problem Formulation

Sensor selection for multitarget tracking entails a scenario where there are a large number of densely distributed sensors with fixed and known positions. Each sensor has limited field of view (FOV) but a relatively large communication range. The data allocated by each sensor are transmitted to a fusion center to perform centralized information fusion in order to track all targets in the surveillance area. Sensor selection is required to balance between tracking accuracy and the network workload which is highly dependent on the number of activated sensors. In this paper, we consider two-dimensional coordinate tracking as a particular interest that is demonstrated in Figure1.Figure 1
Multitarget tracking in a static sensor network. Yellow region is the FOV of a sensor. Only five of all sensors can detect targets. Target 1 and target 2 can be detected by two sensors, while target 3 is only detected by one sensor. The fusion center is not shown here but exists.As illustrated in Figure1, multitarget tracking in a sensor network requires multisensor fusion for multitarget state estimation and appropriate sensor selection approach to guarantee tracking accuracy with minimum number of activated sensors. Assume that target moves according to the nearly constant velocity model given by(1)
x
k
=
F
x
k
-
1
+
G
v
k
,where x
k
=
[
p
x
,
k
,
v
x
,
k
,
p
y
,
k
,
v
y
,
k
]
T and p
x
,
k, p
y
,
k are planar position and v
x
,
k, v
y
,
k are planar velocity, respectively, along x-, y-coordinate. Also F
=
I
2
⊗
1
Δ
0
1; G
=
I
2
⊗
Δ
2
/
2
Δ. I
2 is 2
×
2 identity matrix and ⊗ denotes Kronecker product. Δ is the sampling period, and v
k
~
N
(
0
,
Q
k
) is a 2
×
1 i.i.d. Gaussian noise. Assume process noise is time invariant and identical for both v
x
,
k and v
y
,
k; then Q
=
σ
v
2
I
2, where σ
v is the standard deviation.The observation of sensorl originated from target with state x
k is a noisy vector z
k
l, and the measurement model is given by(2)
z
k
l
=
h
l
x
k
+
w
l
,where w
l is zero mean Gaussian noise N
(
w
l
;
0
,
R
l
) and h
l
(
·
) is dependent on the position of the lth sensor s
l
=
[
x
l
y
l
]
T and the type of sensor l.


## 3. Bayesian Multitarget Filtering

This section provides the basic concepts and notations of RFS based Bayesian framework and presents the multisensor information fusion approach thereafter. Section3.1 gives a general description of the RFS based multisensor multitarget Bayesian filtering. Due to the fact that the multisensor multitarget Bayesian filtering is intractable and computationally unacceptable, we provide the sequential multi-Bernoulli filter as a feasible approach. For clarity, we first introduce the multi-Bernoulli filter in Section 3.2 and then provide the SMC implementation of the sequential multi-Bernoulli filter in Section 3.3.


### 3.1. Multisensor Multitarget Bayesian Framework

Stochastic filtering in Bayesian framework has been developed for decades [10]. Under the assumption of linear model and Gaussian distribution, Kalman filter was first derived in [11] and has been widely used for target tracking since then. With respect to the multisource multitarget tracking case, RFS based Bayesian framework provides an elegant solution that outweighs the data association approach significantly.A random finite set is a random variable that takes value as an unordered finite set. The randomness of an RFS refers to two aspects: the set cardinality (the number of elements in the set) is random; each element in the set is also a random variable. The probabilistic description of RFS has been studied regarding different types of probability distributions such as multi-Bernoulli (or Bernoulli) RFS, i.i.d. (short for independent identically distributed) cluster RFS, and Poisson RFS [12]. Let X
k and Z
k, respectively, denote the multitarget state set and the observation set of multiple sensors,(3)
X
k
=
x
k
,
1
,
…
,
x
k
,
N
k
,
Z
k
=
Z
k
1
,
…
,
Z
k
m
,where m is the total number of sensors and Z
k
l
=
z
k
l
,
…
,
z
k
,
M
k
,
l
l for l
=
1
,
…
,
m. N
k is the time-varying cardinality of targets while M
k
,
l is the cardinality of the measurement set generated by sensor l.Using the RFS representation, the movement of multiobject can be described using two parts: an RFS for survival targets from previous time stepS
k and an RFS for spontaneous birth targets at current time Γ
k. Thus, at time k we have the predicted RFS X
k
=
S
k
∪
Γ
k. The RFS for measurements Z
k
l of sensor l can be represented as a union of two parts: target-generated measurements Θ
k
l and clutter K
k
l; thus, Z
k
l
=
Θ
k
l
∪
K
k
l.Given any specific type of RFS, the Bayesian framework for optimal estimation via RFS is in the same form as the classical Bayesian filtering given as follows:(4)
f
k
∣
k
-
1
X
∣
Z
1
:
k
-
1
=
∫
‍
f
k
∣
k
-
1
X
∣
X
′
f
k
-
1
X
′
∣
Z
1
:
k
-
1
δ
X
′
,
(5)
f
k
X
Z
1
:
k
=
f
k
Z
k
X
f
k
∣
k
-
1
X
Z
1
:
k
-
1
∫
‍
f
k
Z
k
X
f
k
∣
k
-
1
X
Z
1
:
k
-
1
δ
X
,which represent the prediction and update process of Bayesian recursion via set integrals, respectively. Notice that (4) and (5) are for multisensor multitarget Bayesian filtering and computationally intractable. Meanwhile, the single senor multitarget filtering version given in [13] is a particular case of it, from which the probability hypothesis density (PHD) [14], cardinalized PHD (CPHD) [15], and multi-Bernoulli filter [13] have been derived under different forms of RFS.


### 3.2. Cardinality-Balanced Multi-Bernoulli Filter

Here, we introduce the multi-Bernoulli RFS for multitarget state modelling, which offers a better alternative than the Poisson and i.i.d. cluster RFS in applications with highly nonlinear model and/or nonhomogeneous sensor type [13]. Assume the dimension of the target state is n; then the target state space is denoted by X
⊆
R
n. A multi-Bernoulli RFS X on X is a union of a fixed number of independent Bernoulli RFSs X
(
j
) with existence probability r
(
j
)
∈
(
0,1
) and probability density p
(
j
) (defined on X), j
=
1
,
…
,
M; that is X
=
∪
j
=
1
M
{
X
(
j
)
}.Use a Bernoulli set for modelling the state of a single target; then the multitarget state can be modeled as multi-Bernoulli RFSΞ with probability density given in [13] as follows:(6)
f
x
1
,
…
,
x
n
=
f
∅
·
∑
‍
r
j
,
p
j
∈
Ξ
∏
i
=
1
n
r
j
p
j
x
i
1
-
r
j
,where r
(
j
) and p
(
j
), respectively, represent the existence probability and distribution of the jth target and f
(
∅
)
=
∏
j
=
1
M
‍
(
1
-
r
(
j
)
). It is clear that the multitarget density can be completely specified by multi-Bernoulli parameter set {
(
r
(
j
)
,
p
(
j
)
)
}
j
=
1
M. Hence, let us denote the multitarget density at time k as π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short in the following content. In multi-Bernoulli filter, the probability hypothesis density (also known as the intensity function in [13, 16]), as the first-order moment, is propagated over time as approximations of the full posteriors f
(
X
∣
Z
) [6]. In the following, we refer to “probability hypothesis density” as “density” for short. Let us denote the multitarget posterior density at time k using multi-Bernoulli parameters by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short. The following gives the recursion of the Cardinality-Balanced multi-Bernoulli filter.Prediction. At time k, if the posterior multitarget density is multi-Bernoulli given by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k and the density of new births is also multi-Bernoulli π
Γ
,
k
+
1
=
{
(
r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
)
}
j
=
1
M
Γ
,
k
+
1, then the predicted density is given by(7)
π
k
+
1
∣
k
=
r
k
+
1
∣
k
j
,
p
k
+
1
∣
k
j
j
=
1
M
k
∪
r
Γ
,
k
+
1
j
,
p
Γ
,
k
+
1
j
j
=
1
M
Γ
,
k
+
1
,where for survival targets(8)
r
k
+
1
∣
k
j
=
r
k
j
·
p
k
j
,
p
S
,
k
,
p
k
+
1
∣
k
j
x
=
f
k
+
1
∣
k
x
∣
·
,
p
k
j
p
S
,
k
p
k
j
,
p
S
,
k
,and, for new born targets, r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
(
x
) are prior existence probability and distribution of birth model.Update. At time k
+
1, if the predicted multitarget density is multi-Bernoulli π
k
+
1
∣
k
=
{
(
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
)
}
j
=
1
M
k
+
1
∣
k, the output of corrector is composed of legacy tracks and measurement-updated tracks,(9)
π
k
+
1
=
r
L
,
k
+
1
j
,
p
L
,
k
+
1
j
j
=
1
M
k
+
1
∣
k
∪
r
U
,
k
+
1
z
,
p
U
,
k
+
1
·
∣
z
z
∈
Z
k
,where(10)
r
L
,
k
+
1
j
=
r
k
+
1
∣
k
j
1
-
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
,
p
L
,
k
+
1
j
x
=
p
k
+
1
∣
k
j
x
1
-
p
D
,
k
+
1
x
1
-
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
,
r
U
,
k
+
1
z
=
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
2
×
κ
k
+
1
z
+
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
-
1
,
p
U
,
k
+
1
z
=
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
x
ψ
k
+
1
,
z
x
×
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
-
1
,
ψ
k
+
1
,
z
x
=
g
k
+
1
z
x
p
D
,
k
+
1
x
.
p
S and p
D are probability of survival and detection for target. The inner product ·
,
· is defined between two real valued functions β and γ by β
,
γ
=
∫
‍
β
(
x
)
γ
(
x
)
d
x. Note that without loss of generality we refer to the Cardinality-Balanced multi-Bernoulli filter as “multi-Bernoulli” filter for simplicity in this paper.


### 3.3. Sequential Multi-Bernoulli Filter

Sequential update has been widely used and verified to be a good approximation for information fusion of multiple sensors [3, 17]. Here, we introduce the sequential multi-Bernoulli filter in the SMC implementation which is proposed in [17].Suppose that, at timek, the posterior multitarget density is given as {
r
k
(
j
)
,
p
k
(
j
)
}
j
=
1
M
k, and the distribution of each target is given by a set of weighted particles p
k
(
j
)
(
x
)
=
∑
i
=
1
L
k
(
j
)
‍
ω
i
,
k
(
j
)
·
δ
x
i
,
k
(
j
)
(
x
). Then, the SMC implementation of the sequential multisensor multi-Bernoulli filter is provided in Algorithm 1. We refer the readers to subsection IV-A of [13] for detailed equations.Algorithm 1:SMC sequential multisensor multi-Bernoulli filter.
Prediction:
Input:  π
k
=
{
r
k
(
j
)
,
p
k
(
j
)
}
j
=
1
M
k
(1) for  Survival targets: j
=
1
,
…
,
M
k  do
(2)  for i
=
1
,
…
,
L
k
(
j
), sample x
S
,
i
,
k
+
1
|
k
(
j
), compute weight ω
S
,
i
,
k
+
1
|
k
(
j
) and normalization ω
~
S
,
i
,
k
+
1
|
k
(
j
)
(3)  compute r
S
,
k
+
1
|
k
(
j
)
,
p
S
,
k
+
1
∣
k
(
j
)
(
x
)
(4) end for
(5) for newborn targets: j
=
1
,
…
,
M
Γ
,
k
+
1  do
(6)  for i
=
1
,
…
,
L
Γ
,
k
(
j
), sample x
Γ
,
i
,
k
+
1
(
j
), compute weight ω
Γ
,
i
,
k
+
1
(
j
) and normalization ω
~
Γ
,
i
,
k
+
1
(
j
)
(7)  compute r
Γ
,
k
+
1
(
j
), p
Γ
,
k
+
1
(
j
)
(
x
)
(8) end for
Output:  π
k
+
1
|
k
=
{
(
r
k
+
1
|
k
(
j
)
,
p
k
+
1
|
k
(
j
)
)
}
j
=
1
M
k
∪
{
(
r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
)
}
j
=
1
M
Γ
,
k
+
1
Update:
Input:  π
k
+
1
|
k
=
{
r
k
+
1
|
k
(
j
)
,
p
k
+
1
|
k
(
j
)
}
j
=
1
M
k
+
1
|
k, candidate sensors located at s
l with observation set Z
k
+
1
l
(
l
=
1
,
…
,
m
′
)
(1) for  l
=
1
,
…
,
m
′  do
(2)  for  Legacy targets: j
=
1
,
…
,
M
k
+
1
|
k  do
(3)   compute ω
L
,
i
,
k
+
1
(
j
)
,
l and normalization ω
~
L
,
i
,
k
+
1
(
j
)
,
l, compute pseudo-likelihood ϱ
L
,
k
+
1
(
j
)
,
l
(4)   compute r
L
,
k
+
1
(
j
)
,
l, p
L
,
k
+
1
(
j
)
,
l
(
x
)
(5)  end for
(6)  for  z
∈
Z
k
+
1
l  do
(7)   for  Measurement-updated targets: j
=
1
,
…
,
M
k
+
1
|
k  do
(8)    compute ω
U
,
i
,
k
+
1
(
j
)
,
l
(
z
) and normalization ω
~
U
,
i
,
k
+
1
(
j
)
,
l
(
z
), compute pseudo-likelihood ϱ
U
,
k
+
1
(
j
)
,
l
(
z
)
(9)    compute r
U
,
k
+
1
(
j
)
,
l
(
z
), p
U
,
k
+
1
(
j
)
,
l
(
x
;
z
)
(10)  end for
(11)  end for
(12)  π
k
+
1
|
k
=
{
(
r
L
,
k
+
1
(
j
)
,
l
,
p
L
,
k
+
1
(
j
)
,
l
)
}
j
=
1
M
k
+
1
|
k
∪
{
(
r
U
,
k
+
1
l
(
z
)
,
p
U
,
k
+
1
l
(
⋅
|
z
)
)
}
z
∈
Z
k
+
1
l
(13) end for
Output:  π
k
+
1
=
{
(
r
L
,
k
+
1
(
j
)
,
m
′
,
p
L
,
k
+
1
(
j
)
,
m
′
)
}
j
=
1
M
k
+
1
|
k
∪
{
(
r
U
,
k
+
1
m
′
(
z
)
,
p
U
,
k
+
1
m
′
(
⋅
|
z
)
)
}
z
∈
Z
k
+
1
m
′The superscript(
j
)
,
l in Algorithm 1 represents the predicted jth Bernoulli set updated with the lth sensor. To avoid the infinite growth of multi-Bernoulli set number, those with existence probability less than a predefined threshold (e.g., 0.001) are removed. Meanwhile, the particle number is limited between L
m
i
n and L
m
a
x, in case that sampling is not enough or resampling reallocates too many particles. The number of particles for each Bernoulli set is proportional to each target existence r
k
(
j
) during the resampling step. With a given existence threshold 0.75, those sets with r
k
(
j
) over 0.75 are true tracks while the others are not.


## 3.1. Multisensor Multitarget Bayesian Framework

Stochastic filtering in Bayesian framework has been developed for decades [10]. Under the assumption of linear model and Gaussian distribution, Kalman filter was first derived in [11] and has been widely used for target tracking since then. With respect to the multisource multitarget tracking case, RFS based Bayesian framework provides an elegant solution that outweighs the data association approach significantly.A random finite set is a random variable that takes value as an unordered finite set. The randomness of an RFS refers to two aspects: the set cardinality (the number of elements in the set) is random; each element in the set is also a random variable. The probabilistic description of RFS has been studied regarding different types of probability distributions such as multi-Bernoulli (or Bernoulli) RFS, i.i.d. (short for independent identically distributed) cluster RFS, and Poisson RFS [12]. Let X
k and Z
k, respectively, denote the multitarget state set and the observation set of multiple sensors,(3)
X
k
=
x
k
,
1
,
…
,
x
k
,
N
k
,
Z
k
=
Z
k
1
,
…
,
Z
k
m
,where m is the total number of sensors and Z
k
l
=
z
k
l
,
…
,
z
k
,
M
k
,
l
l for l
=
1
,
…
,
m. N
k is the time-varying cardinality of targets while M
k
,
l is the cardinality of the measurement set generated by sensor l.Using the RFS representation, the movement of multiobject can be described using two parts: an RFS for survival targets from previous time stepS
k and an RFS for spontaneous birth targets at current time Γ
k. Thus, at time k we have the predicted RFS X
k
=
S
k
∪
Γ
k. The RFS for measurements Z
k
l of sensor l can be represented as a union of two parts: target-generated measurements Θ
k
l and clutter K
k
l; thus, Z
k
l
=
Θ
k
l
∪
K
k
l.Given any specific type of RFS, the Bayesian framework for optimal estimation via RFS is in the same form as the classical Bayesian filtering given as follows:(4)
f
k
∣
k
-
1
X
∣
Z
1
:
k
-
1
=
∫
‍
f
k
∣
k
-
1
X
∣
X
′
f
k
-
1
X
′
∣
Z
1
:
k
-
1
δ
X
′
,
(5)
f
k
X
Z
1
:
k
=
f
k
Z
k
X
f
k
∣
k
-
1
X
Z
1
:
k
-
1
∫
‍
f
k
Z
k
X
f
k
∣
k
-
1
X
Z
1
:
k
-
1
δ
X
,which represent the prediction and update process of Bayesian recursion via set integrals, respectively. Notice that (4) and (5) are for multisensor multitarget Bayesian filtering and computationally intractable. Meanwhile, the single senor multitarget filtering version given in [13] is a particular case of it, from which the probability hypothesis density (PHD) [14], cardinalized PHD (CPHD) [15], and multi-Bernoulli filter [13] have been derived under different forms of RFS.


## 3.2. Cardinality-Balanced Multi-Bernoulli Filter

Here, we introduce the multi-Bernoulli RFS for multitarget state modelling, which offers a better alternative than the Poisson and i.i.d. cluster RFS in applications with highly nonlinear model and/or nonhomogeneous sensor type [13]. Assume the dimension of the target state is n; then the target state space is denoted by X
⊆
R
n. A multi-Bernoulli RFS X on X is a union of a fixed number of independent Bernoulli RFSs X
(
j
) with existence probability r
(
j
)
∈
(
0,1
) and probability density p
(
j
) (defined on X), j
=
1
,
…
,
M; that is X
=
∪
j
=
1
M
{
X
(
j
)
}.Use a Bernoulli set for modelling the state of a single target; then the multitarget state can be modeled as multi-Bernoulli RFSΞ with probability density given in [13] as follows:(6)
f
x
1
,
…
,
x
n
=
f
∅
·
∑
‍
r
j
,
p
j
∈
Ξ
∏
i
=
1
n
r
j
p
j
x
i
1
-
r
j
,where r
(
j
) and p
(
j
), respectively, represent the existence probability and distribution of the jth target and f
(
∅
)
=
∏
j
=
1
M
‍
(
1
-
r
(
j
)
). It is clear that the multitarget density can be completely specified by multi-Bernoulli parameter set {
(
r
(
j
)
,
p
(
j
)
)
}
j
=
1
M. Hence, let us denote the multitarget density at time k as π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short in the following content. In multi-Bernoulli filter, the probability hypothesis density (also known as the intensity function in [13, 16]), as the first-order moment, is propagated over time as approximations of the full posteriors f
(
X
∣
Z
) [6]. In the following, we refer to “probability hypothesis density” as “density” for short. Let us denote the multitarget posterior density at time k using multi-Bernoulli parameters by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short. The following gives the recursion of the Cardinality-Balanced multi-Bernoulli filter.Prediction. At time k, if the posterior multitarget density is multi-Bernoulli given by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k and the density of new births is also multi-Bernoulli π
Γ
,
k
+
1
=
{
(
r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
)
}
j
=
1
M
Γ
,
k
+
1, then the predicted density is given by(7)
π
k
+
1
∣
k
=
r
k
+
1
∣
k
j
,
p
k
+
1
∣
k
j
j
=
1
M
k
∪
r
Γ
,
k
+
1
j
,
p
Γ
,
k
+
1
j
j
=
1
M
Γ
,
k
+
1
,where for survival targets(8)
r
k
+
1
∣
k
j
=
r
k
j
·
p
k
j
,
p
S
,
k
,
p
k
+
1
∣
k
j
x
=
f
k
+
1
∣
k
x
∣
·
,
p
k
j
p
S
,
k
p
k
j
,
p
S
,
k
,and, for new born targets, r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
(
x
) are prior existence probability and distribution of birth model.Update. At time k
+
1, if the predicted multitarget density is multi-Bernoulli π
k
+
1
∣
k
=
{
(
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
)
}
j
=
1
M
k
+
1
∣
k, the output of corrector is composed of legacy tracks and measurement-updated tracks,(9)
π
k
+
1
=
r
L
,
k
+
1
j
,
p
L
,
k
+
1
j
j
=
1
M
k
+
1
∣
k
∪
r
U
,
k
+
1
z
,
p
U
,
k
+
1
·
∣
z
z
∈
Z
k
,where(10)
r
L
,
k
+
1
j
=
r
k
+
1
∣
k
j
1
-
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
,
p
L
,
k
+
1
j
x
=
p
k
+
1
∣
k
j
x
1
-
p
D
,
k
+
1
x
1
-
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
,
r
U
,
k
+
1
z
=
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
2
×
κ
k
+
1
z
+
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
p
D
,
k
+
1
-
1
,
p
U
,
k
+
1
z
=
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
x
ψ
k
+
1
,
z
x
×
∑
j
=
1
M
k
+
1
∣
k
r
k
+
1
∣
k
j
1
-
r
k
+
1
∣
k
j
p
k
+
1
∣
k
j
,
ψ
k
+
1
,
z
-
1
,
ψ
k
+
1
,
z
x
=
g
k
+
1
z
x
p
D
,
k
+
1
x
.
p
S and p
D are probability of survival and detection for target. The inner product ·
,
· is defined between two real valued functions β and γ by β
,
γ
=
∫
‍
β
(
x
)
γ
(
x
)
d
x. Note that without loss of generality we refer to the Cardinality-Balanced multi-Bernoulli filter as “multi-Bernoulli” filter for simplicity in this paper.


## 3.3. Sequential Multi-Bernoulli Filter

Sequential update has been widely used and verified to be a good approximation for information fusion of multiple sensors [3, 17]. Here, we introduce the sequential multi-Bernoulli filter in the SMC implementation which is proposed in [17].Suppose that, at timek, the posterior multitarget density is given as {
r
k
(
j
)
,
p
k
(
j
)
}
j
=
1
M
k, and the distribution of each target is given by a set of weighted particles p
k
(
j
)
(
x
)
=
∑
i
=
1
L
k
(
j
)
‍
ω
i
,
k
(
j
)
·
δ
x
i
,
k
(
j
)
(
x
). Then, the SMC implementation of the sequential multisensor multi-Bernoulli filter is provided in Algorithm 1. We refer the readers to subsection IV-A of [13] for detailed equations.Algorithm 1:SMC sequential multisensor multi-Bernoulli filter.
Prediction:
Input:  π
k
=
{
r
k
(
j
)
,
p
k
(
j
)
}
j
=
1
M
k
(1) for  Survival targets: j
=
1
,
…
,
M
k  do
(2)  for i
=
1
,
…
,
L
k
(
j
), sample x
S
,
i
,
k
+
1
|
k
(
j
), compute weight ω
S
,
i
,
k
+
1
|
k
(
j
) and normalization ω
~
S
,
i
,
k
+
1
|
k
(
j
)
(3)  compute r
S
,
k
+
1
|
k
(
j
)
,
p
S
,
k
+
1
∣
k
(
j
)
(
x
)
(4) end for
(5) for newborn targets: j
=
1
,
…
,
M
Γ
,
k
+
1  do
(6)  for i
=
1
,
…
,
L
Γ
,
k
(
j
), sample x
Γ
,
i
,
k
+
1
(
j
), compute weight ω
Γ
,
i
,
k
+
1
(
j
) and normalization ω
~
Γ
,
i
,
k
+
1
(
j
)
(7)  compute r
Γ
,
k
+
1
(
j
), p
Γ
,
k
+
1
(
j
)
(
x
)
(8) end for
Output:  π
k
+
1
|
k
=
{
(
r
k
+
1
|
k
(
j
)
,
p
k
+
1
|
k
(
j
)
)
}
j
=
1
M
k
∪
{
(
r
Γ
,
k
+
1
(
j
)
,
p
Γ
,
k
+
1
(
j
)
)
}
j
=
1
M
Γ
,
k
+
1
Update:
Input:  π
k
+
1
|
k
=
{
r
k
+
1
|
k
(
j
)
,
p
k
+
1
|
k
(
j
)
}
j
=
1
M
k
+
1
|
k, candidate sensors located at s
l with observation set Z
k
+
1
l
(
l
=
1
,
…
,
m
′
)
(1) for  l
=
1
,
…
,
m
′  do
(2)  for  Legacy targets: j
=
1
,
…
,
M
k
+
1
|
k  do
(3)   compute ω
L
,
i
,
k
+
1
(
j
)
,
l and normalization ω
~
L
,
i
,
k
+
1
(
j
)
,
l, compute pseudo-likelihood ϱ
L
,
k
+
1
(
j
)
,
l
(4)   compute r
L
,
k
+
1
(
j
)
,
l, p
L
,
k
+
1
(
j
)
,
l
(
x
)
(5)  end for
(6)  for  z
∈
Z
k
+
1
l  do
(7)   for  Measurement-updated targets: j
=
1
,
…
,
M
k
+
1
|
k  do
(8)    compute ω
U
,
i
,
k
+
1
(
j
)
,
l
(
z
) and normalization ω
~
U
,
i
,
k
+
1
(
j
)
,
l
(
z
), compute pseudo-likelihood ϱ
U
,
k
+
1
(
j
)
,
l
(
z
)
(9)    compute r
U
,
k
+
1
(
j
)
,
l
(
z
), p
U
,
k
+
1
(
j
)
,
l
(
x
;
z
)
(10)  end for
(11)  end for
(12)  π
k
+
1
|
k
=
{
(
r
L
,
k
+
1
(
j
)
,
l
,
p
L
,
k
+
1
(
j
)
,
l
)
}
j
=
1
M
k
+
1
|
k
∪
{
(
r
U
,
k
+
1
l
(
z
)
,
p
U
,
k
+
1
l
(
⋅
|
z
)
)
}
z
∈
Z
k
+
1
l
(13) end for
Output:  π
k
+
1
=
{
(
r
L
,
k
+
1
(
j
)
,
m
′
,
p
L
,
k
+
1
(
j
)
,
m
′
)
}
j
=
1
M
k
+
1
|
k
∪
{
(
r
U
,
k
+
1
m
′
(
z
)
,
p
U
,
k
+
1
m
′
(
⋅
|
z
)
)
}
z
∈
Z
k
+
1
m
′The superscript(
j
)
,
l in Algorithm 1 represents the predicted jth Bernoulli set updated with the lth sensor. To avoid the infinite growth of multi-Bernoulli set number, those with existence probability less than a predefined threshold (e.g., 0.001) are removed. Meanwhile, the particle number is limited between L
m
i
n and L
m
a
x, in case that sampling is not enough or resampling reallocates too many particles. The number of particles for each Bernoulli set is proportional to each target existence r
k
(
j
) during the resampling step. With a given existence threshold 0.75, those sets with r
k
(
j
) over 0.75 are true tracks while the others are not.


## 4. Sensor Selection with Spatial Nonmaximum Suppression

In this section, we first illustrate the definition of hypothesis information gain in Section4.1. Then, Section 4.2 provides the description of the spatial nonmaximum suppression. Then, we propose a greedy implementation of sensor selection in Section 4.3 given the proposed spatial nonmaximum suppression.


### 4.1. Hypothesis Information Gain

Before applying nonmaximum suppression, we need to define a proper measure of each sensor to determine whether the sensor is informative for target state estimation. Here, we first introduce the HIG for individual target and then give the HIG for multitarget directly. The HIG for individual target is defined with respect to target density, which is the product of existence probability and state distribution given in [6] as follows:(11)
π
k
j
=
r
k
j
·
p
k
j
,for Bernoulli set {
r
k
(
j
)
,
p
k
(
j
)
}. Notice that the density is the product of a scalar and a probability distribution. The Rényi divergence is adopted to measure the information gain between target prior and posterior distribution. Hence, the HIG of sensor l with respect to target j is defined as(12)
H
k
+
1
j
,
l
=
r
k
+
1
j
,
l
r
k
+
1
∣
k
j
E
R
l
p
k
+
1
j
,
l
∥
p
k
+
1
∣
k
j
,where E
[
·
] is the expectation operator and R
l
(
p
k
+
1
(
j
)
,
l
∥
p
k
+
1
∣
k
(
j
)
) is the Rényi divergence between the prior and posterior distribution of target j denoted by R
(
j
)
,
l for short. Given the Bayesian recursion, R
(
j
)
,
l is given in [8] as follows:(13)
R
j
,
l
=
1
α
-
1
l
o
g
∫
‍
g
k
+
1
l
z
∣
x
k
+
1
j
α
p
k
+
1
∣
k
j
x
d
x
k
+
1
j
p
z
∣
Z
k
α
,where p
(
z
∣
Z
k
)
=
∫
‍
g
k
+
1
l
(
z
∣
x
k
+
1
(
j
)
)
p
k
+
1
∣
k
(
j
)
(
x
)
d
x
k
+
1
(
j
) and α is a parameter that determines how much we emphasize the tails of two distributions in the metric, and the Rényi divergence becomes the Kullback-Leibler discrimination and Hellinger affinity, respectively, when α
→
1 and α
=
0.5 [18].To compute (12), we only generate one future measurement z for sensor l based on the predicted state, assuming no clutter or unity detection rate as illustrated in [9]. Thus, r
k
+
1
(
j
)
,
l
=
r
U
,
k
+
1
(
j
)
,
l
(
z
) for unity detection rate. Assume that at time k the predicted multitarget density π
k
+
1
∣
k
=
{
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
}
j
=
1
M
k
+
1
∣
k is given in SMC form where p
k
+
1
∣
k
j
(
x
)
=
∑
i
=
1
L
k
+
1
∣
k
j
‍
ω
i
,
k
+
1
∣
k
j
·
δ
x
i
,
k
+
1
∣
k
j
x; substitute SMC form of p
k
+
1
∣
k
(
j
)
(
x
) into (13); then, we obtain(14)
H
k
+
1
j
,
l
=
r
U
,
k
+
1
j
,
l
z
α
-
1
r
k
+
1
∣
k
j
×
l
o
g
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
x
k
+
1
j
α
.Notice that the HIG of sensorl described by (14) is defined for each individual target. Consequently, we also define another HIG of sensor l given by(15)
H
k
+
1
l
=
∑
j
=
1
M
k
+
1
∣
k
‍
H
k
+
1
j
,
l
=
∑
j
=
1
M
k
+
1
∣
k
r
U
,
k
+
1
j
,
l
z
α
-
1
r
k
+
1
∣
k
j
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
‍
×
l
o
g
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
,to capture the information gain of sensor l with regard to multitarget state. In this paper, we name H
k
+
1
(
j
)
,
l in (14) as “individual HIG” and H
k
+
1
l in (15) as “sum HIG.” The benefits of using the HIG is twofold: firstly, maximizing the measurement-updated existence probability tends to avoid losing targets; secondly, maximum Rényi divergence between the predicted and updated distribution obtains more information from future measurements and makes target state estimation more accurate.


### 4.2. Spatial Nonmaximum Suppression

Nonmaximum suppression (NMS) plays a very important role in computer vision field especially in the object detection process, which aims to pick real objects with local maxima and suppress those that are outliers. Nonmaximum suppression is first proposed in an edge detection context [19] and then widely used in many detectors, such as points [20, 21], edges [22], and objects [23, 24]. The original version of NMS is one-dimensional (1D) [19] and then extended to isotropic NMS to locate two-dimensional (2D) feature points from an image [20]. In this paper, we use the underlying rationale of the isotropic NMS and propose SNMS method in order to eliminate sensors with low HIG.Since the SNMS method here is different from that in computer vision area, we first present a general description of the SNMS algorithm for spatially distributed valued-points. At timek, assume that there are n valued-points with known positions u
i and values V
k
i for i
=
1
,
…
,
m, given a predefined suppression gate ϵ; then, the ith point is suppressed if there is any point j that satisfies(16)
u
i
-
u
j
≤
ϵ
,
V
k
j
≥
V
k
i
,where · is the Euclidean distance between point i and point j. It is clear that the spatial nonmaximum suppression is trivial and easy to implement in a straightforward way. Algorithm 2 provides the straightforward implementation of proposed SNMS method. It is obvious that in each loop 3 comparisons are made. Hence, the computational complexity of the straightforward implementation is O
(
n
) for the best case and O
(
n
2
) for the worst case and O
(
n
(
n
+
1
)
/
2
) on average. This is because the inner loop in Algorithm 2 may stop at any i for i
=
1
,
…
,
n with equal probability p
(
i
)
=
1
/
n; thus, the expected number of comparisons E
c
m
p is(17)
E
c
m
p
=
3
×
n
×
∑
i
=
1
n
‍
l
·
1
n
=
3
n
n
+
1
2
.Algorithm 2:Spatial nonmaximum suppression (straightforward).
Input:  V
k
i
,
s
i
(
i
=
1
,
…
,
n
), suppression gate ϵ
(1) S
k
=
{
s
1
,
…
,
s
n
}
(2) for  i
=
1
,
…
,
n  do
(3)  for  j
=
1
,
…
,
n  do
(4)   if  i
≠
j and |
|
s
i
-
s
j
|
|
≤
ϵ and V
k
j
≥
V
k
i  then
(5)    remove s
i from  S
k; continue;
(6)   end if
(7)  end for
(8) end for
Output:  S
kTo reduce the computational cost when the number of valued-points is extremely large, we introduce the greedy implementation of SNMS by assuming suppressed points will no longer suppress other points. The greedy implementation has been proven to be as effective as the straightforward implementation and much more efficient in the literature [23, 24]. The greedy SMC implementation of SNMS is given in Algorithm 3. It is clear that the greedy implementation requires sorting the point first by its associated value, and this procedure can be achieved at O
(
n
l
o
g
n
) using proper sorting algorithm.Algorithm 3:Spatial nonmaximum suppression (greedy).
Input:  V
k
i
,
s
i
(
l
=
1
,
…
,
n
), suppression gate ϵ
(1) S
k
=
{
s
1
,
…
,
s
n
}, S
k
s
c
l
t
=
∅
(2) sort S
k given V
k
i in descending order
(3) while  S
k
≠
∅  do
(4)  take the first point s
first from S
k
(5)  for  i
=
1
,
…
,
n  do
(6)   if  |
|
s
i
-
s
first
|
|
≤
ϵ  then
(7)    remove s
i from S
k;
(8)    S
k
s
c
l
t
=
S
k
s
c
l
t
∪
{
s
first
}
(9)   end if
(10) end for
(11) end while
Output:  S
k
s
c
l
t


### 4.3. Implementation of SNMS Sensor Selection

Given the proposed SNMS above, the sensor selection via SNMS is direct given the location and HIG of each sensor. Assume the sensor network containsm candidate sensors with fixed and known position s
l
(
l
=
1
,
…
,
m
), given the multitarget prediction π
k
+
1
∣
k
=
{
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
}
j
=
1
M
k
+
1
∣
k at time k
+
1; we present the SMC implementation of sensor selection via SNMS in Algorithms 4 and 5, respectively, for individual HIG and sum HIG. In the following content, we refer to the sensor selection using individual HIG as “sensor selection I” and using sum HIG as “sensor selection II”. The u
n
i
q
u
e
(
) function in Algorithm 4 eliminates repeatedly chosen sensor to ensure that each sensor will only be updated once in the sequential multi-Bernoulli filter.Algorithm 4:Sensor selection I (individual HIG).
Input: π
k
+
1
|
k
,
s
l
(
l
=
1
,
…
,
m
)
(1) for  j
=
1
,
…
,
M
k
+
1
|
k  do
(2)  for  l
=
1
,
…
,
m  do
(3)   predict z
=
∑
i
=
1
L
k
+
1
|
k
(
j
)
‍
ω
i
,
k
+
1
|
k
(
j
)
h
l
(
x
i
,
k
+
1
|
k
(
j
)
)
(4)   compute H
k
+
1
(
j
)
,
l given by (14)
(5)  end for
(6)  greedy SNMS to obtain S
k
+
1
(
j
)
(7) end for
Output:  S
k
+
1
=
u
n
i
q
u
e
(
∪
j
=
1
M
S
k
+
1
(
j
)
)Algorithm 5:Sensor selection II (sum HIG).
Input:  π
k
+
1
|
k
,
s
l
(
l
=
1
,
…
,
m
)
(1)for
j
=
1
,
…
,
M
k
+
1
|
k  do
(2)for  l
=
1
,
…
,
m  do
(3)   predictz
=
∑
i
=
1
L
k
+
1
|
k
(
j
)
‍
ω
i
,
k
+
1
|
k
(
j
)
h
l
(
x
i
,
k
+
1
|
k
(
j
)
)
(4)   computeH
k
+
1
(
j
)
,
l given by (14)
(5)end for
(6)  computeH
k
+
1
l given by (15)
(7)  greedy SNMS to obtainS
k
+
1
(8)end for
Output:  S
k
+
1In order to track targets accurately, we need to select sensors with higher HIG and eliminate those with lower HIG. The underlying rationale of SNMS sensor selection is that sensors with higher value of HIG are very likely to have some targets near them since each sensor has only limited FOV. Thus, sensors should be selected by seeking the local maxima given their spatial locations. A question may rise here:why not use thresholding method to eliminate sensors with lower HIG. Actually, the thresholding technique is inappropriate in this case for two reasons: firstly, it mayover select sensors that observe the same targets; secondly, sensors that are not so informative but actually useful may be eliminated wrongly.


## 4.1. Hypothesis Information Gain

Before applying nonmaximum suppression, we need to define a proper measure of each sensor to determine whether the sensor is informative for target state estimation. Here, we first introduce the HIG for individual target and then give the HIG for multitarget directly. The HIG for individual target is defined with respect to target density, which is the product of existence probability and state distribution given in [6] as follows:(11)
π
k
j
=
r
k
j
·
p
k
j
,for Bernoulli set {
r
k
(
j
)
,
p
k
(
j
)
}. Notice that the density is the product of a scalar and a probability distribution. The Rényi divergence is adopted to measure the information gain between target prior and posterior distribution. Hence, the HIG of sensor l with respect to target j is defined as(12)
H
k
+
1
j
,
l
=
r
k
+
1
j
,
l
r
k
+
1
∣
k
j
E
R
l
p
k
+
1
j
,
l
∥
p
k
+
1
∣
k
j
,where E
[
·
] is the expectation operator and R
l
(
p
k
+
1
(
j
)
,
l
∥
p
k
+
1
∣
k
(
j
)
) is the Rényi divergence between the prior and posterior distribution of target j denoted by R
(
j
)
,
l for short. Given the Bayesian recursion, R
(
j
)
,
l is given in [8] as follows:(13)
R
j
,
l
=
1
α
-
1
l
o
g
∫
‍
g
k
+
1
l
z
∣
x
k
+
1
j
α
p
k
+
1
∣
k
j
x
d
x
k
+
1
j
p
z
∣
Z
k
α
,where p
(
z
∣
Z
k
)
=
∫
‍
g
k
+
1
l
(
z
∣
x
k
+
1
(
j
)
)
p
k
+
1
∣
k
(
j
)
(
x
)
d
x
k
+
1
(
j
) and α is a parameter that determines how much we emphasize the tails of two distributions in the metric, and the Rényi divergence becomes the Kullback-Leibler discrimination and Hellinger affinity, respectively, when α
→
1 and α
=
0.5 [18].To compute (12), we only generate one future measurement z for sensor l based on the predicted state, assuming no clutter or unity detection rate as illustrated in [9]. Thus, r
k
+
1
(
j
)
,
l
=
r
U
,
k
+
1
(
j
)
,
l
(
z
) for unity detection rate. Assume that at time k the predicted multitarget density π
k
+
1
∣
k
=
{
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
}
j
=
1
M
k
+
1
∣
k is given in SMC form where p
k
+
1
∣
k
j
(
x
)
=
∑
i
=
1
L
k
+
1
∣
k
j
‍
ω
i
,
k
+
1
∣
k
j
·
δ
x
i
,
k
+
1
∣
k
j
x; substitute SMC form of p
k
+
1
∣
k
(
j
)
(
x
) into (13); then, we obtain(14)
H
k
+
1
j
,
l
=
r
U
,
k
+
1
j
,
l
z
α
-
1
r
k
+
1
∣
k
j
×
l
o
g
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
x
k
+
1
j
α
.Notice that the HIG of sensorl described by (14) is defined for each individual target. Consequently, we also define another HIG of sensor l given by(15)
H
k
+
1
l
=
∑
j
=
1
M
k
+
1
∣
k
‍
H
k
+
1
j
,
l
=
∑
j
=
1
M
k
+
1
∣
k
r
U
,
k
+
1
j
,
l
z
α
-
1
r
k
+
1
∣
k
j
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
‍
×
l
o
g
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
∑
i
=
1
L
k
+
1
∣
k
j
ω
i
,
k
+
1
∣
k
j
g
k
+
1
l
z
∣
x
k
+
1
j
α
,to capture the information gain of sensor l with regard to multitarget state. In this paper, we name H
k
+
1
(
j
)
,
l in (14) as “individual HIG” and H
k
+
1
l in (15) as “sum HIG.” The benefits of using the HIG is twofold: firstly, maximizing the measurement-updated existence probability tends to avoid losing targets; secondly, maximum Rényi divergence between the predicted and updated distribution obtains more information from future measurements and makes target state estimation more accurate.


## 4.2. Spatial Nonmaximum Suppression

Nonmaximum suppression (NMS) plays a very important role in computer vision field especially in the object detection process, which aims to pick real objects with local maxima and suppress those that are outliers. Nonmaximum suppression is first proposed in an edge detection context [19] and then widely used in many detectors, such as points [20, 21], edges [22], and objects [23, 24]. The original version of NMS is one-dimensional (1D) [19] and then extended to isotropic NMS to locate two-dimensional (2D) feature points from an image [20]. In this paper, we use the underlying rationale of the isotropic NMS and propose SNMS method in order to eliminate sensors with low HIG.Since the SNMS method here is different from that in computer vision area, we first present a general description of the SNMS algorithm for spatially distributed valued-points. At timek, assume that there are n valued-points with known positions u
i and values V
k
i for i
=
1
,
…
,
m, given a predefined suppression gate ϵ; then, the ith point is suppressed if there is any point j that satisfies(16)
u
i
-
u
j
≤
ϵ
,
V
k
j
≥
V
k
i
,where · is the Euclidean distance between point i and point j. It is clear that the spatial nonmaximum suppression is trivial and easy to implement in a straightforward way. Algorithm 2 provides the straightforward implementation of proposed SNMS method. It is obvious that in each loop 3 comparisons are made. Hence, the computational complexity of the straightforward implementation is O
(
n
) for the best case and O
(
n
2
) for the worst case and O
(
n
(
n
+
1
)
/
2
) on average. This is because the inner loop in Algorithm 2 may stop at any i for i
=
1
,
…
,
n with equal probability p
(
i
)
=
1
/
n; thus, the expected number of comparisons E
c
m
p is(17)
E
c
m
p
=
3
×
n
×
∑
i
=
1
n
‍
l
·
1
n
=
3
n
n
+
1
2
.Algorithm 2:Spatial nonmaximum suppression (straightforward).
Input:  V
k
i
,
s
i
(
i
=
1
,
…
,
n
), suppression gate ϵ
(1) S
k
=
{
s
1
,
…
,
s
n
}
(2) for  i
=
1
,
…
,
n  do
(3)  for  j
=
1
,
…
,
n  do
(4)   if  i
≠
j and |
|
s
i
-
s
j
|
|
≤
ϵ and V
k
j
≥
V
k
i  then
(5)    remove s
i from  S
k; continue;
(6)   end if
(7)  end for
(8) end for
Output:  S
kTo reduce the computational cost when the number of valued-points is extremely large, we introduce the greedy implementation of SNMS by assuming suppressed points will no longer suppress other points. The greedy implementation has been proven to be as effective as the straightforward implementation and much more efficient in the literature [23, 24]. The greedy SMC implementation of SNMS is given in Algorithm 3. It is clear that the greedy implementation requires sorting the point first by its associated value, and this procedure can be achieved at O
(
n
l
o
g
n
) using proper sorting algorithm.Algorithm 3:Spatial nonmaximum suppression (greedy).
Input:  V
k
i
,
s
i
(
l
=
1
,
…
,
n
), suppression gate ϵ
(1) S
k
=
{
s
1
,
…
,
s
n
}, S
k
s
c
l
t
=
∅
(2) sort S
k given V
k
i in descending order
(3) while  S
k
≠
∅  do
(4)  take the first point s
first from S
k
(5)  for  i
=
1
,
…
,
n  do
(6)   if  |
|
s
i
-
s
first
|
|
≤
ϵ  then
(7)    remove s
i from S
k;
(8)    S
k
s
c
l
t
=
S
k
s
c
l
t
∪
{
s
first
}
(9)   end if
(10) end for
(11) end while
Output:  S
k
s
c
l
t


## 4.3. Implementation of SNMS Sensor Selection

Given the proposed SNMS above, the sensor selection via SNMS is direct given the location and HIG of each sensor. Assume the sensor network containsm candidate sensors with fixed and known position s
l
(
l
=
1
,
…
,
m
), given the multitarget prediction π
k
+
1
∣
k
=
{
r
k
+
1
∣
k
(
j
)
,
p
k
+
1
∣
k
(
j
)
}
j
=
1
M
k
+
1
∣
k at time k
+
1; we present the SMC implementation of sensor selection via SNMS in Algorithms 4 and 5, respectively, for individual HIG and sum HIG. In the following content, we refer to the sensor selection using individual HIG as “sensor selection I” and using sum HIG as “sensor selection II”. The u
n
i
q
u
e
(
) function in Algorithm 4 eliminates repeatedly chosen sensor to ensure that each sensor will only be updated once in the sequential multi-Bernoulli filter.Algorithm 4:Sensor selection I (individual HIG).
Input: π
k
+
1
|
k
,
s
l
(
l
=
1
,
…
,
m
)
(1) for  j
=
1
,
…
,
M
k
+
1
|
k  do
(2)  for  l
=
1
,
…
,
m  do
(3)   predict z
=
∑
i
=
1
L
k
+
1
|
k
(
j
)
‍
ω
i
,
k
+
1
|
k
(
j
)
h
l
(
x
i
,
k
+
1
|
k
(
j
)
)
(4)   compute H
k
+
1
(
j
)
,
l given by (14)
(5)  end for
(6)  greedy SNMS to obtain S
k
+
1
(
j
)
(7) end for
Output:  S
k
+
1
=
u
n
i
q
u
e
(
∪
j
=
1
M
S
k
+
1
(
j
)
)Algorithm 5:Sensor selection II (sum HIG).
Input:  π
k
+
1
|
k
,
s
l
(
l
=
1
,
…
,
m
)
(1)for
j
=
1
,
…
,
M
k
+
1
|
k  do
(2)for  l
=
1
,
…
,
m  do
(3)   predictz
=
∑
i
=
1
L
k
+
1
|
k
(
j
)
‍
ω
i
,
k
+
1
|
k
(
j
)
h
l
(
x
i
,
k
+
1
|
k
(
j
)
)
(4)   computeH
k
+
1
(
j
)
,
l given by (14)
(5)end for
(6)  computeH
k
+
1
l given by (15)
(7)  greedy SNMS to obtainS
k
+
1
(8)end for
Output:  S
k
+
1In order to track targets accurately, we need to select sensors with higher HIG and eliminate those with lower HIG. The underlying rationale of SNMS sensor selection is that sensors with higher value of HIG are very likely to have some targets near them since each sensor has only limited FOV. Thus, sensors should be selected by seeking the local maxima given their spatial locations. A question may rise here:why not use thresholding method to eliminate sensors with lower HIG. Actually, the thresholding technique is inappropriate in this case for two reasons: firstly, it mayover select sensors that observe the same targets; secondly, sensors that are not so informative but actually useful may be eliminated wrongly.


## 5. Simulation

In order to demonstrate the performance of proposed sensor selection approach for multitarget tracking, we present a planar multitarget tracking scenario in a static sensor network, in which9
×
9 sensors were laid out uniformly over a square of size [
-
1000
m
,
1000
m
]
×
[
-
1000
m
,
1000
m
] divided into 200
m
×
200
m blocks. There are unknown and time-varying number of targets observed in clutter while range and bearing measurements from each sensor are available.Targets can appear or disappear in the scene at any time, and survival probabilityp
S
=
0.95 for each existing target. Newborn targets appear spontaneously according to γ
k
=
0.2
N
(
·
;
x
¯
,
Q
). Each target moves according to the constant velocity model given by (1). Four targets are presented for tracking purpose as illustrated in Figure 2.Figure 2
Target tracks. Start/stop positions for each track are shown with•
/
■. Blue dots are sensor locations.The tracker, composed of SNMS sensor selection and the sequential multi-Bernoulli filter, runs for 50 scans with sampling periodΔ
=
1 s. The standard derivation of process noise σ
v
=
1
m
/
s for both v
x
,
k and v
y
,
k. The h
l
(
x
k
) of the measurement model described by (2) is given as(18)
h
l
x
k
=
p
k
-
s
l
a
r
c
t
a
n
y
l
-
p
y
,
k
x
l
-
p
x
,
k
,in which p
k
-
s
l
=
(
x
l
-
p
x
,
k
)
2
+
(
y
l
-
p
y
,
k
)
2. The covariance of measurement noise w
k
l for sensor l at time k is R
k
l
=
d
i
a
g
(
[
(
σ
r
,
k
l
)
2
,
(
σ
ϕ
,
k
l
)
2
]
), where(19)
σ
r
,
k
l
2
=
σ
0
+
β
r
p
k
-
s
l
2
,
σ
ϕ
,
k
l
2
=
σ
1
+
β
ϕ
p
k
-
s
l
,with σ
0
=
1 m, β
r
=
5
×
1
0
-
5 m−1, σ
1
=
π
/
180 rad, β
ϕ
=
1
0
-
5 rad·m−1. The FOV of each sensor is [
-
π
/
2
,
+
π
/
2
]
×
[
0,500
m
] with clutter uniformly distributed over this interval. The clutter rate of each sensor is λ
c
=
5 per scan. The probability of detection of sensor l is modelled by (20)
p
D
l
x
k
=
0.99
p
k
-
s
l
≤
R
1
0.99
-
c
p
k
-
s
l
R
1
<
p
k
-
s
l
≤
R
2
0
p
k
-
s
l
>
R
2
,with R
1
=
200 m, R
2
=
500 m, and c
=
5
×
1
0
-
4.Birth parameters for the sequential multi-Bernoulli filter are configured asx
¯
1
=
[
650
m
;
-
10
m
/
s
;
650
m
;
-
30
m
/
s
]
T, x
¯
2
=
[
-
300
m
;
15
m
/
s
;
300
m
;
-
15
m
/
s
]
T, x
¯
3
=
[
-
650
m
;
15
m
/
s
;
-
450
m
;
0
]
T, x
¯
4
=
[
-
300
m
;
30
m
/
s
;
500
m
;
-
2
m
/
s
]
T, respectively, for each target, and Q
=
d
i
a
g
(
[
(
30
m
,
2
m
/
s
,
30
m
,
2
m
/
s
]
2
) is identical for all four targets. L
m
i
n
=
300 and L
m
a
x
=
1000 are the minimum and maximum particle number of each Bernoulli set for track maintenance.We first compare the tracking performance, respectively, using sensor selection I and sensor selection II. The optimal subpattern assignment (OSPA) metric composed of location error and cardinality error is adopted for tracking performance evaluation [25]. Figure 3 shows the OSPA distance (p
=
1, c
=
300 m) comparison from 500 Monte Carlo runs with SNMS gate ϵ
=
500 m. It can be seen that there is only slight difference of the two selection approaches in tracking performance and either sensor selection approach can provide accurate target tracking. The difference occurs near time k
=
10 and k
=
42 when there are target births and deaths, and sensor selection using sum HIG has a relatively large error than that using individual HIG.Figure 3
Comparison of OSPA distance.To further illustrate the SNMS sensor selection and multitarget tracking procedure, we take sensor selection II as an example and present three consecutive frames of one trial run in Figure4. It can be seen that all three figures contain more than four sensors selected while there are only four real targets even when targets are well separated. This is because in the sequential multi-Bernoulli filter there are birth sets for new born targets and legacy sets for temporally undetected targets, and the SNMS sensor selection will treat all these sets as potential targets and compute the sensor HIG with respect to them. Thus, there are sensors selected to detect potential targets, which guarantee that newborn targets can be recognized and temporally missing target can be picked up again. For example, the sequential multi-Bernoulli filter lost target 3 in Figure 4(b) and picked it up right away in Figure 4(c). Trial runs of both selection approaches have been recorded in videos attached as supplement materials (see Supplementary Material available online at http://dx.doi.org/10.1155/2015/148081) for demonstration.Figure 4
Consecutive tracking frames. Blue dots stand for sensors that cannot detect any targets, and red⋄ are sensors that can detect at least one target. Blue line and ∘ are used for target track and current position, black □ for selected sensor, and red ∗ for target estimation.
(a)
(b)
(c)We also study the impact the SNMS gateϵ has on the sensor selection procedure. We set four different SNMS gates ϵ
i
=
500,1000,1500
,
i
n
f for i
=
1,2
,
3,4 and perform 500 Monte Carlo runs on each gate for the two selection methods. Figure 5 shows the simulation results of OSPA distance, average selected sensor number, and average computation time (algorithms are implemented in MATLAB 2012a on a PC with 8 GB RAM and Intel Core i7-4770k CPU). With respect to sensor selection I, it is clear in Figure 5(a) that the OSPA distance barely changes with the size of SNMS gate while the number of allocated sensors and computation time are big when ϵ
1
=
500, which indicates that ϵ
1
=
500 is too small to suppress less informative sensors. For sensor selection II in Figure 5(b), we can see that with bigger gate the OSPA distance increases and the number of selected sensors decreases to 1 (minimum). This is because there are chances ofoversuppression when using sum HIG in SNMS sensor selection. Taking ϵ
4
=
i
n
f as illustration, sensor selection II actually only picks one most informative sensor among the network and suppress all the other sensors, whereas the most informative sensor can barely detect all targets for a common case. Hence, bigger value of SNMS gate may causeoversuppression in sensor selection II and consequently lose targets. Nevertheless, the advantage of sensor selection II over sensor selection I lies in fewer selected sensors and less computation time given an appropriate SNMS gate.Figure 5
Comparison of different gates.
(a)
Sensor selection I
(b)
Sensor selection IITo sum up, either of the proposed sensor selection approachs is shown effective and efficient in the simulations. Sensor selection with individual HIG will always provide satisfactory tracking performance while the number of selected sensors and computation time are relatively large. On the other hand, sensor selection with sum HIG may suffer fromoversuppression with an inappropriately big SNMS gate but it is much more efficient and can be as effective as sensor selection I when using a proper SNMS gate.


## 6. Conclusion

In this paper, we propose a novel sensor selection method for multitarget tracking in the static sensor network. We provide the SMC implementation of the sequential multi-Bernoulli filter to perform multisensor fusion for multitarget tracking. With respect to sensor selection, we propose the SNMS sensor selection approach by considering the locations and HIGs of sensors and present the SMC implementation of two SNMS sensor selection approaches, respectively, using individual HIG and sum HIG. We thoroughly compare the performance of two sensor selection approaches and analyze the impact the SNMS gate has on either sensor selection approach in simulations. It is shown by the simulation results that either SNMS sensor selection approach is efficient and effective for multitarget tracking in sensor network.Our future work is to consider more challenging measurement model, such as time-difference-of-arrival measurement or Doppler measurement, which is less informative than the bearing and range sensor used in the simulation. Besides, we will also analyze the proposed sensor selection method in the heterogeneous sensor network.



---

*Source: 148081-2015-06-29.xml*