Bayesian Smoothing for the Extended Object Random Matrix Model

TL;DR

This paper develops closed-form Bayesian smoothing formulas for two variants of the extended object random matrix model, enhancing tracking accuracy in extended object tracking applications.

Contribution

It introduces novel Bayesian smoothing expressions for both conditional and factorized models in the extended object random matrix framework.

Findings

Bayesian smoothers improve tracking performance.

Closed-form solutions are derived for both model types.

Simulation results compare the effectiveness of different smoothers.

Abstract

The random matrix model is popular in extended object tracking, due to its relative simplicity and versatility. In this model, the extended object state consists of a kinematic vector for the position and motion parameters (velocity, etc), and an extent matrix. Two versions of the model can be found in literature, one where the state density is modelled by a conditional density, and one where the state density is modelled by a factorized density. In this paper, we present closed form Bayesian smoothing expression for both the conditional and the factorised model. In a simulation study, we compare the performance of different versions of the smoother.

Tables

Table 1. TABLE II: Conditional model: prediction

	$m_{k+1|k}$	$=\left(F_k\otimes {𝐈}_d\right)m_{k|k}$
	$P_{k+1|k}$	$=F_k{P}_{k|k}{F}_k^T+D$
	$v_{k+1|k}$	$=d+1+{\left(1+\frac{v_{k|k}-2d-2}{n}\right)}^{-1}\left(v_{k|k}-d-1\right)$
	$V_{k+1|k}$	$={\left(1+\frac{v_{k|k}-d-1}{n-d-1}\right)}^{-1}A{V}_{k|k}{A}^T$

Table 2. TABLE III: Conditional model: update

Table 3. TABLE IV: Conditional model: smoothing

Table 4. TABLE V: Factorized model: prediction

Table 5. TABLE VI: Factorized model: update

Table 6. TABLE VII: Factorized model: smoothing

Equations

$$p({\xi }_{k+1}|{\mathbf{Z}}_{1:k})$$

$$p({\xi }_{k+1}|{\mathbf{Z}}_{1:k+1})$$

$$p({\xi }_k|{\mathbf{Z}}_{1:K})$$

$$p({\xi }_k|{\mathbf{Z}}_{1:\ell })=$$

$$=$$

$$\times {\mathcal{I}\mathcal{W}}_d\left(X_k;v_{k|\ell },V_{k|\ell }\right),$$

$$p({\xi }_k|{\mathbf{Z}}_{1:\ell })=$$

$$=$$

$$p({\xi }_{k+1}|{\xi }_k)=$$

$$=$$

$$p({\mathbf{Z}}_k|{\xi }_k)\propto \prod _{\mathbf{z}\in {\mathbf{Z}}_k}p(\mathbf{z}|{\mathbf{x}}_k,X_k)$$

$$p({\xi }_k|{\mathbf{Z}}_{1:K})=$$

$$=$$

$$\times {\mathcal{I}\mathcal{W}}_d\left(X_k;v_{k|K},V_{k|K}\right).$$

$$p\left(X_{k+1}|{\mathbf{x}}_k,X_k\right)=p\left(X_{k+1}|X_k\right).$$

$$p({\mathbf{x}}_k|X_k)$$

$$p({\mathbf{x}}_{k+1}|X_{k+1})$$

$$p\left({\mathbf{x}}_{k+1}|X_{k+1},{\mathbf{x}}_k\right)$$

$$p({\xi }_{k+1}|{\xi }_k)\approx$$

$$=$$

$$\times {\mathcal{W}}_d\left(X_{k+1};n_k,\frac{X_k}{n_k}\right)$$

$$p({\xi }_{k+1}|{\xi }_k)\approx$$

$$\times {\mathcal{W}}_d\left(X_{k+1};n_k,\frac{A{X}_k{A}^{\mathrm{T}}}{n_k}\right)$$

$$p(\mathbf{z}|{\mathbf{x}}_k,X_k)=\mathcal{N}\left(\mathbf{z};\left(H_k\otimes {\mathbf{I}}_d\right){\mathbf{x}}_k,X_k\right),$$

