Sage Journals: Discover world-class research

Abstract

In this paper, we propose a novel simultaneous registration and fusion approach for tracking. This method is based on a recursive Variational Bayesian (RVB) algorithm, which is the online variant of the Variational Bayesian (VB) approach. Under the Bayesian framework, the states and parameters are recursively estimated. It is shown by simulation that the proposed RVB method has better estimation performance than the conventional approach.

Keywords

Sensor Registration Sensor Fusion Recursive Variational Bayesian

1. Introduction

Sensor fusion is an essential component for a sensor network to improve the overall system performance [1]. It is a method of integrating signals from multiple sources into a single signal or piece of information. Considering the fact that multiple sources have their own coordinate systems, it is necessary to transform the local measurements into a common coordinate and estimate the sensor bias to perform alignment. The fusion centre is then applied to estimate the target states. Therefore, the fusion process relies on an accurate registration of sensors. The lack of proper registration will cause poor performance in the fusion process.

In order to obtain more information, intelligent vehicles should register and fuse numerous data from multi-sensors, such as cameras, lidars and radars. For instance, in a multi-radar multi-target system, each radar range and azimuth bias can cause significant errors in the target location. Therefore, many sensor registration algorithms have been proposed in the literature, such as the least squares (LS) method [2] and the maximum likelihood estimator method [3–4].

Some methods have been proposed to address the problem of registration and fusion simultaneously, which can be performed at the measurement level or track level. In [5], a joint sensor registration and fusion approach at the measurement level was proposed using an expectation-maximization (EM) algorithm. In [6], by using a pseudo-measurement approach, simultaneous registration and fusion was performed at the track level.

However, such studies focus on the batch algorithm, which estimates the parameter based on all data. This method is impractical for real-time processing. Therefore, many online algorithms have been proposed to address this problem. The traditional method is that the sensor bias is included in the augmented state equations, and then the Kalman filter, extended Kalman filter (EKF) or unscented Kalman filter (UKF) is used to estimate the augmented states’ systems. In [7], Okello, et al. formulated the joint registration and fusion at the track level as a Bayesian estimation problem, and proposed the extended Kalman filter (EKF) by augmenting the state vector with sensor bias. In [8], the augmented Kalman filter was proposed to perform the sensor registration.

In this paper, we propose a novel method to recursively perform simultaneous registration and fusion at the measurement level. The RVB [9] is the recursive variant of Variational Bayesian (VB) inference, which is an approximation to the exact Bayesian inference, where the true posterior is approximated using a simpler distribution. VB has low computational costs in comparison with the sampling methods, such as Markov Chain Monte Carlo (MCMC). It has widespread use in machine learning, signal processing and many other fields, such as the state space model [10, 11], time series [12, 13], Mixture Models [14], filter [15 –17], image [18 –20], communication [21, 22], speech recognition [23] and graphical models [24, 25].

The problem of joint sensor registration and fusion at the measurement level has been seen as the estimation of parameters on the Bayesian framework. Combining them with the expectation-maximization (EM) algorithm, states and parameters are formulated via Recursive Variational Bayesian Expectation Maximization (RVBEM). Compared with the augmented Kalman filter [8], the proposed method has a better performance.

The paper is organized as follows. The problem of simultaneous registration and fusion is formulated in Section 2. In Section 3, the RVBEM for simultaneous registration and fusion are developed in detail. Computer simulations are shown in Section 4. Finally, some concluding remarks are given in Section 5.

2. Problem Formulation

In tracking, it is assumed that the target state evolves according to a linear, discrete time model of the form:

x_{t} = F x_{t - 1} + w_{t}

(1)

where $x_{t}$ is the state vector at time instant t, F is the known transition matrix and $w_{t}$ is a zero mean, white Gaussian noise with precision Q (inverse covariance). The measurement equation can be described as:

z_{t} = h (x_{t}) + η + v_{t}

(2)

where $z_{t}$ is the measurement at time instant t, h(.) is the known nonlinear measurement function, $η$ denotes the sensor bias and $v_{t}$ is also a zero mean, white Gaussian noise with precision R (inverse covariance). The measurement noise $v_{t}$ is independent of the process noise.

The state is driven by an initial state x₀, which is assumed to be a Gaussian process with mean $θ_{0}$ and precision $Δ_{0}$ . In the real world, the measurement noise $v_{t}$ is non-stationary with the time varied precision matrix R. Set $ρ = [η, R]$ . Our objective is to obtain a real-time accurate estimate of $ρ$ and the system state $x_{t}$ based on the sequential observations {z₁, z₂……., z_T}, where T is the total number of measurements.

Then, we choose the standard conjugate priors for the model. Firstly, the prior of the initial state x0 and the sensor bias is defined as:

\begin{array}{l} p (x_{0}) = N (θ_{0}, Δ_{0}) \\ p (η) = \prod_{t = 1}^{T} N (η_{t}; u_{t}, Σ_{t}) \end{array}

