# Distributed Multiagent Control Approach for Multitarget Tracking

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

---

## Abstract

In multiagent systems, tracking multiple targets is challenging for two reasons: firstly, it is nontrivial to dynamically deploy networked agents of different types for utility optimization; secondly, information fusion for multitarget tracking is difficult in the presence of uncertainties, such as data association, noise, and clutter. In this paper, we present a novel control approach in distributed manner for multitarget tracking. The control problem is modelled as a partially observed Markov decision process, which is a NP-hard combinatorial optimization problem, by seeking all possible combinations of control commands. To solve this problem efficiently, we assume that the measurement of each agent is independent of other agents’ behavior and provide a suboptimal multiagent control solution by maximizing the local Rényi divergence. In addition, we also provide the SMC implementation of the sequential multi-Bernoulli filter so that each agent can utilize the measurements from neighbouring agents to perform information fusion for accurate multitarget tracking. Numerical studies validate the effectiveness and efficiency of our multiagent control approach for multitarget tracking.


---

## Body

## 1. Introduction

Advances in microelectromechanical systems have significantly boosted the development of multiagent systems in the past two decades. Low-cost agents, for example, robots, unmanned vehicles, or autonomous platforms, with high mobility, various sensing types, and powerful communication ability, are capable of different tasks in complex environment, for example, environmental monitoring, target localization and tracking, and event recognition [1, 2]. Multiagent coordinated control, as the fundamental problems in multiagent system, has received increasing research interest for utility optimization recently [3–5]. An area that benefits greatly from multiagent system is target tracking [6]. However, very little progress has been made in this direction since it is an extremely challenging problem in two aspects: multiagent control and information fusion for multitarget tracking.The multiagent control problem in essence is a decision making issue which is to fulfill a given task in an optimal/suboptimal way. In this paper, we study the multiagent control problem particularly for multitarget tracking including accurate estimation of both target number and the location of each target. In the literature, researchers and practitioners have done extensive work on the agent control problem. Olfati-Saber et al. [6] established a consensus based Kalman filter for distributed single target tracking. In [7], multiagent control is performed by path following of a virtual leader. Olfati-Saber [8] provided flocking algorithms in both theory and applications to handle large number of agents. Unfortunately, none of them fit for the multitarget tracking issue well.Recently, random finite set (RFS) based Bayesian framework has opened doors for multisensor multiobject system and provides elegant mathematical tools to address multitarget tracking problem in multiagent systems [9]. The probability hypothesis density (PHD) [10], cardinalized PHD (CPHD) [11], and multi-Bernoulli filter [12] have been developed as approximations under different posterior assumptions. Gaussian mixture and particle implementation of these filters have been shown to be effective in different tracking applications [13–17]. Using tools from FISST, agent control can pose as a partially observed Markov decision process (POMDP) [18], which has been shown to be effective in a single agent case in recent works [19, 20].In this paper, we will further extend the work in [19, 20] to a more general problem: multiagent control for multitarget tracking, which is far more difficult than the single agent control case. We first model the multiagent coordinated control problem as a 1-step look-ahead POMDP since the multiple-step look-ahead is computationally intractable and show it is a NP combinatorial optimization problem when seeking all possible combinations of admissible control commands. Hence, we propose a suboptimal solution under the assumption that the measurements of each agent are independent of other agents’ behavior. The multiagent coordinated control is decoupled into distributed control of each agent by maximizing the local Rényi divergence between prior and posterior multitarget probability density. Besides, we present the sequential Monte Carlo (SMC) implementation of the sequential multi-Bernoulli filter for each agent to utilize measurements from neighbouring agents. Numerical simulations demonstrate the effectiveness and efficiency of our approach.The remainder of this paper is organized as follows. Section2 presents some preliminary knowledge of RFS based Bayesian framework to lay a foundation for the rest of this paper. In Section 3, we illustrate the distributed agent control approach and present its implementation in detail. Section 4 provides the information fusion scheme for multitarget tracking for each agent to utilize measurements from other agents. Section 5 provides numerical results that verify the proposed agent control and multitarget tracking approach.


## 2. Preliminary Knowledge of RFS Based Bayesian Framework

This section provides the basic concepts and notations of RFS based Bayesian framework. Section2.1 gives a general description of RFS and how to model multitarget by multi-Bernoulli RFS. Then, the RFS based Bayesian filtering is provided thereafter in Section 2.2.


### 2.1. Multi-Bernoulli RFS

A random finite set (RFS) is a random variable that takes values as unordered finite sets. The randomness of an RFS refers to two aspects: the set cardinality (number of elements of the set) is random; each element in the set is also a random variable. The probabilistic description of RFS has been studied regarding various types of probability distributions such as multi-Bernoulli (or Bernoulli) RFS, IID (short for independent identically distributed) cluster RFS, and Poisson RFS [21]. Here, we introduce the multi-Bernoulli RFS for multitarget state modelling, which offers a better alternative than the Poisson RFS and IID cluster RFS in applications with highly nonlinear model and/or nonhomogeneous sensor type [12].Assume the dimension of target state isn; 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 a single target; then, the multitarget state can be modeled as multi-Bernoulli RFSΞ with probability density given in [12] as follows: (1)
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 by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short in the following content.


### 2.2. Multisource Multiobject Bayesian Framework