$$p({\xi }_k|{\mathbf{Z}}_{1:K})=$$

$$=$$

$$\times {\mathcal{I}\mathcal{W}}_d\left(X_k;v_{k|K},V_{k|K}\right).$$

$$p\left({\mathbf{x}}_{k+1}|X_{k+1},{\mathbf{x}}_k\right)=p\left({\mathbf{x}}_{k+1}|{\mathbf{x}}_k\right)$$

$$p({\xi }_{k+1}|{\xi }_k)=$$

$$\times {\mathcal{W}}_d\left(X_{k+1};n_k,\frac{M({\mathbf{x}}_k)X_k{M}^{\mathrm{T}}({\mathbf{x}}_k)}{n_k}\right)$$

$$p(\mathbf{z}|{\mathbf{x}}_k,X_k)=\mathcal{N}\left(\mathbf{z};{\tilde{H}}_k{\mathbf{x}}_k,\rho X_k+R_k\right),$$

$$\displaystyle\begin{array}[]{rcl}m_{k+1|k}&=&f_{k}(m_{k|k})\\ P_{k+1|k}&=&\tilde{F}_{k}P_{k|k}\tilde{F}_{k}^{\scriptscriptstyle\mathrm{T}}+Q\\ v_{k+1|k}&=&d+1+\eta^{-1}(v_{k|k}-d-1)\\ V_{k+1|k}&=&\eta^{-1}\left(1-\frac{d+1}{s}\right)\left(1-\frac{d+1}{n}\right)C_{2}\\ \eta&=&1+(v_{k|k}-2d-2)\left(\frac{1}{s}+\frac{1}{n}-\frac{d+1}{ns}\right)\\ s&=&\frac{d+1}{d}\operatorname{Tr}\left\{C_{1}C_{2}\left(C_{1}C_{2}-\mathbf{I}_{d}\right)^{-1}\right\}\\ \tilde{F}_{k}&=&\left.\nabla_{\mathbf{x}}f_{k}(\mathbf{x})\right|_{\mathbf{x}=m_{k|k}}\\ C_{1}&=&\mathbb{E}_{k|k}\left[\left(M_{\mathbf{x}}V_{k|k}M_{\mathbf{x}}^{\scriptscriptstyle\mathrm{T}}\right)^{-1}\right]\\ C_{2}&=&\mathbb{E}_{k|k}\left[M_{\mathbf{x}}V_{k|k}M_{\mathbf{x}}^{\scriptscriptstyle\mathrm{T}}\right]\end{array}$$

$$\displaystyle\begin{array}[]{rcl}m_{k|K}&=&m_{k|k}+G_{k}\left(m_{k+1|K}-m_{k+1|k}\right)\\ P_{k|K}&=&P_{k|k}-G\left(P_{k+1|k}-P_{k+1|K}\right)G^{\scriptscriptstyle\mathrm{T}}\\ v_{k|K}&=&v_{k|k}+\eta_{2}^{-1}\left(g-\frac{2(d+1)^{2}}{h+d+1}\right)\\ V_{k|K}&=&V_{k|k}+\eta_{3}^{-1}C_{4}\\ G&=&P_{k|k}\tilde{F}_{k}^{\scriptscriptstyle\mathrm{T}}P_{k+1|k}^{-1}\\ W&=&V_{k+1|K}-V_{k+1|k}\\ w&=&v_{k+1|K}-v_{k+1|k}\\ g&=&\eta_{1}^{-1}\left(w-\frac{2(d+1)^{2}}{n}\right)\\ h&=&\frac{d+1}{d}\operatorname{Tr}\left\{C_{3}C_{4}\left(C_{3}C_{4}-\mathbf{I}_{d}\right)^{-1}\right\}\\ \eta_{1}&=&1+\frac{w-3(d+1)}{n}\\ \eta_{2}&=&1+\frac{g-3d-3}{h+d+1}\\ \eta_{3}&=&1+\frac{g-d-1}{h-d-1}\\ C_{3}&=&\mathbb{E}_{k|K}\left[\left(M_{\mathbf{x}}^{-1}W(M_{\mathbf{x}}^{-1})^{\scriptscriptstyle\mathrm{T}}\right)^{-1}\right]\\ &=&\mathbb{E}_{k|K}\left[M_{\mathbf{x}}^{\scriptscriptstyle\mathrm{T}}W^{-1}M_{\mathbf{x}}\right]\\ C_{4}&=&\mathbb{E}_{k|K}\left[M_{\mathbf{x}}^{-1}W(M_{\mathbf{x}}^{-1})^{\scriptscriptstyle\mathrm{T}}\right]\\ &=&\mathbb{E}_{k|K}\left[\left(M_{\mathbf{x}}^{\scriptscriptstyle\mathrm{T}}W^{-1}M_{\mathbf{x}}\right)^{-1}\right]\end{array}$$