(3)

The Gaussian density is defined as:

N (x; u, Σ) = {(2 π)}^{- \frac{d}{2}} {| Σ |}^{\frac{1}{2}} e^{- \frac{1}{2} t r (Σ (x - u) {(x - u)}^{T})}

(4)

where u and Σ is the mean and precision of the normal distributed, respectively, and d is the number of the dimension of x.

The prior over the measurement noise precision matrix R, which is assumed to be diagonal, has Gamma prior:

p (R) = \prod_{t = 1}^{T} \prod_{j = 1}^{m} G a m m a (R_{t}_{_{(j, j)}}; a_{t}, b_{t})

(5)

where m is the number of dimension of R, $R_{(j, j)}$ denotes the element of the jth row and the jth column of a matrix R.

The Gamma density is defined as:

G a m m a (x; a, b) = \frac{b^{a}}{Γ (a)} x^{a - 1} e^{- b x}

(6)

where a and b is the shape and inverse scale of the gamma distributed, respectively.

Therefore, the overall joint density as:

\begin{array}{l} P (X, Z, ρ) = \prod_{t = 1}^{T} p (z_{t} | x_{t}, R_{t}) \times \\ \times \prod_{t = 1}^{T} p (x_{t} | x_{t - 1}, Q, η_{t}) p (R) p (η) p (x_{0}) \end{array}

(7)

where X={x₁, x₂…… x_T} is the hidden state and Z={z₁, z₂……., z_T} denotes the data.

3. RVBEM for Simultaneous Registration and Fusion Recursively

Based on the Bayes’ theorem:

p (X, ρ | Z) = \frac{p (Z, X, ρ)}{p (Z)} = \frac{p (Z | X, ρ) p (X, ρ)}{p (Z)}

(8)

where $p (X, ρ | Z)$ is the posterior distribution of the central interest in Bayesian learning. The equation fully describes our knowledge regarding the parameter of the model. The likelihood of the data is $p (Z | X, ρ)$ , $p (X, ρ)$ denotes the prior distribution and $p (Z)$ represents the evidence or marginal likelihood for the measurements. In most cases, the posterior is difficult to estimate, as it involves an intractable integral. The VB or ensemble learning can circumvent this problem. A general tutorial on the VB can be found in [26, 27]. In what follows, we describe the key feature.

The log of marginal likelihood can be written as:

\begin{array}{l} \log p (Z) = \int q (X, ρ | Z) \log p (Z) d X d ρ \\ = \int q (X, ρ | Z) \log \frac{p (Z, X, ρ)}{p (X, ρ | Z)} d X d ρ \\ = \int q (X, ρ | Z) \log \frac{p (Z, X, ρ)}{p (X, ρ | Z)} \frac{q (X, ρ | Z)}{q (X, ρ | Z)} d X d ρ \\ = \int q (X, ρ | Z) \log \frac{p (Z, X, ρ)}{q (X, ρ | Z)} d X d ρ + \\ + \int q (X, ρ | Z) \log \frac{q (X, ρ | Z)}{p (X, ρ | Z)} d X d ρ \end{array}

(9)

where $q (X, ρ | Z)$ is the approximation of the posterior density. The equation can be rewritten as:

\log p (Z) = F + K L

(10)

where F is the free energy:

F = \int q (X, ρ | Z) \log \frac{p (Z, X, ρ)}{q (X, ρ | Z)} d X d ρ

(11)

and KL is the Kullback-Liebler divergence between the approximate posterior $q (X, ρ | Z)$

and the true posterior $p (X, ρ | Z)$ , given as:

K L = \int q (X, ρ | Z) \log \frac{q (X, ρ | Z)}{p (X, ρ | Z)} d X d ρ

(12)

The goal of VB is to minimize the KL-divergence. However, the KL-divergence cannot be directly minimized, because the true posterior $p (X, ρ | Z)$ is unknown. Since the KL-divergence is always positive and the free energy and KL-divergence must equal the log marginal likelihood, we can achieve the same aim by maximizing the free energy instead. To make the integrals in free energy tractable, the mean field approximation is proposed. Thus:

q (X, ρ | Z) = q (X | Z) \prod_{i} q (ρ_{i} | Z)

(13)

where $ρ_{i}$ is the group of parameters. It assumes that the parameters in the separate group are independent. The computation of $q (ρ_{i})$ by the maximization of $q (ρ_{i})$ over free energy is given as:

\log q (ρ_{i}) \propto \int q (ρ_{i}) \log (p (X, Z, ρ)) d ρ_{i}

(14)

q (ρ_{i} | Z) = \frac{\exp [I (ρ_{i})]}{\int \exp [I (ρ_{i})] d ρ_{i}}

(15)

where

