PARAMETRIC DATA ASSOCIATION PRIOR FOR MULTI-TARGET

TRACKING BASED ON RAO-BLACKWELLIZED MONTE CARLO

DATA ASSOCIATION

Oliver Greß and Stefan Posch

Institute of Computer Science, Martin Luther University Halle-Wittenberg,

Von-Seckendorff-Platz 1, Halle (Saale), Germany

Keywords:

Probabilistic Tracking, Multiple Targets, Data Association, Monte Carlo.

Abstract:

Association of observations to underlying targets is a crucial task in probabilistic tracking of multiple targets.

The Rao-Blackwellized Monte Carlo Data Association (RBMCDA) framework circumvents the combinatorial

explosion by approximating the joint distribution of targets and association variables by Monte Carlo samples

in the space of association variables. We present a parametric data association prior distribution required by

RBMCDA, which models the formation of observations. To sample from this distribution an efﬁcient algo-

rithm is developed. The Interacting Multiple Models (IMM) ﬁlter is integrated into the RBMCDA framework

to model the changing dynamics of targets aiming at tracking small particles in microscopy images. The

proposed method is evaluated in a proof of concept and evaluated using synthetic data.

1 INTRODUCTION

Probabilistic tracking techniques are beneﬁcial to

track multiple targets with fast motion and similar ap-

pearance. Such methods employ a dynamic model to

incorporate prior information. In detail, the state dis-

tribution of targets is predicted in each time step us-

ing the dynamic model and previous estimates. Sub-

sequently this estimate is corrected using an obser-

vation of the target, if available. The major issue in

multi-target tracking is the combinatorial problem to

associate observations to underlying targets.

In (S

arkk

a et al., 2007) the Rao-Blackwellized

Monte Carlo Data Association (RBMCDA) frame-

work was presented for tracking targets with linear

dynamics. RBMCDA circumvents the combinato-

rial problem of data association by sampling in the

space of association variables, while distributions of

targets’ states are kept in analytical form. It also

allows to integrate the Interacting Multiple Models

(IMM) ﬁlter (Bar-Shalom and Blair, 2000) to account

for changing dynamics of targets, which is required

by the application we are aiming at, namely tracking

of small particles in ﬂuorescence microscopy images.

We propose a parametric data association prior

required for sampling of associations, which models

the formation of observations, because of the lack of

a generic deﬁnition in (S

arkk

a et al., 2007).

The main contribution of this work is (i) a rigorous

deﬁnition of the parametric data association prior and

(ii) the development of an efﬁcient algorithm to carry

out the sampling procedure. We give a proof of con-

cept of the proposed method evaluating performance

of tracking on synthetic data.

The remainder of this work is structured as fol-

lows: Sec. 2 gives an overview of related work con-

cerned with the problem of tracking multiple targets.

An introduction to probabilistic multi-target tracking

and the RBMCDA framework is given in Sec. 3. In

Sec. 4 the parametric data association prior is devel-

oped as well as an efﬁcient algorithm for sampling.

Employed models, generation of data, methods for

evaluation, and the complexity are described in Sec. 5.

The performance is discussed in Sec. 6 and some con-

cluding remarks are given in Sec. 7.

2 RELATED WORK

Various approaches exist for tracking multiple targets.

Multi-Hypothesis Tracking (MHT) (Reid, 1979) pro-

vides an optimal solution to multi-target tracking for

linear systems with Gaussian noise in theory. It fur-

ther allows for the generation of new targets, extinc-

tion of targets as well as classifying observations as

387

Greß O. and Posch S..

PARAMETRIC DATA ASSOCIATION PRIOR FOR MULTI-TARGET TRACKING BASED ON RAO-BLACKWELLIZED MONTE CARLO DATA ASSOCIA-

TION.

DOI: 10.5220/0003863103870394

In Proceedings of the International Conference on Computer Vision Theory and Applications (VISAPP-2012), pages 387-394

ISBN: 978-989-8565-04-4

 2012 SCITEPRESS (Science and Technology Publications, Lda.)

clutter. In detail, MHT maintains a tree of hypotheses

for all possible associations of observations, which re-

sults in an exponentially growing number of hypothe-

ses. However, for practical application this number

must be reduced rigorously by heuristic strategies like

tree pruning or gating of observations.

The Joint Probabilistic Data Association (JPDA)

ﬁlter (Bar-Shalom and Fortmann, 1987; Bar-Shalom

and Blair, 2000) updates target state distributions by

observations in a gate around the target. The up-

dated state distribution results in a Gaussian mix-

ture which is subsequently approximated by a sin-

gle Gaussian to prevent exponential growth. JPDA

was developed to provide robust estimates in environ-