$$C_1\approx$$

$$+\sum _{i=1}^{n_x}\sum _{j=1}^{n_x}{\frac{{\mathrm{d}}^2{\left(M_{\mathbf{x}}V{M}_{\mathbf{x}}^{\mathrm{T}}\right)}^{-1}}{\mathrm{d}{\mathbf{x}}^{[i]}\mathrm{d}{\mathbf{x}}^{[j]}}|}_{\mathbf{x}=m}{P}^{[i,j]}$$

$$\frac{\mathrm{d}{M}_{\mathbf{x}}V{M}_{\mathbf{x}}^{\mathrm{T}}}{\mathrm{d}{\mathbf{x}}^{[j]}}=$$

$$\frac{{\mathrm{d}}^2{M}_{\mathbf{x}}V{M}_{\mathbf{x}}^{\mathrm{T}}}{\mathrm{d}{\mathbf{x}}^{[i]}\mathrm{d}{\mathbf{x}}^{[j]}}=$$

$$+\frac{\mathrm{d}{M}_{\mathbf{x}}}{\mathrm{d}{\mathbf{x}}^{[i]}}V\frac{\mathrm{d}{M}_{\mathbf{x}}^{\mathrm{T}}}{\mathrm{d}{\mathbf{x}}^{[j]}}+M_{\mathbf{x}}V\frac{{\mathrm{d}}^2{M}_{\mathbf{x}}^{\mathrm{T}}}{\mathrm{d}{\mathbf{x}}^{[i]}\mathrm{d}{\mathbf{x}}^{[j]}}$$

$$\frac{{\mathrm{d}}^2{N}_{\mathbf{x}}^{-1}}{\mathrm{d}{\mathbf{x}}^{[i]}\mathrm{d}{\mathbf{x}}^{[j]}}=$$

Full text

I Introduction

Multiple Object Tracking (mot) denotes the process of successively determining the number and states of multiple dynamic objects based on noisy sensor measurements. Tracking is a key technology for many technical applications in areas such as robotics, surveillance, autonomous driving, automation, medicine, and sensor networks. Extended object tracking is defined as mot where each object generates multiple measurements per time step and the measurements are spatially structured on the object, see [1].

Extended object tracking is applicable in many different scenarios, e.g., environment perception for autonomous vehicles using camera, lidar and automotive radar. The multiple measurements per object and time step create a possiblity to estimate the object extent, in addition to the position and the kinematic properties such as velocity and heading. This estimation requires an object state space model, including modelling of the object dynamics and the measurement process. Extended obejct models include the Random Matrix model [2, 3], the Random Hypersurface model [4], and Gaussian Process models [5]. A comprehensive overview of extended object tracking can be found in [1].

In this paper we focus on the Random Matrix model, also known as the Gaussian inverse Wishart (giw) model. The random matrix model was originally proposed by Koch [2], and is an example of a spatial model. In this model the shape of the object is assumed to be elliptic. The ellipse shape is simple but still versatile, and the random matrix model has been integrated into many different multiple extended object tracking frameworks [6, 7, 8, 9, 10, 11, 12, 13, 14, 15]. Indeed, the random matrix model is applicable in many real scenarios, e.g., pedestrian tracking using video [6, 7, 8] or lidar [9] and tracking of boats and ships using marine radar [16, 17, 15, 18, 19, 20, 14].