I (ρ_{i}) = \int q (ρ_{i} | Z) \log (p (X, Z, ρ)) d ρ_{i}

(16)

and $ρ_{i}$ refer to all the parameters not in the ith group. The equation (16) can be rewritten in terms of an expectation:

I (ρ_{i}) = {〈 \log (p (X, Z, ρ)) 〉}_{q_{ρ_{i}}}

(17)

where $< >_{x}$ is the expectation of the expression taken with respect to x.

This leads to a set of coupled update rules. Iterated application of these leads to the desired maximization. The resulting VB algorithm follows an EM update procedure. In the Variational Bayesian expectation (VBE)-step, the variational posteriors over hidden variables are calculated. In the Variational Bayesian Maximization (VBM)-step, the variational posteriors over parameters are updated. Each step is guaranteed to increase the lower bound of the marginal likelihood. We choose every size of the batch from a single data point. The RVBE and RVBM are iterated until convergence, which occurs rapidly. After convergence, the prior of parameters for the next batch is set as the current posterior.

The optimal components of the approximation posterior are defined as:

\begin{array}{l} q (η) = \prod_{t = 1}^{T} N ({\hat{η}}_{t}; {\hat{u}}_{t}, {\hat{Σ}}_{t}) \\ q (R) = \prod_{t = 1}^{T} \prod_{j = 1}^{m} G a m m a ({\hat{R}}_{t}_{_{(j, j)}}; {\hat{a}}_{t_{(j)}}, {\hat{b}}_{t_{(j)}}) \\ q (x_{t}) = N (θ_{t}, Δ_{t}) \end{array}

(18)

Note that the form of the posterior distribution is the same as in the prior distribution. The prior distributions are chosen to conjugate with the likelihood [26]. In what follows, the states and parameters of distributions are updated.

E-step

In this step, we update the distribution of hidden variables using the variational forward algorithm [11].

Forward pass: the predictive density:

\begin{array}{l} α_{t}^{*} (x_{t}) = p (x_{t} | z_{1 : t - 1}) \\ \propto \int d x_{t - 1} p (x_{t - 1} | z_{1 : t - 1}) p (x_{t} | x_{t - 1}) \\ = \int d x_{t - 1} α_{t - 1} (x_{t - 1}) p (x_{t} | x_{t - 1}) \\ = \int d x_{t - 1} N (x_{t - 1}; θ_{t - 1}, Δ_{t - 1}) N (x_{t}; F x_{t - 1}, Q) \\ = N (x_{t}; m_{t}, Ω_{t}) \end{array}

(19)

where $α_{t - 1} (x_{t - 1})$ is the previous hidden state density. Set $f = α_{t - 1} (x_{t - 1}) p (x_{t} | x_{t - 1})$ , then