ments with many clutter observations (Bar-Shalom

and Fortmann, 1987). However, JPDA is limited to

a ﬁxed and known number of targets. JPDA was

combined with IMM ﬁltering (Bar-Shalom and Blair,

2000) to account for changing dynamics of targets and

extended to overcome the limitation to a ﬁxed number

of targets (Musicki and Evans, 2002).

Sequential Monte Carlo (SMC) methods approx-

imate state distributions not in closed-form as the

above methods but by a ﬁnite set of samples (Doucet

et al., 2001). SMC methods are capable to deal with

nonlinear dynamics and non-Gaussian noise. SMC

method are also available for multi-target tracking,

but a large number of samples is necessary to avoid

overly sparse representation of the high dimensional

joint distribution of targets’ states. Numerous exten-

sions exist, e.g. to achieve IMM- and JPDA-like be-

havior (Blom and Bloem, 2003). Monte Carlo meth-

ods avoid the combinatorial problem of data associ-

ation by sampling of association variables. In (Hue

et al., 2002) the joint distribution of state and asso-

ciation variables is completely approximated by sam-

ples. To reduce the number of samples (S

arkk

a et al.,

2007) introduce Rao-Blackwellization of state vari-

ables, i.e. their distributions are kept in analytical

form while associations are still sampled.

Targets and observations often are assumed to be

statistically independent for efﬁciency reasons. It thus

is reasonable to model dependencies of targets only

locally. In (Khan et al., 2005) the joint distribution

of otherwise independent targets is penalized using

Markov random ﬁelds, if targets approach. In (Smal

et al., 2008) state samples are re-clustered to recon-

sider their afﬁliation to targets in microscopy images.

3 MULTI-TARGET TRACKING

AND THE RBMCDA

FRAMEWORK

The scenario considered for tracking of multiple tar-

gets is as follows: At each time step multiple targets

might be present and their unknown number may vary

with time because targets may (dis-)appear in the ﬁeld

of view. Individual targets are identiﬁed by unique

target IDs which are preserved across time steps. The

set of target IDs of existing targets at time t is de-

noted with N

. The state of the target with ID n is

denoted by x

. Let X

:= {x

}

n∈N

describe the set

of the states of the N

= |X

| = |N

| existing targets

at time t. Targets may “die” if they disappear and

thus their ID and state have to be removed from the

set of states and IDs. On the other hand new targets

may be “born” in which case their state and ID must

be inserted into the corresponding sets. With X

−

, N

−

and N

−

we refer to the above sets and variables be-

fore target birth and death at time t and with X

, N

and N

thereafter. Further on, multiple observations

are detected at each time step, comprising observa-

tions from existing and newborn targets, but also false

alarms. False alarms are referred to as clutter obser-

vations and do not originate from targets but are due

to noise. Furthermore, a target may not be detected in

all time steps. Let

:= (

, . . . ,

) be the vector of

observations at time t, with an arbitrary but ﬁxed

order. Note that N

and

are implicitly deﬁned by

the number of elements in X

and

respectively. Let

:= (c

, . . . , c

) denote a vector of hidden associ-

ation variables, where c

describes association of the

m-th observation

to a target n (c

= n, n ∈ N

)

or clutter (c

= 0). Note that realizations of random

variables are marked with ˜.

This multi target tracking problem can be solved

by the iterative estimation of the joint distribution

of targets’ states analogous to the predict&update-

scheme well known from single-target tracking (see

e.g. (Doucet et al., 2001)). In fact, for the multi-target

case we seek for an estimate of the joint distribution of

targets’ states and association variables given all ob-

servations up to the current time t as stated in Eq. (1).

The factorization in Eq. (1) demonstrates the com-

binatorial problem inherent to multi-target tracking:

The ﬁrst term corresponds to the estimate of targets’

state distribution for ﬁxed association variables, i.e.

we know which observation to use for updating which

target. However, estimation of the complete distri-

bution is prohibitive for practical application because

the number of valid associations C

1:t

to be considered

explodes with time and increasing number of targets

VISAPP 2012 - International Conference on Computer Vision Theory and Applications

388

and observations.

p(X

, C

1:t

) = p(X

1:t

, C

1:t

) ·p(C

1:t

) (1)

For this reason (S

arkk

a et al., 2007) propose to

approximate Eq. (1) by a ﬁnite set of samples

{

(s)

1:t

}

s=1,...,S

in the space of association variables

whereas the distributions of targets’ states are main-

tained in analytical form, e.g. by Gaussian distribu-

tions.

p(X

, C

1:t

) ≈

∑

s=1

p(X

1:t

(s)

1:t

) ·δ(C

−

(s)

), (2)