The focus of this paper is on Bayesian smoothing for the random matrix model. A preliminary version of this work was presented in [21]. This paper is a significant extension of [21], and presents the following contributions:

• •

Closed form smoothing expressions for the conditional giw model from [2].

• •

Closed from smoothing expressions for the factorized giw model from [3].

• •

A simulation study that compares the derived smoothers to both prediction and filtering.

As a minor contribution, a closed form expression for the random matrix prediction from [22] is presented.

The rest of the paper is organized as follows. A problem formulation is given in the next section. In Section III a review of the random matrix model is given. Smoothing for the conditional giw model is presented in Section IV; smoothing for the factorised giw model is presented in Section V. Results from a simulation study are presented in Section VI. Concluding remarks are given in Section VII.

II Problem formulation

Let $\xi_{k}$ denote the extended object state at time $k$, let $\mathbf{Z}_{k}$ denote the set of measurements at time step $k$, and let $\mathbf{Z}_{1:k}$ denote the sets of measurements from time $1$ up to, and including, time $k$. Bayesian extended object filtering builds upon two steps, the Chapman-Kolmogorov prediction

[TABLE]

where $p(\xi_{k+1}|\xi_{k})$ is the transition density, and the Bayes update

[TABLE]

where $p(\mathbf{Z}_{1:k+1}|\xi_{k+1})$ is the measurement likelihood. The focus of this paper is on Bayesian extended object smoothing,

[TABLE]

where $K$ is the final time step. For the random matrix model, the Chapman-Kolmogorov prediction and Bayes update have been covered extensively in previous litterature, see, e.g., [2, 3, 22, 24] for the prediction, and, e.g., [2, 3, 25, 26, 27, 24] for the update. In this paper, we focus on Bayesian extended object smoothing. In previous literature, smoothing is only discussed briefly in [2, Sec. 3.F], and complete details are not given.

Bayesian filtering and smoothing for the random matrix model is an example of assumed density filtering: the functional form of the state density is to be preserved in the prediction and the update. It is therefore necessary that Bayesian smoothing also preserves the functional form of the extended object state density. Two different assumed state densities can be found in the literature: the conditional Gaussian inverse Wishart [2], and the factorized Gaussian inverse Wishart [3].

The problem considered in this paper is to use the Bayesian smoothing equation (3) to compute the smoothing giw parameters for both the conditional model and the factorized model.

III Review of random matrix model

In this section we give a brief review of the random matrix model; a longer review can be found in [1, Sec. 3.A].

In the random matrix model [2, 3], the extended object state is a tuple $\xi_{k}=\left(\mathbf{x}_{k},X_{k}\right)\in\mathbb{R}^{n_{x}}\times\mathbb{S}_{++}^{d}$. The vector $\mathbf{x}_{k}\in\mathbb{R}^{n_{x}}$ represents the object’s position and its motion properties, such as velocity, acceleration, and turn-rate. The matrix $X_{k}\in\mathbb{S}_{++}^{d}$ represents the object’s extent, where $d$ is the dimension of the object; $d=2$ for tracking with 2D position and $d=3$ for tracking with 3D position. The matrix $X_{k}$ is modelled as being symmetric and positive definite, which means that the object shape is approximated by an ellipse.

In the literature, there are two alternative models for the extended object state density, the conditional and the factorised. In the conditional model, first presented in [2], the following state density is used,

[TABLE]

where $m_{k|\ell}\in\mathbb{R}^{n_{x}}$, $P_{k|\ell}\in\mathbb{S}_{+}^{s}$, $v_{k|\ell}>2d$, $V_{k|\ell}\in\mathbb{S}_{++}^{d}$, and $s=\frac{n_{x}}{d}$. In this model, the random vector $\mathbf{x}$ consists of a $d$-dimensional spatial component (the position) and its derivatives (velocity, acceleration, etc.), see [2, Sec. 3]. Thus, $s-1=\frac{n_{x}}{d}-1$ describes up to which derivative the kinematics are described, see [2, Sec. 3].

In the factorised model, first presented in [3], the following state density is used,