Stochastic filtering in Bayesian framework has developed in decades [22]. Under the assumption of linear model and Gaussian distribution, Kalman filter is derived [23] and has been widely used for tracking since then. To extend the standard Bayesian framework to a multisource multitarget version, we need the help of RFS modelling. Let X
k and Z
k denote the state set and observation set, respectively, as follows: (2)
X
k
=
x
k
,
1
,
…
,
x
k
,
N
k
,
Z
k
=
Z
k
1
,
…
,
Z
k
W
,where W is the total number of agents (we treat each agent as a single sensor) and Z
k
i
=
z
k
i
,
…
,
z
k
,
M
k
,
i
i for i
=
1
,
…
,
W. N
k is the time-varying cardinality of targets, while M
k
,
i is the cardinality of the measurement set generated by agent i.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
i of agent i can be represented as a union of two parts: target-generated measurements Θ
k
i and clutter K
k
i; thus, Z
k
i
=
Θ
k
i
∪
K
k
i.Given the specific type of RFS, the Bayesian framework for optimal estimation via RFS which is the same form as the classical Bayesian filtering is given as follows:(3)
f
k
∣
k
-
1
X
∣
Z
1
:
k
-
1
=
∫
‍
f
k
∣
k
-
1
X
∣
X
′
f
k
-
1
X
′
∣
Z
1
:
k
-
1
δ
X
′
,
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. Under different assumptions of the RFS type, the PHD, CPHD, and multi-Bernoulli filter have been derived from (3) by finite set statistics [9]. Notice that the integrals in (3) are FISST set integrals (see [9–11] for details).


## 2.1. Multi-Bernoulli RFS

A random finite set (RFS) is a random variable that takes values as unordered finite sets. The randomness of an RFS refers to two aspects: the set cardinality (number of elements of the set) is random; each element in the set is also a random variable. The probabilistic description of RFS has been studied regarding various types of probability distributions such as multi-Bernoulli (or Bernoulli) RFS, IID (short for independent identically distributed) cluster RFS, and Poisson RFS [21]. Here, we introduce the multi-Bernoulli RFS for multitarget state modelling, which offers a better alternative than the Poisson RFS and IID cluster RFS in applications with highly nonlinear model and/or nonhomogeneous sensor type [12].Assume the dimension of target state isn; 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 a single target; then, the multitarget state can be modeled as multi-Bernoulli RFSΞ with probability density given in [12] as follows: (1)
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 by π
k
=
{
(
r
k
(
j
)
,
p
k
(
j
)
)
}
j
=
1
M
k for short in the following content.


## 2.2. Multisource Multiobject Bayesian Framework

Stochastic filtering in Bayesian framework has developed in decades [22]. Under the assumption of linear model and Gaussian distribution, Kalman filter is derived [23] and has been widely used for tracking since then. To extend the standard Bayesian framework to a multisource multitarget version, we need the help of RFS modelling. Let X
k and Z
k denote the state set and observation set, respectively, as follows: (2)
X
k
=
x
k
,
1
,
…
,
x
k
,
N
k
,
Z
k
=
Z
k
1
,
…
,
Z
k
W
,where W is the total number of agents (we treat each agent as a single sensor) and Z
k
i
=
z
k
i
,
…
,
z
k
,
M
k
,
i
i for i
=
1
,
…
,
W. N
k is the time-varying cardinality of targets, while M
k
,
i is the cardinality of the measurement set generated by agent i.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
i of agent i can be represented as a union of two parts: target-generated measurements Θ
k
i and clutter K
k
i; thus, Z
k
i
=
Θ
k
i
∪
K
k
i.Given the specific type of RFS, the Bayesian framework for optimal estimation via RFS which is the same form as the classical Bayesian filtering is given as follows:(3)
f
k
∣
k
-
1
X
∣
Z
1
:
k
-
1
=
∫
‍
f
k
∣
k
-
1
X
∣
X
′
f
k
-
1
X
′
∣
Z
1
:
k
-
1
δ
X
′
,
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. Under different assumptions of the RFS type, the PHD, CPHD, and multi-Bernoulli filter have been derived from (3) by finite set statistics [9]. Notice that the integrals in (3) are FISST set integrals (see [9–11] for details).


## 3. POMDP Based Distributed Multiagent Control Approach

In this section, we first illustrate multiagent control in the framework of a POMDP in Section3.1. Section 3.2 provides maximizing expected Rényi divergence between prior and posterior multitarget density as the objective function for the control scheme. By assuming the measurement of each agent is independent of other agents’ behavior, we propose a distributed agent control approach by maximizing the local Rényi divergence in Section 3.3.


### 3.1. POMDP Based Multiagent Control

We begin with the notations of using a POMDP as the solution for multiagent control. At timek, denote the control command of agent i by u
k
i
∈
U
k
i, where U
k
i is the set of all admissible control commands for agent i. Let U
k
∈
U
k denote the multiagent control command, and U
k is the set of all possible control command combinations. Then, U
k
=
{
u
k
1
,
…
,
u
k
i
,
…
,
u
k
W
} for i
=
1
,
…
,
W and W is the total number of agents.DefineD
(
V
,
f
,
Z
) as the objective function dependent on multiagent control command V, multitarget density f, and the associated measurement set Z when control command V was applied. The aim of multiagent control is to find the optimal multiagent control command U
k by maximizing/minimizing the statistical expectation of predefined objective function D
(
V
,
f
,
Z
) as (4)
U
k
=
a
r
g
m
a
x
m
i
n
V
∈
U
k
E
×
D
V
,
f
k
-
1
X
k
-
1
∣
Z
0
:
k
-
1
,
U
0
:
k
-
1
,
Z
k
V
,where f
k
-
1
(
X
k
-
1
∣
Z
0
:
k
-
1
,
U
0
:
k
-
1
) represents the multitarget posterior density after applying a sequence of multiagent control commands U
0
,
…
,
U
k
-
1. Remark that the general formulation of POMDP is a p-step future decision process of which the computational cost would grow exponentially with the number of future steps. In this paper, we only consider a one-step future decision described by (4) as an approximation.