\begin{array}{l} f = α_{t - 1} (x_{t - 1}) p (x_{t} | x_{t - 1}) \\ \propto - \frac{1}{2} {(x_{t - 1} - θ_{t - 1})}^{T} Δ_{t - 1} (x_{t - 1} - θ_{t - 1}) - \frac{1}{2} {(x_{t} - F x_{t - 1})}^{T} Q (x_{t} - F x_{t - 1}) \\ = - \frac{1}{2} (x_{t - 1}^{T} Δ_{t - 1} x_{t - 1} - x_{t - 1}^{T} Δ_{t - 1} θ_{t - 1} - θ_{t - 1}^{T} Δ_{t - 1} x_{t - 1} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1} + \\ x_{t}^{T} Q x_{t} - x_{t}^{T} Q F x_{t - 1} - x_{t - 1}^{T} F^{T} Q x_{t} + x_{t - 1}^{T} F^{T} Q F x_{t - 1}) \\ = - \frac{1}{2} [x_{t - 1}^{T} (Δ_{t - 1} + F^{T} Q F) x_{t - 1} - x_{t - 1}^{T} (Δ_{t - 1} θ_{t - 1} + F^{T} Q x_{t}) - \\ (θ_{t - 1}^{T} Δ_{t - 1} + x_{t}^{T} Q F) x_{t - 1} + x_{t}^{T} Q x_{t} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \end{array}

(20)

Set:

\begin{array}{l} Β_{t} = Δ_{t - 1} + F^{T} Q F \\ ω_{t} = Β_{t}^{- 1} (Δ_{t - 1} θ_{t - 1} + F^{T} Q x_{t}) \end{array}

(21)

Then the predictive density $α_{t}^{*} (x_{t})$ becomes:

\begin{array}{l} α_{t}^{*} (x_{t}) \propto \int d x_{t - 1} - \frac{1}{2} [{(x_{t - 1} - ω_{t})}^{T} Β_{t} (x_{t - 1} - ω_{t}) - ω_{t}^{T} Β_{t} ω_{t} + x_{t}^{T} Q x_{t} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \\ \propto - \frac{1}{2} [- ω_{t}^{T} Β_{t} ω_{t} + x_{t}^{T} Q x_{t} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \\ = - \frac{1}{2} [- (θ_{t - 1}^{T} Δ_{t - 1} + x_{t}^{T} Q F) Β_{t}^{- 1} Β_{t} Β_{t}^{- 1} (Δ_{t - 1} θ_{t - 1} + F^{T} Q x_{t}) + \\ + x_{t}^{T} Q x_{t} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \\ = - \frac{1}{2} [- θ_{t - 1}^{T} Δ_{t - 1} Β_{t}^{- 1} Δ_{t - 1} θ_{t - 1} - θ_{t - 1}^{T} Δ_{t - 1} Β_{t}^{- 1} F^{T} Q x_{t} - x_{t}^{T} Q F Β_{t}^{- 1} Δ_{t - 1} θ_{t - 1} \\ - x_{t}^{T} Q F Β_{t}^{- 1} F^{T} Q x_{t} + x_{t}^{T} Q x_{t} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \\ = - \frac{1}{2} [x_{t}^{T} (Q + Q F Β_{t}^{- 1} F^{T} Q) x_{t} - θ_{t - 1}^{T} Δ_{t - 1} Β_{t}^{- 1} F^{T} Q x_{t} - x_{t}^{T} Q F Β_{t}^{- 1} Δ_{t - 1} θ_{t - 1} \\ - θ_{t - 1}^{T} Δ_{t - 1} Β_{t}^{- 1} Δ_{t - 1} θ_{t - 1} + θ_{t - 1}^{T} Δ_{t - 1} θ_{t - 1}] \end{array}

(22)

Therefore, the mean and covariance of the predictive density are defined as:

\begin{array}{l} Ω_{t} = Q + Q F Β_{t}^{- 1} F^{T} Q \\ m_{t} = Ω_{t}^{- 1} (Q F Β_{t}^{- 1} Δ_{t - 1} θ_{t - 1}) \end{array}

(23)

Forward pass: the updated density:

α_{t} (x_{t}) = p (x_{t} | z_{1 : t}) \propto α_{t}^{*} (x_{t}) 〈 p (z_{t} | x_{t}) 〉

(24)

Taylor-expand h(x) at the predictive estimate of state:

h (x_{t}) \approx h (m_{t}) + H_{m_{t}} (x_{t} - m_{t})

(25)

where H is the Jacobian

H_{m_{t}} = \frac{d h (x)}{d x} |_{x = m_{t}}

(26)

\begin{array}{l} α_{t} (x_{t}) \propto α_{t}^{*} (x_{t}) 〈 p (z_{t} | x_{t}) 〉 \\ \propto - \frac{1}{2} [{(x_{t} - m_{t})}^{T} Ω_{t} (x_{t} - m_{t}) + 〈 {(z_{t} - η_{t} - h (m_{t}) - H_{m_{t}} (x_{t} - m_{t}))}^{T} \\ R_{t} (z_{t} - η_{t} - h (m_{t}) - H_{m_{t}} (x_{t} - m_{t})) 〉] \\ \propto - \frac{1}{2} [x_{t}^{T} Ω_{t} x_{t} - x_{t}^{T} Ω_{t} m_{t} - m_{t}^{T} Ω_{t} x_{t} + m_{t}^{T} Ω_{t} m_{t} - x_{t}^{T} H_{m_{t}}^{T} 〈 R_{t} 〉 \\ 〈 (z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t}) 〉 - \\ 〈 {(z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t})}^{T} 〉 〈 R_{t} 〉 H_{m_{t}} x_{t} + x_{t}^{T} H_{m_{t}}^{T} 〈 R_{t} 〉 H_{m_{t}} x_{t} \\ + 〈 {(z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t})}^{T} R_{t} (z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t}) 〉 \\ \propto - \frac{1}{2} [x_{t}^{T} (Ω_{t} + H_{m_{t}}^{T} 〈 R_{t} 〉 H_{m_{t}}) x_{t} - x_{t}^{T} \\ (Ω_{t} m_{t} + H_{m_{t}}^{T} 〈 R_{t} 〉 〈 (z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t}) 〉) \\ - (m_{t}^{T} Ω_{t} + 〈 {(z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t})}^{T} 〉 〈 R_{t} 〉 H_{m_{t}}) x_{t} \\ + m_{t}^{T} Ω_{t} m_{t} + 〈 {(z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t})}^{T} \\ 〈 R_{t} 〉 (z_{t} - η_{t} - h (m_{t}) + H_{m_{t}} m_{t}) 〉] \end{array}

(27)

Then, we get the mean and covariance of the updated density:

\begin{array}{l} Δ_{t} = Ω_{t} + H_{m_{t}}^{T} {\hat{R}}_{t} H_{m_{t}} \\ θ_{t} = Δ_{t}^{- 1} (Ω_{t} m_{t} + H_{m_{t}}^{T} {\hat{R}}_{t} (z_{t} - {\hat{η}}_{t} - h (m_{t}) + H_{m_{t}} m_{t})) \end{array}