[TABLE]

where $m_{k|\ell}\in\mathbb{R}^{n_{x}}$, $P_{k|\ell}\in\mathbb{S}_{+}^{n_{x}}$, $v_{k|\ell}>2d$, and $V_{k|\ell}\in\mathbb{S}_{++}^{d}$. In this model, the random vector $\mathbf{x}$ consists of a $d$-dimensional spatial component (the position) and additional motion parameters; note that, in contrast to the conditional model, here the motion parameters are not restricted to being derivatives of the spatial component, and non-linear dynamics can be modelled, see further in [3].

The random matrix transition density can expressed as

[TABLE]

where the last equality follows from a Markov assumption, see [2]. The random matrix measurement likelihood can be expressed on a general form as

[TABLE]

Note that the modelling of the extended object measurement set cardinality is outside the scope of this work, see [1, Sec. 2.C] for an overview of different models for the number of measurements.

IV Conditional model smoothing

In the conditional model, we have conditional Gaussian inverse Wishart densities, cf. (4b), and under assumed density filtering we seek a smoothed density of the same form, i.e.,

[TABLE]

IV-A Assumptions and modelling

The following assumptions are made for the conditional giw model, see [2, Sec. 2].

Assumption 1

The time evolution of the extent state is assumed independent of the kinematic state,

[TABLE]

$\square$**

Assumption 2

The extent changes slowly with time, $X_{k+1}\approx X_{k}$, such that for the kinematic state, conditioned on the extent state, the following holds,

[TABLE]

$\square$**

The validity of Assumptions 1 and 2 is discussed in [2].

In the conditional random matrix model, the transition density (6c) is Gaussian-Wishart, see [2, Sec. 3.A/B],

[TABLE]

In the remainder of the paper, we consider this generalised transition density. The measurement model is [2, Sec. 3.D]

[TABLE]

where the $1\times s$ matrix $H_{k}$ is the measurement model.

IV-B Prediction, update, and smoothing

The prediction and the update for the conditional model are reproduced in Table II and in Table III, respectively. The smoothing is given in the following theorem.

Theorem 1

Let the densities $p(\xi_{k}|\mathbf{Z}_{1:k})$, $p(\xi_{k+1}|\mathbf{Z}_{1:K})$ and $p(\xi_{k+1}|\mathbf{Z}_{1:k})$ be conditional Gaussian inverse Wishart (4b), and let the transition density be Gaussian Wishart (13c). The smoothed density $p(\xi_{k}|\mathbf{Z}_{1:K})$, see (3), is conditional Gaussian inverse Wishart, see with parameters $(m_{k|K},P_{k|K},v_{k|K},V_{k|K})$ given in Table IV. $\square$

The proof of Theorem 1 is given in Appendix -B.

V Factorized model

For the random matrix model in [3], we have factorised Gaussian inverse Wishart densities, cf. (5), and under assumed density filtering we seek a smoothed density of the same form, i.e.,

[TABLE]

V-A Assumptions, approximations

The following assumption is made for the factorized giw model, see [3].

Assumption 3

The time evolution of the kinematic state is independent of the extent state,

[TABLE]

$\square$**

The validity of this assumption is discussed in [3, 22]. The transition density is Gaussian Wishart [22],

[TABLE]

where the function $f_{k}(\cdot):\mathbb{R}^{n_{x}}\rightarrow\mathbb{R}^{n_{x}}$ is the motion model, the $n_{x}\times n_{x}$ matrix $Q_{k}$ is the process noise covariance, the degrees of freedom $n_{k}\geq d$ govern the uncertainty of the time evolution of the extent, and the function $M(\cdot):\mathbb{R}^{n_{x}}\rightarrow\mathbb{I}^{d\times d}$ describes how the extent changes over time due to the object motion. For example, $M(\cdot)$ can be a rotation matrix. In what follows, we write $M_{\mathbf{x}}=M(\mathbf{x})$ for brevity.

The measurement model is

[TABLE]