### 3.2. Global Objective Function

As shown in (4), the objective function plays a crucial role in POMDP based multiagent control problem. Information theoretic method is a typical objective function for sensor management. Here, we propose maximizing the information gain of multitarget prior and posterior density as the objective function for tracking. The Rényi divergence, also known as alpha divergence, measures the information gain between any two probability densities. The objective function for multitarget tracking is given as follows: (5)
R
U
k
=
1
α
-
1
×
l
o
g
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
H
H
H
i
H
H
×
f
k
+
1
∣
k
X
k
+
1
∣
Z
1
:
k
δ
X
k
+
1
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
×
p
Z
k
+
1
∣
Z
1
:
k
,
U
k
α
-
1
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
,where p
(
Z
k
+
1
∣
Z
1
:
k
,
U
k
)
=
∫
‍
g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
,
U
k
)
f
k
+
1
∣
k
(
X
k
+
1
∣
Z
1
:
k
)
δ
X and α is a parameter that determines how much we emphasize the tails of two densities in the metric. Notice that the Rényi divergence becomes the Kullback-Leibler discrimination and Hellinger affinity, respectively, when α
→
1 and α
=
0.5 [24].To compute the expectation of (5), we introduce the SMC implementation of the objective function. At time k, assume that the multitarget predicted density f
k
+
1
∣
k
(
X
k
∣
Z
1
:
k
) is given in SMC form; that is, f
k
+
1
∣
k
(
X
k
∣
Z
1
:
k
)
=
∑
j
=
1
S
‍
ω
j
·
δ
X
k
+
1
j
(
X
), where each multitarget particle X
k
+
1
j
=
{
x
k
+
1,1
j
,
…
,
x
k
+
1
,
n
k
j
j
}, n
k
+
1
j
=
X
k
+
1
j is the number of targets, x
k
+
1
,
i
j for i
=
1
,
…
,
n
k
+
1
j represent target position in the state space, and ω
j is the weight associated to particle X
k
+
1
j. Notice that X
k
+
1
j is a particle sampled from a RFS, which accounts for the randomness of both the cardinality and target positions of particle X
k
+
1
j. Given the Bayesian prediction and update equations (3), we obtain (6)
E
R
U
k
,
Z
k
+
1
=
1
α
-
1
l
o
g
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
α
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
α
,where the multisensor multitarget likelihood function dependent on multiagent control g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
j
,
U
k
) is given by (7)
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
=
∏
i
=
1
W
‍
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
U
kand the single-sensor multitarget likelihood function dependent on multiagent control g
k
+
1
i
(
Z
k
+
1
i
∣
X
k
+
1
j
,
U
k
) is given in [9] as follows: (8)
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
U
k
=
f
c
Z
k
+
1
i
1
-
p
D
i
n
k
+
1
j
×
∑
θ
‍
∏
l
:
θ
l
>
0
‍
p
D
i
g
k
+
1
i
z
θ
l
i
∣
x
k
+
1
,
l
j
,
U
k
1
-
p
D
i
λ
c
z
θ
l
i
,where g
k
+
1
i
(
z
θ
(
l
)
i
∣
x
k
+
1
,
l
j
,
U
k
) is the standard single-sensor single-target likelihood function described by the measurement model of agent i; f
c
(
Z
k
+
1
i
) is the probability density of clutter RFS and f
c
(
Z
k
+
1
i
)
=
e
-
λ
∏
z
i
∈
Z
k
+
1
i
‍
λ
c
(
z
i
) for Poisson RFS; λ is the clutter rate, while c
(
·
) is the probability distribution of clutter; l
:
θ
(
l
)
>
0 represents all possible associations between the particle set X
k
+
1
j and the measurement set Z
k
+
1
i; that is, θ
:
{
1
,
…
,
n
k
+
1
j
}
→
{
0,1
,
…
,
M
k
,
i
}. Notice that (6) is a multisensor multitarget case, whereas the derivation is similar to the single sensor case in [19]. Thus, we directly omit the tedious proof here.


### 3.3. Distributed Agent Control