(28)

M-step

The detailed derivations of the updated parameters are given in this step.

Update for η

\begin{array}{l} \log q (η_{t}) \propto {〈 \log p (z_{t} {| x_{t}, R}_{t}, η_{t}) 〉}_{q (x) q (R_{t})} + \log p (η_{t}) \\ \propto 〈 - \frac{m}{2} \log 2 π + \frac{1}{2} \log R_{t} - \frac{1}{2} \\ {〈 t r (R_{t} (z_{t} - η_{t} - h (x_{t})) {(z_{t} - η - h (x_{t}))}^{T}) 〉}_{q (x) q (R_{t})} \\ - \frac{1}{2} \log 2 π + \frac{1}{2} \log Σ_{t} - \frac{1}{2} t r [Σ_{t} (η_{t} - u_{t}) {(η_{t} - u_{t})}^{T}] \\ \propto - \frac{m T}{2} \log 2 π + \frac{T}{2} {〈 \log R_{t} 〉}_{q (R_{t})} - \\ - \frac{1}{2} \sum_{t = 1}^{T} 〈 t r ({〈 R_{t} 〉}_{q (R_{t})} (z_{t} - η_{t} - h (θ_{t}) - H_{θ_{t}} (x_{k} - θ_{t})) \\ {\times {(z_{t} - η_{t} - h (θ_{t}) - H_{θ_{t}} (x_{k} - θ_{t}))}^{T}) 〉}_{q (x)} - \frac{1}{2} \log 2 π + \frac{1}{2} \log Σ_{t} \\ - \frac{1}{2} t r [Σ_{t} (η_{t} - u_{t}) {(η_{t} - u_{t})}^{T}] \\ \propto - \frac{m T + 1}{2} \log 2 π - \frac{T}{2} {〈 \log R_{t} 〉}_{q (R_{t})} + \frac{1}{2} \log Σ_{t} \\ - \frac{1}{2} \sum_{t = 1}^{T} t r ({〈 R_{t} 〉}_{q (R_{t})} [(z_{t} - h (θ_{t})) {(z_{t} - h (θ_{t}))}^{T} - η_{t} {(z_{t} - h (θ_{t}))}^{T} \\ - (z_{t} - h (θ_{t})) η_{t}^{T} + η_{t} η_{t}^{T} + H_{θ_{t}} Δ_{t} H_{θ_{t}}^{T}]) - \\ - \frac{1}{2} t r [Σ_{t} (η_{t} η_{t}^{T} - η_{t} u_{t}^{T} - u_{t} η_{t}^{T} + u_{t} u_{t}^{T})] \end{array}

(29)

The form of log posterior density can be rewritten as:

\log q (η) = - \frac{1}{2} \log 2 π + \frac{1}{2} \log {\hat{Σ}}_{t} - \frac{1}{2} t r [{\hat{Σ}}_{t} (η_{t} - {\hat{u}}_{t}) {(η_{t} - {\hat{u}}_{t})}^{T}]

(30)

Comparing the coefficients with equations (29) and (30), we get:

\begin{array}{l} {\hat{Σ}}_{t} = Σ_{t} + {\hat{R}}_{t} \\ {\hat{u}}_{t} = {\hat{Σ}}_{t}^{- 1} {\hat{R}}_{t} (z_{t} - h (θ_{t}) + Σ_{t} u_{t}) \end{array}

(31)

So: ${\hat{η}}_{t} = {\hat{u}}_{t}$

The prior of the next batch is set as:

\begin{array}{l} Σ_{t + 1} = {\hat{Σ}}_{t} \\ u_{t + 1} = {\hat{u}}_{t} \end{array}

(32)

Update for R:

\begin{array}{l} \log q (R_{t}) \propto {〈 \log p (z_{t} {| x_{t}, R}_{t}, η_{t}) 〉}_{q (x) q (η_{t})} + \log p (R_{t}) \\ \propto 〈 \frac{1}{2} \log R_{t} - \frac{1}{2} {〈 t r (R_{t} (z_{t} - η_{t} - h (x_{t})) {(z_{t} - η_{t} - h (x_{t}))}^{T}) 〉}_{q (x) q (η_{t})} \\ + \sum_{j = 1}^{m} (a_{t}_{_{(j)}} - 1) \log R_{t}_{_{(j, j)}} - b_{t}_{_{(j)}} R_{t}_{_{(j, j)}} \\ \propto \frac{T}{2} \log R_{t} - \frac{1}{2} 〈 t r (R_{t} (z_{t} - η_{t} - h (θ_{t}) - H_{θ_{t}} (x_{k} - θ_{t})) \\ {\times {(z_{t} - η_{t} - h (θ_{t}) - H_{θ_{t}} (x_{k} - θ_{t}))}^{T}) 〉}_{q (x) q (η_{t})} + \\ + \sum_{j = 1}^{m} (a_{t}_{_{(j)}} - 1) \log R_{t}_{_{(j, j)}} - b_{t}_{_{(j)}} R_{t}_{_{(j, j)}} \\ \propto \frac{1}{2} \log R_{t} - \frac{1}{2} t r (R_{t} [(z_{t} - h (θ_{t})) {(z_{t} - h (θ_{t}))}^{T} - {\hat{u}}_{t} {(z_{t} - h (θ_{t}))}^{T} \\ - (z_{t} - h (θ_{t})) {\hat{u}}_{t}^{T} + {\hat{u}}_{t} {\hat{u}}_{t}^{T} + {\hat{Σ}}_{t}^{- 1} + H_{θ_{t}} Δ_{t} H_{θ_{t}}^{T}]) + \\ + \sum_{j = 1}^{m} (a_{t}_{_{(j)}} - 1) \log R_{t}_{_{(j, j)}} - b_{t}_{_{(j)}} R_{t}_{_{(j, j)}} \end{array}

(33)

The form of the log posterior density can be rewritten as:

\log q (R_{t}) = \sum_{j = 1}^{m} ({\hat{a}}_{t}_{_{(j)}} - 1) \log R_{t}_{_{(j, j)}} - {\hat{b}}_{t}_{_{(j)}} R_{t}_{_{(j, j)}}

(34)

Comparing the coefficients with equations (33) and (34), we get:

\begin{array}{l} {\hat{a}}_{t}_{_{(j)}} = a_{t}_{_{(j)}} + \frac{1}{2} \\ {\hat{b}}_{t}_{_{(j)}} = b_{t}_{_{(j)}} + ((z_{t} - h (θ_{t}) - {\hat{u}}_{t}) {(z_{t} - h (θ_{t}) - {\hat{u}}_{t})}^{T} + \\ + {\hat{Σ}}_{t}^{- 1} + H_{θ_{t}} Δ_{t} H_{θ_{t}}^{T})_{(j, j)} \end{array}

(35)

Therefore, ${\hat{R}}_{t} = d i a g [\frac{a_{t_{(1)}}}{b_{t}_{_{(1)}}}, \frac{a_{t_{(2)}}}{b_{t}_{_{(2)}}} \dots \frac{a_{t_{(m)}}}{b_{t}_{_{(m)}}}]$

The prior of the next batch is set as:

\begin{array}{l} a_{t}_{_{(j)}} = {\hat{a}}_{t}_{_{(j)}} \\ b_{t}_{_{(j)}} = {\hat{b}}_{t}_{_{(j)}} \end{array}

(36)

In [26], Beal et al. pointed out that the underfitting may occur in this online algorithm, which is due to excessive self-pruning of the parameters by the VB algorithm. An annealing variant of the RVBEM algorithm by making use of an inverse temperature parameter is proposed. We adopt the following update for the RVBM step instead of eqs. (31)–(35).

\begin{array}{l} {\hat{Σ}}_{t} = Σ_{t} + β {\hat{R}}_{t} \\ {\hat{u}}_{t} = {\hat{Σ}}_{t}^{- 1} {\hat{R}}_{t} (β (z_{t} - h (θ_{t})) + Σ_{t} u_{t}) \end{array}

(37)

\begin{array}{l} {\hat{a}}_{t}_{_{(j)}} = a_{t}_{_{(j)}} + \frac{1}{2} β \\ {\hat{b}}_{t}_{_{(j)}} = b_{t}_{_{(j)}} + β ((z_{t} - h (θ_{t}) - {\hat{u}}_{t}) {(z_{t} - h (θ_{t}) - {\hat{u}}_{t})}^{T} + \\ + {\hat{Σ}}_{t}^{- 1} + H_{θ_{t}} Δ_{t} H_{θ_{t}}^{T})_{(j, j)} \end{array}

(38)

The pseudo-code of RVBEM is given as follows:

Initialize

For t=1,2……T (observations)

2.1 RVBEM

RVB-E step: Update the parameters of the approximate distributions’ hidden states using the variational forward algorithm.

RVB-M step: Update the parameters of the approximate posteriors over the parameter as in (37) and (38).

2.2 After convergence, the prior of the next batch is set to the current posterior as in (32) and (36).

4. Simulations

In this section, we used the proposed algorithm to estimate registration and fusion simultaneously online via computer simulations. Two radars are used to track one target with a constant velocity motion. The positions of Sensor 1 and Sensor 2 are at (0 km, 0 km) and (200 km, 0km), respectively. The target speed is 0.5 km per second at the x-coordinate and 0 km per second at the y-coordinate. The initial position of the target is at (50 km, 30 km). The sensors measure the target every two seconds. The registration error of Sensor 1 is (0.2 km, −0.15 rad), and the registration error of Sensor 2 is (−0.2 km, 0.15 rad). The process noise covariance is set as diag(0.0001km², 0.0001rad², 0.0001km2, 0.0001rad²). The number of time instants is 150.