where the $d\times n_{x}$ matrix $\tilde{H}_{k}$ is the measurement model, $\rho>0$ is a scaling factor, and $R_{k}\in\mathbb{S}_{+}^{d}$ is the measurement noise covariance. The scaling factor $\rho$ and the noise covariance $R_{k}$ were added to better model scenarios where the sensor noise is large in relation to the size of the extended object, see discussion in [3, Sec. 3]. In this paper, to enable a straightforward comparison to the conditional model, which assumes that the sensor noise is small in comparison to the size of the extended object, we focus on the case $\rho=1$ and $R_{k}=\mathbf{0}_{d\times d}$.

V-B Prediction, update, and smoothing

The prediction and the update for the conditional model are reproduced in Table V and in Table VI, respectively. The smoothing is given in the following theorem.

Theorem 2

Let the densities $p(\xi_{k}|\mathbf{Z}_{1:k})$, $p(\xi_{k+1}|\mathbf{Z}_{1:K})$ and $p(\xi_{k+1}|\mathbf{Z}_{1:k})$ be factorised Gaussian inverse Wishart (5), and let the transition density be Gaussian Wishart (17). The smoothed density $p(\xi_{k}|\mathbf{Z}_{1:K})$, see (3), is factorised Gaussian inverse Wishart, see with parameters $(m_{k|K},P_{k|K},v_{k|K},V_{k|K})$ given in Table VII. $\square$

The proof of Theorem 2 is given in Appendix -C.

V-C Expected value approximation

Note that both the prediction and the smoothing require expected values, see $C_{1}$ and $C_{2}$ in Table V, and $C_{3}$ and $C_{4}$ in Table VII. For a Gaussian distributed vector $\mathbf{x}\sim\mathcal{N}(m,P)$, the expected value of $M_{\mathbf{x}}VM_{\mathbf{x}}^{\scriptscriptstyle\mathrm{T}}$ can be approximated using third order Taylor expansion,

[TABLE]

where $\mathbf{x}^{[i]}$ is the $i$th element of $\mathbf{x}$, $P^{[i,j]}$ is the $i,j$th element of $P$. The necessary differentiations are

[TABLE]

The expected values $C_{2}$, $C_{3}$ and $C_{4}$ can be approximated analogously.

VI Simulation study

In this section we present the results of a simulation study. In all simulations, the dimension of the extent is $d=2$.

VI-A Implemented smoothers

Three different smoothers were implemented.

VI-A1 Conditional giw model with constant velocity motion model (CCV)

The state vector contains 2D Cartesian position and velocity, $\mathbf{x}_{k}=[p_{k}^{x},\ p_{k}^{y},\ v_{k}^{x},\ v_{k}^{y}]^{\scriptscriptstyle\mathrm{T}}$, $n_{x}=4$ and $s=2$. The following models are used,

[TABLE]

$A=\mathbf{I}_{d}$, and $n_{k}=100$, where $T$ is the sampling time.

VI-A2 Factorized giw model with constant velocity motion model (FCV)

The state vector contains 2D Cartesian position and velocity, $\mathbf{x}_{k}=[p_{k}^{x},\ p_{k}^{y},\ v_{k}^{x},\ v_{k}^{y}]^{\scriptscriptstyle\mathrm{T}}$, and $n_{x}=4$. The following models are used,

[TABLE]

$n_{k}=100$, and $M_{\mathbf{x}}=\mathbf{I}_{d}$.

VI-A3 Factorized giw model with coordinated turn motion model (FCT)

The state vector contains 2D Cartesian position and velocity, as well as turn-rate, $\mathbf{x}_{k}=[p_{k}^{x},\ p_{k}^{y},\ v_{k}^{x},\ v_{k}^{y},\ \omega_{k}]^{\scriptscriptstyle\mathrm{T}}$, and $n_{x}=5$. The following models are used,

[TABLE]

and $n_{k}=\infty$. For the matrix transformation function $M_{\mathbf{x}}$ we have the following,

[TABLE]

VI-B Simulated scenarios

We focused on two types of scenarios: in the first the true tracks were generated by a constant velocity model; in the second the true tracks were generated by a coordinated turn model. This allows us to test the different smoothers both when their respective motion models match the true model, and when there is motion model mis-match.