with

(s)

sampled from p(C

1:t

(s)

1:t−1

) (3)

For each sample

(s)

1:t

a separate joint distribution

of its target states p(X

1:t

(s)

1:t

) is maintained.

This combination of sampled associations and corre-

sponding state distributions will be referred to as a

RBMCDA sample. Estimating the distribution of tar-

gets’ states follows the predict&update-scheme men-

tioned above with an intermediate sampling step to

yield realizations

(s)

before update.

According to (S

arkk

a et al., 2007) the individual

association variables in

(s)

are sampled sequentially

by decomposition of Eq. (3) using Bayes’ rule result-

ing in Eq. (4). The term for the m-th association vari-

able in Eq. (4) is further decomposed (Eq. (5)) into

the predicted likelihood of the corresponding obser-

vation for the given association and the data associa-

tion prior, which does not include information about

the actual observation.

p(C

1:t

, C

1:t−1

) = p(c

1:t

, C

1:t−1

)

∏

m=2

p(c

1:m−1

1:t

, C

1:t−1

) (4)

p(c

1:m−1

1:t

, C

1:t−1

)

∝p(

1:t−1

, C

1:t−1

)p(c

1:m−1

, N

−

) (5)

The formulation in Eq. (5) requires the assumption

that c

is independent from observations that are not

yet associated, i.e.

m+1:

. This assumption is al-

lowed by the scenario considered in (S

arkk

a et al.,

2007) and exploited to achieve efﬁcient sampling of

association variables in Sec. 4, but available informa-

tion clearly remains unused in our case.

In (S

arkk

a et al., 2007) only an informal descrip-

tion is provided of how to compute data association

prior probabilities. Therefore, we present a rigorous

formulation of a data association prior in 4.1 and

an efﬁcient sampling algorithm using this prior in 4.2.

4 PARAMETRIC DATA

ASSOCIATION PRIOR

4.1 Model of the Data Association Prior

In the following we derive the data association prior

from the joint distribution of C

when only the num-

ber of existing targets N

−

and their IDs N

−

is known.

This joint distribution models the formation of a set

of M

observations composed of k

observations from

the N

−

existing targets, b

observations from new-

born targets and u

= M

−k

−b

clutter observa-

tions. Existing targets are assumed to be detected

with probability P

. The number of newborn tar-

gets, which are all presume to be detected, is assumed

to be distributed according to ν(b

) and the number

of clutter observations according to the distribution

µ(u

). The probability of a valid set of association

variables C

is composed of the probability to receive

= k

+ b

+ u

observations and the probability

of the speciﬁc conﬁguration of association variables,

given that the number of existing targets is known.

The joint data association prior (6) is proportional to

this distribution for a ﬁxed number of observations

. It is deﬁned to be 0 if the cardinality of C

different from M

. Otherwise it is proportional to the

described distribution:

p(C

, N

−

) ∝ p(C

−

) =





−





−k





−



(1 −P

)

−

−k

ν(b

) µ(u

) (6)

The fraction in Eq. (6) accounts for the number of

possible conﬁgurations of C

for M

= k

+ b

+ u

observations. k

out of M

observations originate

from existing targets with





! different alterna-

tives to assign target IDs. A total of b

of the re-

maining (M

−k

) observations stem from newborn

targets with



−k



alternatives to choose b

observa-

tions. u

= M

−k

−b

observations are due to clutter.

The data association prior of an association c

re-

quired in Eq. (5) is proportional to Eq. (6) marginal-

ized over all valid conﬁgurations c

m+1:M

given ˜c

1:m−1

denoted by {c

m+1:M

|˜c

1:m−1

p(c

|˜c

1:m−1

, M

, N

−

) =

p( ˜c

1:m−1

, c

, N

−

)

p( ˜c

1:m−1

, N

−

)

∝

∑

m+1:M

|˜c

1:m−1

}

p( ˜c

1:m−1

, c

m:M

−

) (7)

=: q(c

|˜c

1:m−1

, M

, N

−

) (8)

PARAMETRIC DATA ASSOCIATION PRIOR FOR MULTI-TARGET TRACKING BASED ON

RAO-BLACKWELLIZED MONTE CARLO DATA ASSOCIATION

389

4.2 Sampling Algorithm

4.2.1 Computation of Data Association Prior

Three different cases have to be considered. The ﬁrst

case is the association of an observation to clutter, i.e.

= 0. The second case is the association to one of

the existing targets that is c

= n

exist

. Association to a

newborn target as the third case is given by c

= n

new

The marginalization Eq. (7) can be expressed via a

sum over possible numbers of associations to existing

targets k in C

given ˜c