Even though multiagent control can be described as a POMDP given by (4) and (5), it is still intractable to achieve the global optimum of defined objective function for two reasons: firstly, searching all admissible control command combinations is a NP-hard combinatorial optimization problem [25]; secondly, global optimization requires the existence of a centralized fusion center that receives information from all agents, which is unrealistic for most large-scale multiagent systems. Hence, we propose a distributed sensor control method that compromises the local optimum of the global objective function, which is computationally tractable and convenient to implement.Instead of computing the global information gain, we consider finding the optimal command for each individual agent. Similar to (5), the Rényi divergence for the ith agent is given by (9)
E
R
i
u
k
i
,
Z
k
+
1
=
1
α
-
1
l
o
g
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
α
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
α
.Assume that measurement set generated by one agent is independent of other agents’ behaviour; then, the multisensor multitarget likelihood function g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
) in (9) can be written as follows: (10)
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
=
∏
l
=
1
W
‍
g
k
+
1
l
Z
k
+
1
l
∣
X
k
+
1
j
,
u
k
i
=
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
u
k
i
∏
l
=
1
,
l
≠
i
W
‍
g
k
+
1
l
Z
k
+
1
l
∣
X
k
+
1
j
,where g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
,
u
k
i
)
=
g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
) for l
≠
i. The future measurement Z
k
+
1 is generated assuming no clutter and unity detection rate as illustrated in [20]. For agent i, Z
k
+
1
i is predicted based on the multitarget state and possible control command u
k
i, while Z
k
+
1
l for l
≠
i is predicted based on the multitarget state and current location of agent l. Therefore, we can find the optimal control command for each individual agent and then combine them together to form the total control command set, which is given as follows: (11)
u
k
i
=
a
r
g
m
a
x
u
∈
V
k
E
R
i
u
,
Z
k
+
1
,
(12)
U
k
=
⋃
i
=
1
W
u
k
i
.Remark 1.
Equation (9) is exactly a particular case of (6) by assuming all agents except agent i stay still at current time. The proposed distributed control approach is a suboptimal solution of (6) by seeking its local optimum, and the local optimum means that u
k
i is optimal given all other agents keep still.When computing the optimal command of agenti, the distributed agent control given by (11) approximates the likelihood function g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
) via using current locations of agent l for l
≠
i. Hence, we refer to (9) as the “local Rényi divergence” in this paper.The proposed distributed control approach can significantly reduce the computational cost to perform real-time agent control. For illustration, assume that all agents have the same number of possible control command denoted bym; the computation complexity of our control approach is O
(
W
·
m
) which is much smaller than O
(
m
W
) by searching all control combinations, especially in large-scale multiagent system.Generally speaking, each agent can only communicate with its neighbouringW
k
′ agents and obtain their current locations s
k
l as well as their measurement model g
l
(
z
∣
x
) for l
=
1
,
…
,
W
k
′. Here, we assume that the information received by each agent is accurate without any input saturation. For a more challenging case that there is input saturation described in [26, 27], the topic is beyond the scope of this paper. The neighbouring relationship may change over time due to the relative movement of agents. Algorithm 1 provides the SMC implementation of proposed distributed control for each agent by maximizing the local Rényi divergence.Algorithm 1:Control law of agent i.
Input: predicted multi-target density π
k
+
1
∣
k, neighbouring agents s
k
l and g
l
(
z
∣
x
) for l
=
1
,
…
,
W
k
′
(1)   Sampling the multi-object state space{
X
k
+
1
∣
k
j
}
j
=
1
S
~
π
k
+
1
∣
k;
(2)   Compute predicted estimation: target numbern
^
k
+
1
∣
k and state X
^
k
+
1
∣
k
=
∪
j
=
1
n
^
k
+
1
∣
k
{
x
k
+
1
∣
k
j
};
(3)for  l
=
1
,
…
,
W
k
′  do
(4)    PredictZ
k
+
1
l based on X
^
k
+
1
∣
k and g
l
(
z
∣
x
);
(5)    Compute single-sensor multi-target likelihood given by (8) to obtain {
g
(
Z
k
+
1
l
∣
X
k
+
1
∣
k
j
)
}
j
=
1
S;
(6)end for
(7)for each u
k
i
∈
U
k
i  do
(8)    PredictZ
k
+
1
i based on u
k
i and g
i
(
z
∣
x
);
(9)    Compute single-sensor multi-target likelihood given by (8) to obtain {
g
(
Z
k
+
1
i
∣
X
k
+
1
∣
k
j
,
u
k
i
)
}
j
=
1
S;
(10)  Compute multi-sensor multi-target likelihood given by (10) to obtain {
g
(
Z
k
+
1
∣
X
k
+
1
∣
k
j
,
u
k
i
)
}
j
=
1
S;
(11)  Compute the valueρ
(
u
k
i
)
=
E
[
R
i
] given by (9);
(12)end for
Output:  u
^
k
i
=
arg
max
⁡
ρ
(
u
k
i
)


## 3.1. POMDP Based Multiagent Control

We begin with the notations of using a POMDP as the solution for multiagent control. At timek, denote the control command of agent i by u
k
i
∈
U
k
i, where U
k
i is the set of all admissible control commands for agent i. Let U
k
∈
U
k denote the multiagent control command, and U
k is the set of all possible control command combinations. Then, U
k
=
{
u
k
1
,
…
,
u
k
i
,
…
,
u
k
W
} for i
=
1
,
…
,
W and W is the total number of agents.DefineD
(
V
,
f
,
Z
) as the objective function dependent on multiagent control command V, multitarget density f, and the associated measurement set Z when control command V was applied. The aim of multiagent control is to find the optimal multiagent control command U
k by maximizing/minimizing the statistical expectation of predefined objective function D
(
V
,
f
,
Z
) as (4)
U
k
=
a
r
g
m
a
x
m
i
n
V
∈
U
k
E
×
D
V
,
f
k
-
1
X
k
-
1
∣
Z
0
:
k
-
1
,
U
0
:
k
-
1
,
Z
k
V
,where f
k
-
1
(
X
k
-
1
∣
Z
0
:
k
-
1
,
U
0
:
k
-
1
) represents the multitarget posterior density after applying a sequence of multiagent control commands U
0
,
…
,
U
k
-
1. Remark that the general formulation of POMDP is a p-step future decision process of which the computational cost would grow exponentially with the number of future steps. In this paper, we only consider a one-step future decision described by (4) as an approximation.


## 3.2. Global Objective Function