The CV tracks were generated using the CV model (21a) and (21b) with $\sigma_{a}=1$; the extent’s major axis was simulated to be aligned with the velocity vector of the extended object. The CT tracks were generated using the CT model (23a) and (23b) with $\sigma_{a}=1$ and $\sigma_{\omega}=\pi/180$; the extent’s major axis was simulated to be aligned with the velocity vector of the extended object.

For both motion models, in each time step $k$, a detection process was simulated by first sampling a probability of detection $p_{\rm D}$, and, if the object is detected, sampling $N_{z}$ detections using a Gaussian likelihood. We simulated two combinations: $(p_{\rm D},\ N_{z})=(0.25,\ 10)$ and $(p_{\rm D},\ N_{z})=(0.75,\ 10)$.

VI-C Performance evaluation

For performance evaluation of extended object estimates with ellipsoidal extents, a comparison study has shown that among six compared performance measures, the Gaussian Wassterstein Distance (gwd) metric is the best choice [28]. The gwd is defined as [29]

[TABLE]

where the expected extended object state is

[TABLE]

and $\mathbf{p}_{k}$ is the position part of the extended object state vector $\mathbf{x}_{k}$.

VI-D Results

We show results for estimates $\hat{\xi}_{k|\ell}$ for $\ell\in\{k-1,\ k,\ K\}$, i.e., prediction, filtering and smoothing. Results for true tracks generated by a CV model are shown in Figure 1; results for true tracks generated by a CT model are shown in Figure 2.

We see that in all cases, as expected, the smoothing errors are smaller than the filtering errors, which are smaller than the prediction errors. This confirms that the derived smoothers work as they should. It is also in accordance with expectations that the performance is worse when the probability of detection is lower. Perhaps counter-intuitive is that for CV true motion, the FCT smoother performs best despite the modelling error in the motion model. We believe that this is due, at least in part, to the standard assumption that the orientation of the extent ellipse is aligned with the velocity vector. The motion noice on the velocity vector introduces rotations on the extent ellipse, and the motion model used in FCT captures these rotations better.

VII Conclusions and future work

This paper presented Bayesian smoothing for the random matrix model used in extended object tracking. Two variants of Gaussian inverse Wishart state densities exist in the literature, a conditional and a factorised; closed form Bayesian smoothing was derived for both of them. The derived smoothers were implemented and tested in a simulation scenario. In future work the smoothers will be used with real data, e.g., data from camera, lidar or radar.

-A Preliminary results

In this appendix we present some preliminary results that are used in the proofs of Theorems 1 and 2. The first two Lemmas regard products, and ratios, of inverse Wishart pdfs, respectively.

Lemma 1

The product of two inverse Wishart pdfs is proportional to an inverse Wishart pdf,

[TABLE]

$\square$**

Proof: This follows from the definition of the inverse Wishart pdf.

Lemma 2

The fraction of two inverse Wishart pdfs is proportional to an inverse Wishart pdf,

[TABLE]

$\square$**

Proof: This follows from the definition of the inverse Wishart pdf.

The following two Lemmas are related to the use Wishart transition densities and inverse Wishart state densities for the extent matrix.

Lemma 3

For Wishart and inverse Wishart pdfs, the following holds,

[TABLE]

$\square$**

Proof: This follows from the definitions of the Wishart pdf and the inverse Wishart pdf.

Lemma 4

[TABLE]

Proof: see [23, Prob. 5.33].

Lemma 5

[TABLE]

Proof: see [22, Thm. 3].

For approximations of densities, the Kullback-Leibler divergence (kl-div) is often minimised to find the optimal approximation. The following Lemma is about approximation by a factorised density.

Lemma 6

For two random variables $\mathbf{x}$ and $X$, with joint density $p(\mathbf{x},X)$, the factorised density $q^{\star}(\mathbf{x})q^{\star}(X)$ that minimises the kl-div to $p(\mathbf{x},X)$,

[TABLE]

is given by the marginals,

[TABLE]

$\square$**

Proof: This is a previously know result that follows from the definition of the Kullback-Leibler divergence.