1:m−1

are already assigned (see

e.g. Eq (13)). Conditioned on k, an inner sum runs

over the possible number of newborn target associa-

tions b in C

given ˜c

1:m−1

. All conﬁgurations of C

with identical k and b have the same probability, c.f.

Eq. (6).

For brevity of notation we use additional deﬁni-

tions in the following. M := M

denotes the number

of current observations and N := N

−

= N

t−1

is the

number of targets that already exist before time step t.

Case 1: Association to Clutter (c

= 0). We con-

sider all feasible numbers k of associations to exist-

ing targets in C

given ˜c

1:m−1

, if the m-th observa-

tion is associated to clutter. The minimum number is

m−1

, which is deﬁned as the number of associations

to existing targets present in ˜c

1:m−1

. The maximum

number equals k

max

:= min(N, M −m + k

m−1

) since a

maximum of M −m additional associations are possi-

ble in ˜c

m+1:M

and we have a total of N existing target

available. Likewise we consider the number of asso-

ciations to newborn targets in C

given ˜c

1:m−1

. The

minimum number corresponds to the number of pre-

vious associations to newborn targets in ˜c

1:m−1

de-

noted by b

m−1

. The maximum number is b

max

M −m −k + k

m−1

as a maximum of M −m−k

observations may additionally be associated to clut-

ter. The sum over valid conﬁgurations of associa-

tions c

m+1:M

given k

m−1

, b

m−1

, ˜c

1:m−1

and c

= 0 re-

quired in Eq. (7) can be expressed by the number of

these conﬁgurations and is given in Eq. (13). There

are



M−m

k−k

m−1



alternatives to select (k −k

m−1

) out of

(M −m) observations for additional association with

existing targets.



N−k

m−1

k−k

m−1



(k −k

m−1

)! alternatives to

associate (k −k

m−1

) of the (N −k

m−1

) non-associated

existing targets exclusively and



M−m−k+k

m−1

b−b

m−1



alter-

natives to associate (b −b

m−1

) out of the remaining

(M −m −k + k

m−1

) observations to newborn targets.

To highlight relations between the different asso-

ciation cases the following deﬁnitions will be used:

d(m, k

m−1

;M, N) :=