As shown in (4), the objective function plays a crucial role in POMDP based multiagent control problem. Information theoretic method is a typical objective function for sensor management. Here, we propose maximizing the information gain of multitarget prior and posterior density as the objective function for tracking. The Rényi divergence, also known as alpha divergence, measures the information gain between any two probability densities. The objective function for multitarget tracking is given as follows: (5)
R
U
k
=
1
α
-
1
×
l
o
g
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
H
H
H
i
H
H
×
f
k
+
1
∣
k
X
k
+
1
∣
Z
1
:
k
δ
X
k
+
1
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
×
p
Z
k
+
1
∣
Z
1
:
k
,
U
k
α
-
1
∫
g
k
+
1
Z
k
+
1
∣
X
k
+
1
,
U
k
α
,where p
(
Z
k
+
1
∣
Z
1
:
k
,
U
k
)
=
∫
‍
g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
,
U
k
)
f
k
+
1
∣
k
(
X
k
+
1
∣
Z
1
:
k
)
δ
X and α is a parameter that determines how much we emphasize the tails of two densities in the metric. Notice that the Rényi divergence becomes the Kullback-Leibler discrimination and Hellinger affinity, respectively, when α
→
1 and α
=
0.5 [24].To compute the expectation of (5), we introduce the SMC implementation of the objective function. At time k, assume that the multitarget predicted density f
k
+
1
∣
k
(
X
k
∣
Z
1
:
k
) is given in SMC form; that is, f
k
+
1
∣
k
(
X
k
∣
Z
1
:
k
)
=
∑
j
=
1
S
‍
ω
j
·
δ
X
k
+
1
j
(
X
), where each multitarget particle X
k
+
1
j
=
{
x
k
+
1,1
j
,
…
,
x
k
+
1
,
n
k
j
j
}, n
k
+
1
j
=
X
k
+
1
j is the number of targets, x
k
+
1
,
i
j for i
=
1
,
…
,
n
k
+
1
j represent target position in the state space, and ω
j is the weight associated to particle X
k
+
1
j. Notice that X
k
+
1
j is a particle sampled from a RFS, which accounts for the randomness of both the cardinality and target positions of particle X
k
+
1
j. Given the Bayesian prediction and update equations (3), we obtain (6)
E
R
U
k
,
Z
k
+
1
=
1
α
-
1
l
o
g
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
α
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
α
,where the multisensor multitarget likelihood function dependent on multiagent control g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
j
,
U
k
) is given by (7)
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
U
k
=
∏
i
=
1
W
‍
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
U
kand the single-sensor multitarget likelihood function dependent on multiagent control g
k
+
1
i
(
Z
k
+
1
i
∣
X
k
+
1
j
,
U
k
) is given in [9] as follows: (8)
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
U
k
=
f
c
Z
k
+
1
i
1
-
p
D
i
n
k
+
1
j
×
∑
θ
‍
∏
l
:
θ
l
>
0
‍
p
D
i
g
k
+
1
i
z
θ
l
i
∣
x
k
+
1
,
l
j
,
U
k
1
-
p
D
i
λ
c
z
θ
l
i
,where g
k
+
1
i
(
z
θ
(
l
)
i
∣
x
k
+
1
,
l
j
,
U
k
) is the standard single-sensor single-target likelihood function described by the measurement model of agent i; f
c
(
Z
k
+
1
i
) is the probability density of clutter RFS and f
c
(
Z
k
+
1
i
)
=
e
-
λ
∏
z
i
∈
Z
k
+
1
i
‍
λ
c
(
z
i
) for Poisson RFS; λ is the clutter rate, while c
(
·
) is the probability distribution of clutter; l
:
θ
(
l
)
>
0 represents all possible associations between the particle set X
k
+
1
j and the measurement set Z
k
+
1
i; that is, θ
:
{
1
,
…
,
n
k
+
1
j
}
→
{
0,1
,
…
,
M
k
,
i
}. Notice that (6) is a multisensor multitarget case, whereas the derivation is similar to the single sensor case in [19]. Thus, we directly omit the tedious proof here.


## 3.3. Distributed Agent Control