Approximation of matrix valued densities as Wishart, or inverse Wishart, densities, by minimisation of the kl-div, is presented in [22, Thm. 1] and [30, Thms. 3 & 4]. The kl-div minimisation leads to matching of the expected logarithm of the determinant of the extent matrix, as well as either the expected extent matrix, or the expected inverse extent matrix. A simpler closed form approximation is obtained if instead the expected random matrix and the expected inverse are match, i.e., not matching the expected log-determinant; this is shown in the following four Lemmas.

Lemma 7

By matching the expected values $E[X]$ and $E[X^{-1}]$, the inverse Wishart density $\mathcal{IW}_{d}\left(X\,;\,v,V\right)$ can be approximated by a Wishart density $\mathcal{W}_{d}\left(X\,;\,w,W\right)$ with parameters

[TABLE]

$\square$**

Proof: this follows from the definitions of the expected values, see [23].

Lemma 8

By matching the expected values $E[X]$ and $E[X^{-1}]$, the Wishart density $\mathcal{IW}_{d}\left(X\,;\,w,W\right)$ can be approximated by an inverse Wishart density $\mathcal{W}_{d}\left(X\,;\,v,V\right)$ with parameters

[TABLE]

$\square$**

Proof: this follows from the definitions of the expected values, see [23].

Lemma 9

By matching the expected values $E[X]$ and $E[X^{-1}]$, the generalized Beta type 2 density $\mathcal{GB}_{d}^{II}\left(X;\frac{a}{2},\frac{b}{2},A,\mathbf{0}_{d\times d}\right)$ can be approximated by a Wishart density $\mathcal{W}_{d}\left(X\,;\,w,W\right)$ with parameters

[TABLE]

$\square$**

Proof: this follows from the definitions of the expected values, see [23].

Lemma 10

By matching the expected values $E[X]$ and $E[X^{-1}]$, the generalized Beta type 2 density $\mathcal{GB}_{d}^{II}\left(X;\frac{a}{2},\frac{b}{2},A,\mathbf{0}_{d\times d}\right)$ can be approximated by an inverse Wishart density $\mathcal{W}_{d}\left(X\,;\,v,V\right)$ with parameters

[TABLE]

$\square$**

Proof: this follows from the definitions of the expected values, see [23].

-B Conditional model smoothing

For conditional densities (4a) and the transition density (13a), under Assumptions 1 and 2, the Bayesian smoothing (3) leads to a conditional smoothed density

[TABLE]

The proof of (44) is given in (45).

We get the following smoothed conditional giw density

[TABLE]

with the parameters given in Table IV. The proof of (46a) is simple; the details follow the proof of the RTS-smoother, see, e.g., [31, Thm. 8.2]. The proof of (46b) is given in (47).

-C Factorized model smoothing

In the factorized case, there is no known analytical solution that gives a smoothed density of the desired factorized form; therefore approximations are necessary. Factorized density approximations are common in so called variational inference, and the factors are typically found by minimising the Kullback-Leibler divergence, see, e.g., [32, Ch. 10]. To find a factorised smoothed density, we apply Lemma 6 to the smoothed joint density $p(\xi_{k}|\mathbf{Z}_{1:K})$, given in (3), and obtain the following two smoothing equations,

[TABLE]

and

[TABLE]

The proof of the marginalisation (48) is given in (50), and the proof of the marginalisation (49) is given in (51).

For the kinematic vector, we have that for Gaussian densities $p(\mathbf{x}_{k+1}|\mathbf{Z}_{1:K})$ and $p(\mathbf{x}_{k+1}|\mathbf{Z}_{1:k})$, see (5), and a Gaussian transition density $p(\mathbf{x}_{k+1}|\mathbf{x}_{k})$, see (17), the smoothed kinematic state density is Gaussian with parameters given by the standard RTS-smoothing backwards step, given in, e.g., [31, Thm. 8.2]. We get the result in Table VII.

For the extent matrix, the smoothing (49) does not have an analytical solution, and approximations are necessary. The derivation of the result in Table VII is given in (52).