Figure 1.

Estimation of the range bias at Sensor 1

Figure 2.

Estimation of the azimuth bias at Sensor 1

In practice, we set the priors of gamma distribution as $a_{1} = 10^{- 3}, b_{1} = 10^{- 3}$ , and the priors of the sensor bias as $μ_{1} = 0$ , $Σ = I_{4}$ . Here $, I_{4}$ represents the identity matrix of 4*4. The inverse temperature parameter $β$ is set as 1-exp(−4) [17]. The RVBEM method stops in every batch when the difference in the parameters’ estimates in consecutive iterations is 0.001. To illustrate the effectiveness of the proposed RVB methods, we choose the Augmented Kalman Filter [8] for comparison.

Case 1. We choose the measurement noise covariance to be set equal to diag(0.0005km², 0.0005rad²). To illustrate the registration performance of the proposed method, Figures 1 to 4 show the convergence of the estimated registration parameters.

Figure 3.

Estimation of the range bias at Sensor 2

Figure 4.

Estimation of the azimuth bias at Sensor 2

Figure 5.

MSE in the position estimation of x

We evaluate the performance using mean squared error (MSE) over 1000 runs for the estimated target state. Figures 5 to 6 show the MSE of the estimate of the target using the two methods. It is observed that the estimation performance of the proposed method is better than that of the Augmented Kalman Filter.

Figure 6.

MSE in the position estimation of y

The computation is performed on a desktop computer with a Core 2 DUO 2.53 GHz and 2GB RAM. For a single run, the RVB requires 0.188s and the Augmented Kalman filter needs 0.093s. By augmenting the state vector in the augmented Kalman filter, the memory requirement of the augmented Kalman filter is more expensive than that of the RVB.

Figure 7.

Estimation of the variance at Sensor 1 inthe range

Case 2. In order to illustrate the performance of estimating the non-stationary time-varying parameters, we set the measurement noise first, which has variance diag(0.0005km², 0.0005rad²); at time t=1000, the variance is then quickly increased to diag(0.002km², 0.002rad²) and then t=2000 is again quickly decreased back to diag(0.0005km², 0.0005rad²). In addition, the number of measurements taken is 4000. The other parameters are the same as in Case 1. The results of the estimates of non-stationary measurement noise variance are shown in Figures 7 to 10. It is observed that the estimates of the proposed method follow the true variance fairly well.

Figure 8.

Estimation of the variance at Sensor 1 in azimuth

Figure 9.

Estimation of the variance at Sensor 2 in range

Figure 10.

Estimation of the variance at Sensor 2 in azimuth azimuth

5. Conclusion

In this paper, the problem of online joint sensor registration and fusion were seen in terms of estimating parameters of the Bayesian framework. We proposed using a RVBEM algorithm to perform it simultaneously. The posterior distributions of the states and parameters were presented in detail. Our proposed method shows better performance when compared to the augmented Kalman filter, which was validated using the computer simulation. Future work will plan to extend the algorithm to cope with the problem of joint sensor registration and fusion when the target number is unknown and variable.

Footnotes

6. Acknowledgements

This work is jointly supported by the Zhejiang Open Foundation of the Most Important Subjects, the Research Funds of Chongqing Science and Technology Commission (Grant No.cstc2013jcyjA40042) and by the National Natural Science Foundation of China (Grant No. 61301033).

References

Shalom

Y. B.

X. R.

, Multitarget-Multisensor Tracking: Principles and Techniques, Storrs, CT: YBS Publishing, 1995.

Dana

M. P.

, “Registration: A prerequisite for multiple sensor tracking,™ in Multiarget- Multisensor Tracking: Advanced Applications, Bar-Shalom

, Ed. Artech House, Norwood, MA, 1990.

Zhou

Leung

Yip

P.C.

, “An exact maximum likelihood registration algorithm for data fusion,™ IEEE Transactions on Signal Processing, vol. 45, no. 6, pp. 1560–1573, 1997.

Nabaa

Bishop

R. H.

, “Solution to a multisensor tracking problem with sensor registration errors,™ IEEE Transactions on Aerospace and Electronic Systems, vol. 35 no.1, pp. 354–363, 1999.

Huang

Leung

, “An expectation-maximization-based interacting multiple model approach for cooperative driving systems,™ IEEE Transactions on Intelligent Transportation Systems, vol. 6, no. 2, pp. 206–228, 2005.

Huang

Leung

Bosse

, “A *pseudo-measurement approach to simultaneous registration and track fusion,™ IEEE Transactions on Aerospace and Electronic Systems, vol. 48, no. 3, pp. 2315–2331, 2012.