Even though multiagent control can be described as a POMDP given by (4) and (5), it is still intractable to achieve the global optimum of defined objective function for two reasons: firstly, searching all admissible control command combinations is a NP-hard combinatorial optimization problem [25]; secondly, global optimization requires the existence of a centralized fusion center that receives information from all agents, which is unrealistic for most large-scale multiagent systems. Hence, we propose a distributed sensor control method that compromises the local optimum of the global objective function, which is computationally tractable and convenient to implement.Instead of computing the global information gain, we consider finding the optimal command for each individual agent. Similar to (5), the Rényi divergence for the ith agent is given by (9)
E
R
i
u
k
i
,
Z
k
+
1
=
1
α
-
1
l
o
g
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
α
∑
j
=
1
S
‍
ω
j
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
α
.Assume that measurement set generated by one agent is independent of other agents’ behaviour; then, the multisensor multitarget likelihood function g
k
+
1
(
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
) in (9) can be written as follows: (10)
g
k
+
1
Z
k
+
1
∣
X
k
+
1
j
,
u
k
i
=
∏
l
=
1
W
‍
g
k
+
1
l
Z
k
+
1
l
∣
X
k
+
1
j
,
u
k
i
=
g
k
+
1
i
Z
k
+
1
i
∣
X
k
+
1
j
,
u
k
i
∏
l
=
1
,
l
≠
i
W
‍
g
k
+
1
l
Z
k
+
1
l
∣
X
k
+
1
j
,where g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
,
u
k
i
)
=
g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
) for l
≠
i. The future measurement Z
k
+
1 is generated assuming no clutter and unity detection rate as illustrated in [20]. For agent i, Z
k
+
1
i is predicted based on the multitarget state and possible control command u
k
i, while Z
k
+
1
l for l
≠
i is predicted based on the multitarget state and current location of agent l. Therefore, we can find the optimal control command for each individual agent and then combine them together to form the total control command set, which is given as follows: (11)
u
k
i
=
a
r
g
m
a
x
u
∈
V
k
E
R
i
u
,
Z
k
+
1
,
(12)
U
k
=
⋃
i
=
1
W
u
k
i
.Remark 1.
Equation (9) is exactly a particular case of (6) by assuming all agents except agent i stay still at current time. The proposed distributed control approach is a suboptimal solution of (6) by seeking its local optimum, and the local optimum means that u
k
i is optimal given all other agents keep still.When computing the optimal command of agenti, the distributed agent control given by (11) approximates the likelihood function g
k
+
1
l
(
Z
k
+
1
l
∣
X
k
+
1
j
) via using current locations of agent l for l
≠
i. Hence, we refer to (9) as the “local Rényi divergence” in this paper.The proposed distributed control approach can significantly reduce the computational cost to perform real-time agent control. For illustration, assume that all agents have the same number of possible control command denoted bym; the computation complexity of our control approach is O
(
W
·
m
) which is much smaller than O
(
m
W
) by searching all control combinations, especially in large-scale multiagent system.Generally speaking, each agent can only communicate with its neighbouringW
k
′ agents and obtain their current locations s
k
l as well as their measurement model g
l
(
z
∣
x
) for l
=
1
,
…
,
W
k
′. Here, we assume that the information received by each agent is accurate without any input saturation. For a more challenging case that there is input saturation described in [26, 27], the topic is beyond the scope of this paper. The neighbouring relationship may change over time due to the relative movement of agents. Algorithm 1 provides the SMC implementation of proposed distributed control for each agent by maximizing the local Rényi divergence.Algorithm 1:Control law of agent i.
Input: predicted multi-target density π
k
+
1
∣
k, neighbouring agents s
k
l and g
l
(
z
∣
x
) for l
=
1
,
…
,
W
k
′
(1)   Sampling the multi-object state space{
X
k
+
1
∣
k
j
}
j
=
1
S
~
π
k
+
1
∣
k;
(2)   Compute predicted estimation: target numbern
^
k
+
1
∣
k and state X
^
k
+
1
∣
k
=
∪
j
=
1
n
^
k
+
1
∣
k
{
x
k
+
1
∣
k
j
};
(3)for  l
=
1
,
…
,
W
k
′  do
(4)    PredictZ
k
+
1
l based on X
^
k
+
1
∣
k and g
l
(
z
∣
x
);
(5)    Compute single-sensor multi-target likelihood given by (8) to obtain {
g
(
Z
k
+
1
l
∣
X
k
+
1
∣
k
j
)
}
j
=
1
S;
(6)end for
(7)for each u
k
i
∈
U
k
i  do
(8)    PredictZ
k
+
1
i based on u
k
i and g
i
(
z
∣
x
);
(9)    Compute single-sensor multi-target likelihood given by (8) to obtain {
g
(
Z
k
+
1
i
∣
X
k
+
1
∣
k
j
,
u
k
i
)
}
j
=
1
S;
(10)  Compute multi-sensor multi-target likelihood given by (10) to obtain {
g
(
Z
k
+
1
∣
X
k
+
1
∣
k
j
,
u
k
i
)
}
j
=
1
S;
(11)  Compute the valueρ
(
u
k
i
)
=
E
[
R
i
] given by (9);
(12)end for
Output:  u
^
k
i
=
arg
max
⁡
ρ
(
u
k
i
)


## 4. Multisensor Fusion for Multi-Target Tracking

As mentioned before, the multi-Bernoulli filter outweighs the PHD/CPHD filter in the SMC implementation for nonlinear problem since the state extraction in multi-Bernoulli filter is not dependent upon the heuristics in clustering but is dependent only on the Bernoulli parameters. Hence, the multi-Bernoulli filter has been used extensively in computer vision [16, 17], robot SLAM [28], and sensor network [29]. Besides, the state-of-the-art development of the multi-Bernoulli filter offers the power to directly produce tracks of individual targets, which is known as the labelled multi-Bernoulli filter in the community [30].In this section, we first briefly review the cardinality-balanced multi-Bernoulli filter given in [12] in Section 4.1. Then, Section 4.2 provides the sequential update scheme for information fusion of multiple sensors.


### 4.1. 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 (13)
π
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 (14)
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 as (15)
π
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(16)
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
=
∑
‍
M
k
+
1
∣
k
j
=
1
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
+
∑
‍
M
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
-
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
,
ψ
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. 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.


### 4.2. Sequential Multisensor Multi-Bernoulli Filter