(M −m)! (N −k

m−1

M! N!

(9)

f (k, k

m−1

;N) :=

k! ·





(1 −P

)

N−k

(k −k

m−1

(10)

g(b, b

m−1

) :=

b! ·ν(b)

(b −b

m−1

(11)

h(m, r, r

m−1

;M) :=

(M −r)! ·µ(M −r)

(M −r −m + r

m−1

(12)

Here we denote r = k + b and r

m−1

= k

m−1

+ b

m−1

Their properties are examined in detail in 4.2.2 and

are utilized for efﬁcient computation.

Case 2: Association to Existing Target (c

= n

exist

The marginal in Eq. (7) for association to an existing

target differs from the clutter case in the sum limits

for k and is given in Eq. (14). The limits of the clutter

case have to be increased by one, because the number

of associations to existing targets in c

1:m

is (k

m−1

+1)

in this case. This aspect has also to be considered in

the number of possible associations in c

m+1:M

Case 3: Association to Newborn Target (c

new

). For association to a newborn target the sum

limits for b have to be adjusted with regard to the clut-

ter case (see Eq. (15)). Similar to the adjustment of

the k-limits in case 1, the limits of b have to be in-

creased by one, because the number of associations

to newborn targets in c

1:m−1

and c

now equals to

m−1

+ 1).

4.2.2 Iterative Properties of Association Prior

The functions d, f , g and h introduced in equations

(9) to (12) can be computed iteratively and help to

understand the iterative computation of terms used to

sequentially sample the data association prior.

If m is incremented d(m, k

m−1

;M, N) can be easily

updated regardless whether k

= k

m−1

+1 (Case 2) or

= k

m−1

(Case 1 and Case 3):

d(0, 0; M, N) =1 (16)

d(m + 1, k

m−1

;M, N) =

M −m

d(m, k

m−1

;M, N)

d(m, k

m−1

+ 1;M, N) =

N −k

m−1

d(m, k

m−1

;M, N)

The values of f (k, k

m−1

;N) only change if an associ-

ation to an existing target is sampled (Case 2):

f (k, 0; N) =





(1 −P

)

N−k

(17)

f (k, k

m−1

+ 1;N) =(k −k

m−1

) f (k, k

m−1

;N)

Likewise the function g(b, b

m−1

) changes if associa-

tion to a newborn target occurs (Case 3).

g(b, 0) =ν(b) (18)

g(b, b

m−1

+ 1) =(b −b

m−1

) g(b, b

m−1

)

VISAPP 2012 - International Conference on Computer Vision Theory and Applications

390

q(c

= 0|˜c

1:m−1

, M, N

−

) =

max

∑

k=k

m−1

max

∑

b=b

m−1



M−m

k−k

m−1



N−k

m−1

k−k

m−1



(k −k

m−1



M−m−k+k

m−1

b−b

m−1



p(C

|N)

max

∑

k=k

m−1

max

∑

b=b

m−1



M−m

k−k

m−1



N−k

m−1

k−k

m−1



(k −k

m−1



M−m−k+k

m−1

b−b

m−1











M−k







(1 −P

)

N−k

ν(b)µ(M −k −b)

(M −m)! (N −k

m−1

M! N!

max

∑

k=k

m−1

(k −k

m−1





(1 −P

)

N−k



max

∑

b=b

m−1

(b −b

m−1

ν(b)

(M −k −b)!

(M −k −b −m + k

m−1

+ b

m−1

µ(M −k −b)



= d(m, k

m−1

;M, N)

max

∑

k=k

m−1

f (k, k

m−1

;N)

max

∑

b=b

m−1

g(b, b

m−1

) h(m, k + b, k

m−1

+ b

m−1

;M) (13)

q(c

= n

exist

|˜c

1:m−1

, M, N

−

) =

max

∑

k=k

m−1

max

∑

b=b

m−1



M−m

k−k

m−1

−1



N−k

m−1

−1

k−k

m−1

−1



(k −k

m−1

−1)!



M−m−k+k

m−1

b−b

m−1











M−k







(1 −P

)

N−k

ν(b)µ(M −k −b)

= d(m, k

m−1

+ 1; M, N)

max

∑

k=k

m−1

f (k, k

m−1

+ 1; N)

max

∑

b=b

m−1

g(b, b

m−1

)h(m, k + b, k

m−1

+ 1 + b

m−1

;M) (14)

q(c

= n

new

|˜c

1:m−1

, M, N

−

) =

max

∑

k=k

m−1

max

∑

b=b

m−1



M−m

k−k

m−1



N−k

m−1

k−k

m−1



(k −k

m−1



M−m−k+k

m−1

b−b

m−1

−1











M−k







(1 −P

)

N−k

ν(b) ·µ(M −k −b)

=d(m, k

m−1

;M, N)

max

∑

k=k

m−1

f (k, k

m−1

;N)

max

∑

b=b

m−1

g(b, b

m−1

+ 1)h(m, k + b, k

m−1

+ b

m−1

+ 1; M) (15)

If both m and r

m−1

are incremented simultaneously

(Cases 2 and 3) h remains unchanged. If m is incre-

mented but the total number target associations r

m−1

remains unchanged (Case 1) the last recursion is used:

h(1, r, 0; M) =(M −r) µ(M −r) (19)

h(1, r, 1; M) =µ(M −r)

h(m + 1, r, r

m−1

+ 1;M) =h(m, r, r

m−1

;M)

h(m + 1, r, r

m−1

;M) =(M −r −m + r

m−1

)

×h(m, r, r

m−1

;M)

4.2.3 Sampling Algorithm

To sample a complete vector of associations

(s)

for time step t given

(s)

1:t−1

, the individual real-

ization

(s)

˜c

are sampled sequentially for m = 1 to

m = M. An iteration consists of two parts. First the

probabilities for sampling (see Eq. (5)) are computed.

Subsequently these probabilities are used to sample

the m-th association

(s)

˜c

. In the second part vari-

ables and arrays representing functions d, f , g and h

are updated if necessary using the recursions derived

in the last subsection. This sampling is performed

independently for each RBMCDA sample

(s)

turn. For clarity the RBMCDA sample index

(s)

omitted in the following description of the algorithm.

INITIALIZE

d(1, 0; M, N), f (k, 0; M, N), f (k, 1; M, N)

g(b, 0), g(b, 1), h(1, r, 0;M), h(1, r, 1;M)

m−1

← 0 , b

m−1

← 0

← 1 // p(C

1:t

1:t−1

)

C ←{} //

// SEQUENTIAL SAMPLING of C

for m = 1 →M do

// Part I: Compute probabilities and sample ˜c

// compute data association prior probabilities

clutter

← q(c

= 0|˜c

1:m−1

, M, N) // Eq. (13)

PARAMETRIC DATA ASSOCIATION PRIOR FOR MULTI-TARGET TRACKING BASED ON

RAO-BLACKWELLIZED MONTE CARLO DATA ASSOCIATION

391

(a) N

init

= 50. (b) N

init

= 100. (c) N

init

= 200.

Figure 1: Images of the generated sequences at t = 10: Observations are drawn as spots together with their current trajectories.

exist

← q(c

= n

exist

|˜c

1:m−1

, M, N) // Eq. (14)

newborn

← q(c

= n

new

|˜c

1:m−1

, M, N) // Eq. (15)

// compute probabilities for association to clutter

ˆp(c

= 0) ← q

clutter

·p

clutter

)

// ... to each target not associated in ˜c

1:m−1

for all n ∈ N

−

\˜c

1:m−1

ˆp(c

= n) ← q

exist

×p(z

= n, c

1:m−1

, Z

1:t

, C

1:t−1

)

end for

// ... to newborn target

ˆp(c

= n

new

) ← q

newborn

·p

newborn

)

normalize ˆp(c

) → p(c

)

// sample association ˜c

˜c

∼ p(c

)

C ←(

C, ˜c

)

← p

·p( ˜c

)

// Part II: Update variables and arrays

// update d for m ← m + 1

d(m + 1, k

m−1

;M, N) ← Eq. (16)

if ˜c

= 0 then

// update h for m ← m + 1

h(m + 1, r, k

m−1

+ b

m−1

;M) ← Eq. (19)

h(m + 1, r, k

m−1

+ b

m−1

+ 1; M) ← Eq. (19)

else if ˜c

∈ N

−

then

// update f and d for k

m−1

← k

m−1

+ 1

// f (k, k

m−1

+ 1; M, N) remains unaffected

f (k, k

m−1

+ 2; M, N) according to (17)

d(m + 1, k

m−1

+ 1; M, N) ← according to (16)

m−1

← k

m−1

+ 1

else

// update g for b

m−1

← b

m−1

+ 1

// g(b, b

m−1

+ 1) remains unaffected

g(b, b

m−1

+ 2) according to (18)

m−1

← b

m−1

+ 1

end if

end for

return (

C, p

)

5 RESULTS

Model Speciﬁcation. The state of a target com-

prises its current position (x

, y

) and size, deﬁned

as the squareroot of the occupied area, as well as

its position in the previous time step, i.e. x

, y

, x

t−1

, y

t−1

√

area

)

. Observations comprise

position and size only, i.e. z

:= (x

, y

√

area

)

. The

observation of a target is modeled to obtain the current

location and size of the target perturbed with additive

Gaussian noise. The measurement noise covariance

matrix R is deﬁned as a diagonal matrix with vari-

ances for coordinates and size. The dynamic behavior

of targets is modeled by switching linear models with

Gaussian noise. The corresponding noise covariance

matrices Q again are deﬁned as diagonal matrix with

the variances for current and previous coordinates as

well as for size. Two dynamic models are employed:

A random walk model (rw) accounts for stopping and

change of direction by no deterministic change in po-

sition and size of a target, but Gaussian noise only. A

second model for directional motion (ﬂe) extrapolates

the next position from the current and last position

with additional Gaussian noise. See (Genovesio and

Olivo-Marin, 2008) for more details. Switching of

dynamic models between time steps is described by

the probabilities P

rw→rw

to stay in the rw-model and

rw→f le

, P

f le→rw

and P

f le→f le

deﬁned analogously.

We integrate the IMM ﬁlter (Bar-Shalom and Blair,

2000) into the RBMCDA framework, which estimates

the state distribution of a target with changing dynam-

ics by a Gaussian mixture.

Poisson distributions ν(b; λ

) and µ(u; λ

) ac-

count for the number of observations from newborn

targets and clutter observations respectively, con-

trolled by parameters λ

and λ

. The death of targets

follows an exponential distribution of the time to last

association with parameter λ

. The proposed model is

generative and thus may be used to generate synthetic

data.

Synthetic Data Generation. Three time sequences

of 50 time steps were generated for evaluation. Each

sequence was generated for a different number of ini-

tial targets in an image of ﬁxed size to investigate the

behavior of the algorithm regarding different numbers

of targets and different spatial density of observations.

The sequences are referred to by the number of initial

VISAPP 2012 - International Conference on Computer Vision Theory and Applications

392

Table 1: Parameters used for data generation and tracking.

init

= 50 N

init

= 100 N

init

= 200

0.6 1.4 2.25





5 0 0 0 0

0 5 0 0 0

0 0 1.67 0 0

0 0 0 1.67 0

0 0 0 0 1.1





5 0.5 0.97



5 0 0

0 5 0

0 0 1.1



rw→rw

rw→f le

f le→rw

f le→f le

0.7 0.3 0.5 0.5

targets N

init

∈ {50, 100, 200}. Fig. 1 shows images of

the different sequences at t = 10 to give an impres-

sion of how densely observations are situated. All pa-

rameters are ﬁxed for generation of the three differ-

ent sequences except for λ

, which was adjusted for

each sequence to achieve a nearly balanced number of

deaths and births of targets. The parameters for data

generation are presented in Tab. 1. True associations

of all observations in this synthetic data are known

and used for evaluation of tracking results. Tracking

was performed with a number of RBMCDA samples

S ∈ {10, 100, 1000} for each sequence using the cor-

responding parameters from data generation.

Performance Evaluation. Evaluation of

RBMCDA tracking performance is not straight-

forward as target IDs of different RBMCDA samples

do not relate to each other nor to IDs of ground-truth

data. Thus trajectories cannot directly be related

to each other using target IDs or compared to

ground-truth. For the same reason state distributions

estimated within the RBMCDA samples cannot eas-

ily be combined. We therefore consider track graphs

to evaluate tracking results. The nodes of a track

graph are constituted by all observations of a time

sequence. Edges are deﬁned as follows. For each

target we consider all observations associated to this

target. These observations are ordered with respect to

time and neighboring observations in this sequence

are connected by an edge. As a consequence all ob-

servations associated to clutter are not incident to any

edge. In the same way the ground-truth track graph

is constructed using true associations available from

synthetic data. The graph of each RBMCDA sample

is compared to this ground-truth graph as follows.

Edges that are common to both graphs are considered

as true positives (TP), edges that appear only in the

ground-truth graph are counted as false negatives

(FN) and edges only in the graph from the current

sample are regarded as false positives (FP). An

example is given in Fig. 2, where in the ground-truth

track graph observations are horizontally connected,

t=1

t=2

t=3

t=4

t=5 t=6

t=7

1|1

2|1

3|1

4|1 6|1

7|1

2|2

4|2

5|2 7|2

Figure 2: Evaluation of track graphs: nodes represent obser-

vations, solid edges ground-truth trajectories, dashed edges

tracking results. TP are edges present in both graphs, FN

only in ground-truth and FP only in tracking results.

Precision for tracking with S samples

Tracker parameter

Generator parameter

S=10

S=100

S=1000

init

=200

init

=100

init

=50

0.82

0.919

0.95

0.834

0.928

0.959

0.841

0.93

0.961

0.0

0.2

0.4

0.6

0.8

1.0

Precision Tracksegments

(a) Precision P .

Recall for tracking with S samples

Tracker parameter

Generator parameter

S=10

S=100

S=1000

init

=200

init

=100

init

=50

0.762

0.855

0.824

0.782

0.879

0.807

0.899

0.886

0.0

0.2

0.4

0.6

0.8

1.0

Recall Tracksegments

(b) Recall R .

Figure 3: Tracking performance for synthetic data N

init

tracked with different number of RBMCDA samples S.

but observations 2|1 and 2|2 are interchanged in the

track graph of the RBMCDA sample. Furthermore

observations 4|2 and 6|1 are wrongly associated to

clutter. To evaluate the performance of tracking we

use precision P :=

TP+FP

and recall R :=

TP+FN

Note that the track graph is only assessed locally

with respect to time, which may underestimate the

performance. E.g., in Fig. 2 observation 6|1 is

associated to clutter by the tracker resulting in two

FN and one FP, while the association of 4|1 to 7|1 is

obviously correct. Still we prefer this measure for its

clear deﬁnition.

The precision P of tracking results in the different

experiments is presented in Fig. 3(a) and recall R

in Fig. 3(b). Both measures increase for any of the

sequences, when the number of RBMCDA samples is

increased. On the other hand P and R decrease, when

more targets are presents in the domain.

Computation Times. Tracking was conducted on

AMD Opteron 848 processors with clock speed 2.2

GHz using a Java implementation on Sun’s 64-Bit vir-

tual machine. Probabilities are represented by their

logarithms to avoid numerical problems. Their sum-

mation thus is expensive and required extensively by

Eqs. (13), (14) and (15). Computation times are mea-

sured in UNIX’ user CPU time and presented in Tab.

2 for the different sequences and sample sizes.

PARAMETRIC DATA ASSOCIATION PRIOR FOR MULTI-TARGET TRACKING BASED ON

RAO-BLACKWELLIZED MONTE CARLO DATA ASSOCIATION

393

Table 2: Computation time (seconds) of tracking sequences

init

with different number of RBMCDA samples S.

S = 10 S = 100 S = 1000

init

= 50 21.1 167.3 2,083.1

init

= 100 89.8 886.3 11,021.9

init

= 200 384.8 3,714.1 38,478.0

6 DISCUSSION

The proposed method is evaluated for tracking per-

formance in order to prove its concept. Sequence

init

= 50 was generated to pose a tracking problem

of relatively low complexity (see Fig. 1(a)). Targets

are not located densely and noise in position and size

is reasonable. The performance is given in Fig. 2

for varying sample size S. With regard to the pes-

simistic characteristics of the employed measures, the

resulting performance indicates expectedly good per-

formance for tracking with S = 10 RBMCDA sam-

ples. In general R is smaller than P because, e.g.,

wrong associations to clutter have no impact on P ,

however on R . If the number of RBMCDA samples is

increased performance improves with maximum pre-

cision of 0.961 and recall of 0.886 for S = 1000.

The tracking problem gets more difﬁcult as the

number of targets increases in sequence N

init

= 100

(see Fig. 1(b)). Wrong associations become more

likely when targets and their observations approach

and thus both P and R decrease in general, but are

still very satisfactory (see Figs. 3(a) and 3(b)). Again,

the performance improves with sample size S. Maxi-

mum precision is 0.93 with a recall of 0.899.

Even more targets are generated in sequence

init

= 200, posing a very challenging tracking prob-

lem (see Fig. 1(c)). Here, a maximum precision

of 0.841 is achieved where the maximum recall

is 0.807. This is well acceptable taking the complex-

ity of the data into account. It could be speculated

that the performance may be further increased using a

larger set of RBMCDA samples.

In summary, tracking performance shows that the

proposed method is able to successfully track large

numbers of targets in demanding data.

For practical application we recommend to repre-

sent probabilities by their logarithms to track larger

numbers of targets, because sampling would fail as

the small probabilities in Eq. (5) cannot be repre-

sented by double precision any more. Calculation

of the data association prior becomes more expensive

due to the summation of probabilities, but the numer-

ical stability allows to successfully track at least 3600

targets. This is in contrast to about 70 targets for a

conventional representation of probabilities.

7 CONCLUSIONS

In this work we propose a parametric data associa-

tion prior for use with RBMCDA presented in (S

arkk

et al., 2007) to track a varying number of targets. This

prior models the formation of observations from ex-

isting and newborn targets as well as clutter observa-

tions. We developed an efﬁcient algorithm to sample

associations using this prior. In a proof of concept we

show that this sampling procedure allows to success-

fully sample associations in the presence of hundreds

of targets. Computation times are moderate using a

sample size of 1000 and as much as 200 targets on

average. For synthetic data of demanding complex-

ity the performance of sampled associations is well

acceptable. We integrated the IMM ﬁlter to enable

tracking of targets with changing dynamics for an ap-

plication to microscopy image analysis in mind.

REFERENCES

Bar-Shalom, Y. and Blair, W. (2000). Multitarget-

multisensor tracking: Applications and advances. Vol

III. Artech House.

Bar-Shalom, Y. and Fortmann, T. (1987). Tracking and data

association. Academic Press Professional.

Blom, H. and Bloem, E. (2003). Joint IMMPDA particle

ﬁlter. In Information Fusion, Proc 6th International

Conference.

Doucet, A., de Freitas, N., and Gordon, N. (2001). An intro-

duction to sequential monte carlo methods. Sequential

Monte Carlo methods in practice, pages 3–14.

Genovesio, A. and Olivo-Marin, J. (2008). Particle tracking

in 3D+t biological imaging. In Rittscher, J., editor,

Microscopic Image Analysis for Life Science Applica-

tions, pages 223–282. Artech House.

Hue, C., Le Cadre, J., and P

erez, P. (2002). Sequen-

tial monte carlo methods for multiple target track-

ing and data fusion. Signal Processing, IEEE Trans,

50(2):309–325.

Khan, Z., Balch, T., and Dellaert, F. (2005). MCMC-based

particle ﬁltering for tracking a variable number of in-

teracting targets. Pattern Analysis and Machine Intel-

ligence, IEEE Trans, 27:2005.

Musicki, D. and Evans, R. (2002). Joint integrated prob-

abilistic data association-JIPDA. In Information Fu-

sion, 2002. Proc 5th International Conference, vol-

ume 2, pages 1120–1125. IEEE.

Reid, D. (1979). An algorithm for tracking multiple targets.

Automatic Control, IEEE Trans, 24(6):843–854.

arkk

a, S., Vehtari, A., and Lampinen, J. (2007). Rao-

blackwellized particle ﬁlter for multiple target track-

ing. Information Fusion, 8(1):2–15.

Smal, I. et al. (2008). Multiple object tracking in molecu-

lar bioimaging by rao-blackwellized marginal particle

ﬁltering. Medical image analysis, 12(6):764–777.

VISAPP 2012 - International Conference on Computer Vision Theory and Applications

394