Okello

Challa

, “Joint sensor registration and track-to-track fusion for distributed trackers,™ IEEE Transactions on Aerospace and Electronic Systems, vol. 40, no.3, pp. 808–823, 2004.

Leung

Zhou

, “Space-time registration of radar and ESM using unscented Kalman filter,™ IEEE Transactions on Aerospace and Electronic Systems, vol. 40, no. 3, pp. 824–836, 2004.

Smidl

Quinn

, The Variational Bayes Method in Signal Processing. Berlin, Germany: Springer, 2005.

10.

Hao

Leung

Zhongshi

, “State estimation in unknown non-Gaussian measurement noise using Variational Bayesian Technique,™ IEEE Transactions on Aerospace and Electronic Systems, vol. 49, no. 4, pp. 2601–2614, 2013.

11.

Hao

Leung

Zhongshi

, “A variational Bayesian approach to robust sensor fusion based on Student-t distribution,™ Information Sciences, vol. 221, pp. 201–214, 2013.

12.

Roberts

S. J.

Penny

W. D.

, “Variational Bayes for generalized autoregressive models,™ IEEE Transactions on Signal Processing, vol. 50, no. 9, pp. 2245–2257, 2002.

13.

Yan

, “Bayesian nonparametric models for multiway data analysis,™ IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 37, no. 2, pp. 475–487, 2015.

14.

Nasions

Bors

A. G.

, “Variational learning for Gaussian mixture models,™ IEEE Transactions on Systems, Man, and Cybernetics, Part B: Cybernetics, vol. 36, no. 4, pp. 849–862, 2006.

15.

Smidl

Quinn

, “Variational Bayesian Filtering,™ IEEE Transactions on Signal Processing, vol. 56, no. 10, 5020–5030, 2008.

16.

Williams

J. L.

, “An efficient, variational approximation of the best fitting multi-Bernoulli filter,™ IEEE Transactions on Signal Processing, vol. 63, no. 1, 258–273, 2015.

17.

Sarkka

Nummenmaa

, “Recursive noise adaptive Kalman filtering by variational Bayesian approximations,™ IEEE Transactions on Automatic Control, vol. 54, no. 3, pp. 596–600, 2009.

18.

Polatkan

Zhou

Carin

Blei

Daubechies

, “A Bayesian nonparametric approach to image super-resolution,™ IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 37, no. 2, pp. 346–358, 2015.

19.

Tzikas

D. G.

Likas

A. C.

Galatsanos

N. P.

, “Variational Bayesian sparse kernel-based blind image deconvolution with Student's-t priors,™ IEEE Transactions on Image Processing, vol. 18, no. 4, pp. 753–764, 2009.

20.

Zheng

Y. F.

, “Patch-based video processing: A variational Bayesian approach,™ IEEE Transactions on Circuits and Systems for Video Technology, vol. 19, no. 1, pp. 27–40, 2009.

21.

Yik-Chung

Shaoqian

, “The variational inference approach to joint data detection and phase noise estimation in OFDM,™ IEEE Transactions on Wireless Communications, vol. 12, no. 10, pp. 5309–5321, 2013.

22.

Lin

D. D.

Lim

T. J.

, “A variational inference framework for soft-in soft-out detection in multiple-access channels,™ IEEE Transactions on Information Theory, vol. 55, no. 5, pp. 2345–2364, 2009.

23.

Watanabe

Minami

Nakamura

Ueda

, “Variational Bayesian Estimation and Clustering for Speech Recognition,™ IEEE Transactions on Speech and audio processing, vol.12, no.4, pp. 365–381, 2004.

24.

Paisley

Wang

Blei

D. M.

Jordan

M. I.

, “Nested Hierarchical Dirichlet Processes,™ IEEE Transactions on Pattern Analysis and Machine Intelligence, vol.37, No.2, pp.256–270, 2015.

25.

Frey

B. J.

Jojic

, “A Comparison of Algorithm for Inference and Learning in Probabilistic Graphical Models,™ IEEE Transactions on Pattern Analysis and Machine Intelligence, vol.27, No.9, pp. 1392–1416, 2005.

26.

Beal

, Variational algorithms for approximate Bayesian inference, University of London, Ph.D. Thesis, 2003.

27.

Tzikas

D. G.

Likas

A. C.

Galatsanos

N. P.

, “The Variational Approximation for Bayesian Inference,™ IEEE Signal Processing Magazine, vol.25, no. 6, pp.131–145, 200.

Recursive Variational Bayesian Inference to Simultaneous Registration and Fusion

Abstract

Keywords

1. Introduction

2. Problem Formulation

3. RVBEM for Simultaneous Registration and Fusion Recursively

4. Simulations

5. Conclusion

Footnotes

6. Acknowledgements

References