In multisensor multitarget tracking scenario, there is no unified multisensor fusion method which is tractable and computationally acceptable. Sequential update has been widely used and verified to be a good approximation for information fusion of multiple sensors. Here, multitarget tracking in multiagent network is implemented via multi-Bernoulli filter with sequential update scheme based on the SMC implementation.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, we give the SMC implementation of sequential multisensor multi-Bernoulli filter in Algorithm 2. We refer the readers to Section 4.1 of [12] for detailed equations.Algorithm 2: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)    fori
=
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)    computer
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)    fori
=
1
,
…
,
L
Γ
,
k
(
j
), sample x
Γ
,
i
,
k
+
1
(
j
), compute weight ω
Γ
,
i
,
k
+
1
(
j
) and normalization ω
~
Γ
,
i
,
k
+
1
(
j
)
(7)    computer
Γ
,
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; neighbouring agents s
k
l, g
l
(
z
∣
x
) and Z
k
+
1
l for l
=
1
,
…
,
W
k
′
(1)for  l
=
1
,
…
,
W
k
′  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)     computer
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)      computer
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
)
,
W
k
′
,
p
L
,
k
+
1
(
j
)
,
W
k
′
)
}
j
=
1
M
k
+
1
∣
k
∪
{
(
r
U
,
k
+
1
W
k
′
(
z
)
,
p
U
,
k
+
1
W
k
′
(
⋅
∣
z
)
)
}
z
∈
Z
k
+
1
W
k
′The superscript(
j
)
,
l in Algorithm 2 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. Notice that the multi-Bernoulli filter we adopt here cannot produce tracks directly, and it can be replaced by the labelled multi-Bernoulli filter to produce individual tracks at the cost of some extra computation.


## 4.1. 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 (13)
π
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 (14)
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 as (15)
π
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(16)
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
=
∑
‍
M
k
+
1
∣
k
j
=
1
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
+
∑
‍
M
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
-
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
,
ψ
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. 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.


## 4.2. Sequential Multisensor Multi-Bernoulli Filter

In multisensor multitarget tracking scenario, there is no unified multisensor fusion method which is tractable and computationally acceptable. Sequential update has been widely used and verified to be a good approximation for information fusion of multiple sensors. Here, multitarget tracking in multiagent network is implemented via multi-Bernoulli filter with sequential update scheme based on the SMC implementation.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, we give the SMC implementation of sequential multisensor multi-Bernoulli filter in Algorithm 2. We refer the readers to Section 4.1 of [12] for detailed equations.Algorithm 2: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)    fori
=
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)    computer
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)    fori
=
1
,
…
,
L
Γ
,
k
(
j
), sample x
Γ
,
i
,
k
+
1
(
j
), compute weight ω
Γ
,
i
,
k
+
1
(
j
) and normalization ω
~
Γ
,
i
,
k
+
1
(
j
)
(7)    computer
Γ
,
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; neighbouring agents s
k
l, g
l
(
z
∣
x
) and Z
k
+
1
l for l
=
1
,
…
,
W
k
′
(1)for  l
=
1
,
…
,
W
k
′  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)     computer
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)      computer
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
)
,
W
k
′
,
p
L
,
k
+
1
(
j
)
,
W
k
′
)
}
j
=
1
M
k
+
1
∣
k
∪
{
(
r
U
,
k
+
1
W
k
′
(
z
)
,
p
U
,
k
+
1
W
k
′
(
⋅
∣
z
)
)
}
z
∈
Z
k
+
1
W
k
′The superscript(
j
)
,
l in Algorithm 2 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. Notice that the multi-Bernoulli filter we adopt here cannot produce tracks directly, and it can be replaced by the labelled multi-Bernoulli filter to produce individual tracks at the cost of some extra computation.


## 5. Simulation

In order to demonstrate the performance of proposed multiagent control approach for multitarget tracking, we present numerical results for a planar multitarget tracking scenario in multiagent system where three controllable moving observers, two equipped with range-only sensors and one with a bearing-only sensor, are placed in a specified surveillance area of size[
0,1000
m
]
×
[
0,1000
m
] to estimate the number of targets as well as their positions. Each agent shares its current locations, sensor type, and observations with the other agents. It is intuitive that the initial positions of agents have an impact on the agent control and target tracking procedure. However, we are not going to involve the network topology issue in this paper since our approach is supposed to act independently from the network topology which is different from the method described in [31, 32].There are unknown and time-varying numbers of targets observed in clutter for each agent. Assume that targets move according to the nearly constant velocity model given by(17)
x
k
=
F
x
k
-
1
+
G
w
v
,
k
,where x
k
=
[
p
x
,
k
,
v
x
,
k
,
p
y
,
k
,
v
y
,
k
]
T; p
k
=
[
p
x
,
k
,
p
y
,
k
]
T are planar position and v
k
=
[
v
x
,
k
,
v
y
,
k
]
T are planar velocity, respectively, along x-coordinate, y-coordinate. F
=
I
2
⊗
1
T
0
1; G
=
I
2
⊗
T
2
/
2
T. I
2 is 2
×
2 identity matrix and ⊗ denotes Kronecker product. T is the sampling period, and w
v
,
k
~
N
(
0
,
Q
k
) is a 2
×
1 IID 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.Measurement of sensorl originated from target with state x
k is noisy vector of range or bearing measurement, and the measurement model for range-only sensor is given by (18)
z
r
,
k
l
=
p
k
-
s
l
+
w
r
land, for bearing-only sensor, (19)
z
ϕ
,
k
l
=
a
r
c
t
a
n
y
l
-
p
y
,
k
x
l
-
p
x
,
k
+
w
ϕ
l
,where s
k
l
=
[
x
k
l
,
y
k
l
] is the location of agent l and p
k
-
s
k
l
=
(
x
k
l
-
p
x
,
k
)
2
+
(
y
k
l
-
p
y
,
k
)
2 is the Euclidean distance between sensor l and the target. w
r
l is zero mean Gaussian noise N
(
·
;
0
,
(
σ
r
l
)
2
) and w
ϕ
l is also zero mean Gaussian noise N
(
·
;
0
,
(
σ
ϕ
l
)
2
). The standard derivations of w
r
l and w
ϕ
l are, respectively, given by σ
r
l and σ
ϕ
l as follows: (20)
σ
r
,
k
l
=
σ
0
+
β
r
p
k
-
s
l
2
σ
ϕ
,
k
l
=
σ
1
+
β
ϕ
p
k
-
s
l
,with σ
0
=
1
m, β
r
=
5
×
1
0
-
5
m
-
1, σ
1
=
π
/
180
r
a
d, and β
ϕ
=
1
0
-
5
r
a
d
·
m
-
1.Targets can appear or disappear in the scene at any time, and survival probabilityp
S
=
0.95 for each existing target. New born targets appear spontaneously according to γ
k
=
0.2
N
(
·
;
x
-
,
Q
). Three targets are presented as illustrated in Figure 1, where x
-
1
=
[
800
m
;
-
5
m
/
s
;
300
m
;
0
m
/
s
]
T, x
-
2
=
[
650
m
;
0.5
m
/
s
;
800
m
;
-
6
m
/
s
]
T, and x
-
3
=
[
300
m
;
5
m
/
s
;
600
m
;
-
8
m
/
s
]
T, respectively, for each target, and Q
=
d
i
a
g
(
[
50
m
,
2
m
/
s
,
50
m
,
2
m
/
s
]
2
) are identical for all three targets.Figure 1
Target tracks. Start/stop positions for each track are shown with•/■. Colored ⧫ are initial locations of agents.The probability of detection for both range-only sensor and bearing-only sensor is modelled by(21)
p
D
x
k
=
0.99
p
k
-
s
l
≤
R
0
m
a
x
0,0.99
-
c
p
k
-
s
l
p
k
-
s
l
≥
R
0
,with R
0
=
320
m and c
=
0.002
m
-
1. The clutter rate of each sensor λ
c
=
5 per scan. The standard derivation of process noise σ
v
=
1
m
/
s for both v
x
,
k and v
y
,
k. Given the current agent location s
k
l, the set of admissible control commands for each agent is computed as U
k
l
=
{
[
x
k
l
+
i
Δ
R
c
o
s
(
j
Δ
ϕ
)
,
y
k
l
+
i
Δ
R
c
o
s
(
j
Δ
ϕ
)
]
}
i
=
0
,
…
,
N
ϕ
j
=
0
,
…
,
N
R, where Δ
ϕ
=
2
π
/
N
ϕ and Δ
R
=
30
m are for angular and radial step size. N
R
=
2 and N
ϕ
=
8 here.The initial locations of agents are shown in Figure1, and each agent runs for 30 scans with sampling period T
=
1 s. Birth intensity for the multi-Bernoulli filter is approximated using adaptive target birth intensity sampling technique described in [33]. L
m
i
n
=
300 and L
m
a
x
=
500 are the minimum and maximum particle numbers of each Bernoulli set for track maintenance.Since each agent shares its location and observation with other two agents, agents have pretty much the same performance in this scene. Hence, we take agent 3 (bearing-only) as an example to illustrate the tracking procedure. Figure2 presents three key frames of one trial run of the proposed control and tracking approach. It is clear that, with the movement of agents, scattered particles converge to the ground truth locations of targets. This is because each agent is moving forward to a more informative direction so that “ghost” particles that are not generated by actual targets can be quickly eliminated by the multisensor fusion scheme. As a result, the estimation of target positions is getting more accurate over time. Besides, it can be seen from Figure 2(c) that the trajectories of agents are different for their distinctive initial locations and sensor types.Figure 2
One trial result. Green dots are for particles. Colored lines are trajectories of agents. Colored ⧫ are current locations of agents. Red∘ are ground truth of target stop positions, while black ∗ are for target estimation.
(a)
Timek
=
1
(b)
Timek
=
15
(c)
Timek
=
30Due to stochastic nature of our control and tracking approach, we adopt 1000 Monte Carlo runs to evaluate their performance. The optimal subpattern assignment (OSPA) metric composed of location error and cardinality error is used for tracking performance evaluation [34]. Figure 3 shows the OSPA distance of our approach, (c
=
100
m, p
=
1) from 1000 Monte Carlo runs, from which we can see that the OSPA distance converges with the movement of agents over time. The location error is reduced gradually since each agent is obtaining more informative measurements to lower the covariance of position estimation, while the estimation of target number quickly converges to a relatively low value. The average computation time for a single run is only approximately 0.82 s (algorithm is implemented in MATLAB 2012a on a PC with 8 GB RAM and Intel Core i7-4770k CPU). Several runs have been recorded in videos attached as Supplementary Material (available online at http://dx.doi.org/10.1155/2015/903682) for demonstration.Figure 3
OSPA distance.


## 6. Conclusion

In this paper, we propose a novel distributed multiagent control approach by maximizing the local Rényi divergence. The SMC implementation of the sequential multi-Bernoulli filter is provided for each agent to utilize the information from neighbouring agents. Simulation results demonstrate that the proposed approach is capable of distributed multitarget tracking via effective sensor control.Our future work is to use convex relaxation method for seeking global optimal solution for multiagent control and compare with the approach proposed in this paper. We also need to consider more challenging measurement model, such as time-difference-of-arrival measurement or Doppler measurement, which depends on the behavior of multiple agents.



---

*Source: 903682-2015-09-08.xml*