Accurate and Robust Localization in Harsh Environments Based on V2I Communication

1. Introduction

With the arrival of global navigation satellite systems (GNSS), in-car navigation has increasingly become an essential tool for the automotive industry. However, the performance of GNSS is compromised in harsh environments where there is not a line of sight (LOS) to satellites, e.g., tunnels, covered parking areas and dense urban canyons [1]. Hence, in-car navigation requires a localization technology that operates with robustness in such circumstances. The development of vehicular ad-hoc networks (VANETs) provides a promising platform to fulfill this requirement [2].

In VANETs, an on-board unit (OBU) inside the vehicle communicates with other OBUs or with stationary roadside units (RSUs), in vehicle-to-vehicle (V2V) and vehicle-to-infrastructure (V2I) communications, respectively [3]. Cooperation between OBUs can provide good position estimates in V2V communication [4-5]. However, the quick topology changes required by V2V approaches make V2I communication be the preferred option for in-car navigation in harsh environments [6]. In V2I communication, the position of an OBU (the target) can be estimated from range-related measurements taken on the radio-frequency signals transmitted to and from the RSUs (the anchors) [7]. However, the changeable and unpredictable characteristics of the wireless channel in harsh environments make multipath and non-line of sight (NLOS) propagation conditions be predominant [8-9]. Therefore, conventional positioning systems designed for tractable and static signal behavior cannot guarantee an adequate performance.

The position information extracted from the radio-frequency signals varies according to the type of measurement taken. Techniques based on time of arrival (TOA) [9-10] or received signal strength (RSS) [11-12] measurements obtain range-related information, whereas techniques based on angle of arrival (AOA) or time difference of arrival (TDOA) measurements extract information related to directions or difference of distances, respectively [13-14]. AOA and TDOA measurements entail significant costs of antenna-array integration or synchronizing devices. In this chapter, we focus on RSS and TOA measurements that can provide accurate localization with an appropriate complexity.

In the TOA case, measuring the round-trip time avoids the technical difficulty of time synchronization among the nodes.

Range or position estimation is an inference problem where the observations are the RSS and TOA measurements [15-16]. From a Bayesian perspective, determining the posterior distribution of ranges or positions from observations is the optimal approach [17-25]. Then, ranges or positions can be obtained by means of the maximum a posteriori (MAP) or the minimum mean square error (MMSE) estimators.

The optimality of the above mentioned methods depends on the fit between the model assumed for the relationship between measurements and ranges or positions (i.e., the likelihood function) and the actual behavior of the measurements. Tractable and static models for the likelihoods based on Gaussian distributions accurately explain the behavior of measurements only in open areas [26-28]. For harsh environments, several techniques have been developed to address the complex behavior of wireless signal metrics. In the TOA case, the NLOS bias causes range overestimation. Thus, a common procedure is to detect and remove NLOS measurements [29]; other techniques utilize prior knowledge about this NLOS error to subtract it and adjust the measurements to their LOS values [27,30]. In the RSS case, the performance depends on the estimation of the parameters that characterize the propagation channel at each time [12,26]. Certain approaches deal with the dynamic nature of RSS metric through fingerprinting or machine learning [11,21,31]; however, their accuracy is sensitive to fast environmental changes and they do not fuse different signal metrics.

Range and position estimation can be improved by exploiting the relationship among positions in time through Bayesian filtering. Kalman filtering techniques rely on Gaussian models that are not adequate for harsh environments. Different alternative methods based on variations of such filters, as well as on particle filters (PFs), have been proposed: low complexity non-linear/non-parametric adaptive modeling is used for filtering of RSS fingerprints in [11,21]; recursive Bayesian estimation together with multipath and NLOS propagation effects are considered in [22-23]; TOA and RSS data fusion is performed in [32-34]; hybrid information is exploited by particle filtering in [24]; and RSS/TOA Bayesian fusion for multipath and NLOS mitigation are performed in [25]. However, these methods require prior information achieved by arduous training phases or rely on assumptions non-realistic for harsh environments, such as Gaussian and static models.

This chapter presents a framework for adaptive data fusion to handle the difficulties described above, based on non-parametric dynamic modeling of the likelihood. The subsequent usage of a PF leads to the adaptive likelihood particle (ALPA) filter. As we show, the estimation can be carried out without requiring any calibration stage, thus enabling localization capabilities to pre-existing wireless infrastructures, such as VANETs based on V2I communication. The main contributions of this chapter are as follows:

• We present techniques for adaptive and systematic modeling of the relationship between measurements and positions, by means of a dynamic and empirical likelihood function.

• We present a model for Bayesian fusion of TOA and RSS measurements, based on nonlinear and non-Gaussian Bayesian filtering and the likelihoods derived over time.

• We show the suitability of the proposed techniques by experimentation performed using common wireless local area network (WLAN) devices.

• We show the near-optimality of the method by comparing its performance to the posterior Cramér-Rao lower bound (CRLB).

Both empirical and simulation results show that the proposed methods significantly improve the accuracy of conventional approaches with an important reduction on the number of measurements needed.

The structure of the rest of this chapter is as follows: Section II defines the position estimation problem; Section III addresses this problem under a hidden Markov model (HMM) and defines the dynamic and measurements models; Section IV presents the adaptive data fusion technique for likelihood modeling and the recursive Bayesian approach for solving the resulting non-linear and non-Gaussian problem; Section V shows the experimental and simulation results; Section VI includes a discussion on complexity; and finally, Section VII draws the conclusions.

Notations: The notation $p\left(\mathbf{x}\right)$ is the probability density function (pdf) of the random variable $\mathbf{x}$; $f^{\left(m\right)}\left(\bullet \right)$ denotes the $m$th derivative of a real function $f$ evaluated in its argument; $f\left[k\right]$ for $k\in \mathbb{N}$ denotes the value of the function $f$ evaluated in $t_k\in \mathbb{R}$; $\mathbf{X}\left[k\right]$ denotes the set $\{\mathbf{x}\left[i\right],i=1,\cdots ,k\}$; if $M$ is a positive integer, ${\mathcal{M}}^M$ denotes the $M$ Cartesian power of $\left\{1,\dots ,M\right\}$; finally, $\stackrel{-}{\mathbf{z}}$ denotes the sample mean of the components of a vector $\mathbf{z}$.

2. Problem statement

In the following, we consider a two-dimensional scenario where a mobile target (e.g., a car equipped with an OBU) moves freely. To determine its position, the target communicates with several anchors (the RSUs) with known positions. Since the localization system can get measurements in discrete times $\{t_k,k\in \mathbb{N}\}$, we are interested in estimating the sequence $\left\{\mathbf{x}\left[k\right],k\in \mathbb{N}\right\}$ from a sequence of measurements $\left\{\mathbf{z}\left[k\right],k\in \mathbb{N}\right\}$. The entries of vector $\mathbf{x}\left[k\right]$ can be the distances between the target and each anchor or the coordinates of the mobile target’s position. The entries of vector $\mathbf{z}\left[k\right]$ are RSS and TOA measurements.

Next section establishes the probabilistic relationship between vectors $\mathbf{x}\left[k\right]$ and $\mathbf{z}\left[k\right]$ by modeling this problem as an HMM, and defining the state vector, $\mathbf{y}\left[k\right]$ , that consists of vector $\mathbf{x}\left[k\right]$ and several of its derivatives.

3. Hidden Markov model

In addition to the information conveyed by the measurements, the fact that the sequence $\left\{\mathbf{x}\left[k\right],k\in \mathbb{N}\right\}$ is highly correlated in time can likewise be used as another source of information. The position of the target cannot change abruptly in a small lapse of time; hence, we can model the evolution in time of positions or Euclidean distances to each anchor as an analytic function. Being $x\left(t\right)$ a component of the position or the distance to an anchor, we can approximate its value in $t_{k+1}$ by using the $n$th-order Taylor expansion in $t_k$,

$$x\left[k+1\right]\approx x\left[k\right]+x^{\mathrm{\text{'}}}\left[k\right]∆t+x^{\mathrm{\text{'}}\mathrm{\text{'}}}\left[k\right]\frac{{∆t}^2}{2}+\cdots +x^{\left(n\right)}\left[k\right]\frac{{∆t}^n}{n!}$$E1

where $∆t=(t_{k+1}-t_k)\in \mathbb{R}$ is the sampling interval. The error in this approximation is

$x^{\left(n+1\right)}\left({\xi }_0\right)\frac{{∆t}^{n+1}}{\left(n+1\right)!}$

and ${\xi }_0\in \mathbb{R}$ is a point in the interval ${[t}_k,t_{k+1}]$. Therefore, the error in the approximation (1) depends directly on $∆t$, on the smoothness of $x\left(t\right)$ (represented by the $(n+1)$th derivative), and on the order of the approximation.

The correlation in time expressed in (1) implies that $\left\{\mathbf{x}\left[k\right],k\in \mathbb{N}\right\}$ is not a Markov chain, i.e., the current distance or position depends not only on the previous one. However, calling $\mathbf{y}\left[k\right]$ the positional-state consisting of the distance or the position and its $n$ derivatives, we can assume that $\left\{\mathbf{y}\left[k\right],k\in \mathbb{N}\right\}$ is a Markov chain. Moreover, we can likewise assume that, conditioned on $\left\{\mathbf{y}\left[k\right],k\in \mathbb{N}\right\}$, $\left\{\mathbf{z}\left[k\right],k\in \mathbb{N}\right\}$ is a sequence of independent random variables, i.e., given the current positional-state vector, the measurements $\mathbf{z}\left[k\right]$ are independent of all previous and future positional-states and measurements [18]. These assumptions let to build an HMM in which the positional-state vectors $\left\{\mathbf{y}\left[k\right],k\in \mathbb{N}\right\}$ form a non-observable Markov chain, and what is available is the other stochastic process $\left\{\mathbf{z}\left[k\right],k\in \mathbb{N}\right\}$, linked to the Markov chain in that $\mathbf{y}\left[k\right]$ governs the distribution of $\mathbf{z}\left[k\right]$ [35] (see Figure 1).

Figure 1.

Hidden Markov Model for positional-states and measurements evolution. The relationship between $\mathbf{y}\left[k\right]$ and $\mathbf{y}\left[k-1\right]$ and the relationship between $\mathbf{z}\left[k\right]$ and $\mathbf{y}\left[k\right]$ are the only two kinds of dependence.

The conditional independence assumptions reflected in Figure 1 lead to two kinds of dependence between the random variables [36],

• Dynamic model: establishes the relationship between the state vector in time $t_k$ and the state vector in time $t_{k-1}$, i.e., $p\left(\mathbf{y}\left[k\right]\right|\mathbf{y}\left[k-1\right])$.

• Measurements model: establishes the relationship between the measurements and the state vector in each time, i.e., $p\left(\mathbf{z}\left[k\right]\right|\mathbf{y}\left[k\right])$.

Then, the joint distribution of all the random variables involved in the process is given by,

This probabilistic model is a generalization of the maximum likelihood approach in which the estimation is accomplished for a given time instant, $t_k$, neither considering previous nor future positional-states and measurements. In this case, $p(\mathbf{y}\left[k\right],\mathbf{z}\left[k\right])\propto p\left(\mathbf{z}\left[k\right]\right|\mathbf{y}\left[k\right])$

$\begin{array}{l}p\left(Y\left[k\right],Z\left[k\right]\right)=p\left(y\left[1\right]\right)p\left(z\left[1\right]\text{|}y\left[1\right]\right){\displaystyle \prod }_{\text{i}=2}^{\text{k}}p\left(y\left[i\right]\text{|}y\left[i-1\right]\right)p\left(z\left[i\right]\text{|}y\left[i\right]\right)\\ =p(Y\left[k-1\right],Z\left[k-1\right])p\left(y\left[k\right]\text{|}y\left[k-1\right]\right)p\left(z\left[k\right]\text{|}y\left[k\right]\right)\end{array}$E2

The modeling as an HMM shown in (2) makes possible to infer the posterior distribution $p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k\right]\right)$ through a recursive process. In the specific case where the dynamic and measurements model are linear-Gaussian, the posterior distribution is also Gaussian, and the Bayesian inference can be optimally performed by the celebrated Kalman filter (KF) [19]. In the following, we describe these models for harsh propagation environments, showing that the dynamic model can be assumed linear with a wide generality, whereas this assumption for the measurements model yields inaccurate performances.

3.1. Dynamic model

The dynamic model of the positional-state vector can be obtained from the evolution in time given by (1), and by approximating each $m$th derivative, for $m=1,\dots ,n$, by its $\left(n-m\right)$th-order Taylor expansion, as

$$\mathbf{y}\left[k+1\right]={\mathbf{F}}_k\mathbf{y}\left[k\right]+{\mathbf{n}}_d\left[k\right]$$E3

where

$${\mathbf{F}}_k=\left(\begin{array}{ccccc}1& ∆t& \frac{{∆t}^2}{2}& \cdots & \frac{{∆t}^n}{n!}\\ 0& 1& ∆t& \cdots & \frac{{∆t}^{n-1}}{(n-1)!}\\ ⋮& \ddots & \ddots & \ddots & ⋮\\ 0& \cdots & 0& 1& ∆t\\ 0& \cdots & \cdots & 0& 1\end{array}\right)$$E4

is the transition matrix, and ${\mathbf{n}}_d\left[k\right]$ is the error in the approximations. For example, in the case of estimating a one-dimensional parameter, $x\left[k\right]$, the error ${\mathbf{n}}_d\left[k\right]$ is given by

$${\mathbf{n}}_d\left[k\right]=\left(\begin{array}{c}\frac{{∆t}^{n+1}}{\left(n+1\right)!}{x}^{\left(n+1\right)}\left({\xi }_0\right)\\ \frac{{∆t}^n}{n!}{x}^{\left(n+1\right)}\left({\xi }_1\right)\\ ⋮\\ ∆t{x}^{\left(n+1\right)}\left({\xi }_{\mathrm{n}}\right)\end{array}\right)$$E5

where ${\xi }_0,\dots ,{\xi }_{\mathrm{n}}$ are values in the interval ${[t}_k,t_{k+1}]$. The values taken by the $(n+1)$ th derivative of $x\left(t\right)$ in the unknown points ${\xi }_0,\dots ,{\xi }_{\mathrm{n}}$ are modeled as realizations of a random variable that can be assumed to be zero-mean Gaussian variable with a standard deviation ${\sigma }_{d^{(n+1)}}$ [18-19]. Then, we can model the evolution in time of $x\left[k\right]$ as a random walk. Therefore, the dynamic model is a discrete Wiener process velocity (DWPV) model or a discrete Wiener process acceleration (DWPA) model if we use the second- or third-order Taylor expansion, respectively [18]. Hence, the dynamic model, $p\left(\mathbf{y}\left[k\right]|\mathbf{y}\left[k-1\right]\right)$, can be assumed linear-Gaussian.

3.2. Measurements model

The second ingredient to characterize the HMM is the measurements model or likelihood, $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$. This probability distribution relates the measurements to the positional-state. In the case of range-related measurements, we have that $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)=p\left(\mathbf{z}\left[k\right]|d\left[k\right]\right)$, irrespectively of the positional-state used. In the following, we describe realistic models for the relationship between distances and RSS/TOA measurements in concordance with previous essays [10,12].

3.2.1. RSS measurements

In a given specific instant and place, the RSS values are affected by the distance between emitter and receiver. The attenuation caused by the distance between two nodes is known as path-loss and is proportional to this distance raised to a certain exponent, called path-loss exponent [7,12,15,26]. However, the RSS values are likewise affected by a wide range of unpredictable factors, such as multipath propagation (fast fading) and shadowing (slow fading) [37]. By reflecting these factors in the Friis transmission equation for free-space, the relationship between the received signal strength, $P_r$, and the distance, $d\left[k\right]$, is given by [26],

$$P_r=\frac{G_t{G}_r}{4\pi }{P}_t\frac{g^2\gamma }{{\left(d\right[k\left]\right)}^{{\beta }_s}}$$E6

where $P_t$ is the transmitted power, $G_t$ and $G_r$ are the transmitter and receiver gains, respectively, $g$ and $\gamma$ are the parameters of the Rayleigh/Rician and log-normal distributions that model the fast and slow fading, respectively, and ${\beta }_s$ is the path-loss exponent corresponding to the specific propagation environment [37].

By following the procedure described in [26] and taking logarithmic units, we obtain the measurements model for RSS values,

$$z_s\left[k\right]={\alpha }_s-10{\beta }_s{\log }_{10}\left(d\left[k\right]\right)+n_s\left[k\right]$$E7

where $z_s\left[k\right]\in \mathbb{R}$ is the RSS measured value and ${\alpha }_s$ a constant that depends on $P_t$, $G_t$, $G_r$ and the fast and slow fading [12,15,26]. Finally, $n_s\left[k\right]$ is a noise term caused by shadowing that has zero mean in cases where the parameters ${\alpha }_s$ and ${\beta }_s$ fit perfectly the current propagation conditions [12,15,26]. In practice, the value of ${\alpha }_s$ can be previously known [26]. However, in realistic scenarios, the path-loss exponent, ${\beta }_s$ , used to relate RSS values to distances, does not fit exactly the actual propagation conditions [12], and hence, the noise term, $n_s\left[k\right]$, will have a non-zero mean proportional to the logarithm of the distance.

3.2.2. TOA Measurements

The distance between emitter and receiver also affects the time taken by the signal to be propagated from one node to the other. By assuming known the signal speed, we can infer this distance by means of a linear transformation of the TOA values. Due to the technical difficulty of synchronizing devices in a wireless network, techniques that use round-trip time estimation are the most attractive to estimate delays [10,28]. In this case, the processing time at the device that has to transmit the echo causes the relationship between TOA and distance to be affine linear (it has an intercept term). Then, we can model the relationship between the delay, $z_{\tau }\left[k\right]$, measured at time $t_k$, and the distance at that time, $d\left[k\right]$, as,

$$z_{\tau }\left[k\right]={\alpha }_{\tau }+{\beta }_{\tau }d\left[k\right]+n_{\tau }\left[k\right]$$E8

where ${\alpha }_{\tau }$ and ${\beta }_{\tau }$ are constants that can be estimated by a linear regression of measurements previously obtained [28,38-39]. The term $n_{\tau }\left[k\right]$ models the noise that is ussually assumed to be zero-mean and Gaussian in case of LOS propagation. However, in case of NLOS propagation, it is currently not known how to accurately model such error term, where several statistical distributions taking positive values, such as Exponential, Rayleigh, Weibull or Gamma, have been used in the literature [26-28].

From the above discussion, we can notice that in all cases the expected value of the measurements is $\mathrm{{\rm E}}\left\{z\right\}=f\left(d\left[k\right]\right)+b$, where $f$ is a linear or logarithmic function, and $b$ is a systematic error in the model. In addition, we can point out that in harsh environments:

1. the relationship between measurements related to distances and distances is nonlinear and non-Gaussian;

2. such relationship highly depends on the propagation environment that can change rapidly.

These two factors render the linear-Gaussian assumption inadequate for the measurements model, $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$, in harsh environments. Therefore, common inference techniques that use naive and static models may obtain poor results in realistic dense cluttered scenarios.

4. Bayesian adaptive RSS/TOA fusion

Conventional non-Bayesian approaches for parameter estimation are based on maximum-likelihood (ML) estimation (in our case the maximization of $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$). ML commonly assumes tractable models for the likelihood (e.g., Gaussian likelihoods yield a least squares problem), while more intricate models are usually solved by means of expectation-maximization (EM) algorithm [40-41]. In the event that certain prior information about the parameter of interest is available, we can achieve a better estimator by adding this new information. If this prior information is the correlation in time of positional-states, it can be exploited through sequential Bayesian inference. In the following, we briefly describe such estimation process and present the adaptive likelihood particle (ALPA) filter for Bayesian inference based on RSS and TOA non-parametric adaptive likelihoods.

4.1. Bayesian inference

In the above mentioned context, the task is to determine the posterior distribution of positional-states given the measurements, $\mathbf{Z}\left[k\right]$, from the knowledge of the prior, $p\left(\mathbf{y}\left[k\right]\right)$, and the likelihood, $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$, by using the Bayes’ rule [19,42]. The knowledge about the prior distribution, $p\left(\mathbf{y}\left[k\right]\right)$, can come from several avenues, e.g., from environmental knowledge. In this chapter, we use as prior knowledge the positional-states inferred in previous instants over the framework offered by the HMM above explained (see Figure 1). However, any other kind of prior information can be incorporated analogously.

In the case of modeling the positional-state and measurements evolution as an HMM, the expression (2) provides a way to determine the posterior distribution iteratively,

$p\left(\mathbf{Y}\left[1\right]|\mathbf{Z}\left[1\right]\right)=\frac{p\left(\mathbf{y}\left[1\right],\mathbf{z}\left[1\right]\right)}{p\left(\mathbf{z}\left[1\right]\right)}=\frac{p\left(\mathbf{y}\left[1\right]\right)p\left(\mathbf{z}\left[1\right]|\mathbf{y}\left[1\right]\right)}{p\left(\mathbf{z}\left[1\right]\right)}$

and for $k>1$,

$\begin{array}{l}p\left(Y\left[k\right]\text{|}Z\left[k\right]\right)=\frac{p\left(Y\left[\text{k}\right],Z\left[\text{k}\right]\right)}{p\left(Z\left[\text{k}\right]\right)}\\ =\frac{p\left(z\left[k\right]\text{|}y\left[k\right]\right)p\left(y\left[k\right]\text{|}y\left[k-1\right]\right)p\left(Y\left[\text{k}-1\right],Z\left[\text{k}-1\right]\right)}{p\left(Z\left[\text{k}\right]\right)}\\ =\frac{p\left(z\left[k\right]\text{|}y\left[k\right]\right)p\left(y\left[k\right]\text{|}y\left[k-1\right]\right)p\left(Y\left[\text{k}-1\right]|Z\left[\text{k}-1\right]\right)}{p\left(z\left[\text{k}\right]|z\left[\text{k}-1\right]\right)}\\ \end{array}$E9

From the posterior distribution, $p\left(\mathbf{Y}\left[k\right]|\mathbf{Z}\left[k\right]\right)$, we can estimate $\mathbf{y}\left[k\right]$ by,

$$p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k\right]\right)=\int p\left(\mathbf{Y}\left[k\right]|\mathbf{Z}\left[k\right]\right)\mathrm{d}\mathbf{Y}\left[k-1\right]$$E10

leading to a process called filtering.

The positional-state $\mathbf{y}\left[k\right]$ can likewise be estimated by using the measurements until time $t_{k+l}$, leading to a process called smoothing if $l>0$ or prediction if $l<0$.

By replacing (9) in (10) we obtain,

$$p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k\right]\right)=\frac{p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)\int p\left(\mathbf{y}\left[k\right]|\mathbf{y}\left[k-1\right]\right)p\left(\mathbf{Y}\left[\mathrm{k}-1\right]|\mathbf{Z}\left[\mathrm{k}-1\right]\right)\mathrm{d}\mathbf{Y}\left[k-1\right]}{p\left(\mathbf{z}\left[\mathrm{k}\right]|\mathbf{z}\left[\mathrm{k}-1\right]\right)}$$E11

By assuming known the posterior distribution at $t_{k-1}$, $p\left(\mathbf{Y}\left[k-1\right]|\mathbf{Z}\left[k-1\right]\right)$, we can perform the filtering process in two steps [19],

1. Prediction: from the dynamic model we obtain the prediction of the positional-state in time$t_k$, given the measurements until time $t_{k-1}$,

$$p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k-1\right]\right)=\int p\left(\mathbf{y}\left[k\right]|\mathbf{y}\left[k-1\right]\right)p\left(\mathbf{Y}\left[\mathrm{k}-1\right]|\mathbf{Z}\left[\mathrm{k}-1\right]\right)\mathrm{d}\mathbf{Y}\left[k-1\right]$$E12

1. Update: from the measurements model we correct the prediction when a new set of measurements, $\mathbf{z}\left[k\right]$, is available in time $t_k$,

$$p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k\right]\right)=\frac{p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)p\left(\mathbf{y}\left[k\right]|\mathbf{Z}\left[k-1\right]\right)}{p\left(\mathbf{z}\left[\mathrm{k}\right]|\mathbf{z}\left[\mathrm{k}-1\right]\right)}$$E13

and the normalization constant,

$$p\left(\mathbf{z}\left[\mathrm{k}\right]|\mathbf{z}\left[\mathrm{k}-1\right]\right)=\int p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)p\left(\mathbf{y}\left[\mathrm{k}\right]|\mathbf{Z}\left[\mathrm{k}-1\right]\right)\mathrm{d}\mathbf{y}\left[k\right]$$E14

Hence, the objective is to infer the hidden positional-state vector in each time, $\mathbf{y}\left[k\right]$, by using the information achieved by the measurements and the relationship between the variables in time. The Bayesian recursive process given by (12) and (13) avoids the need of reprocessing all the stored data since the posterior distributions are obtained iteratively. Figure 2 graphically explains the evolution of the distributions involved in the filtering process, for the problem of estimating the range between the OBU and an RSU, and for the problem of estimating the position of the OBU when it communicates with three RSUs.

In order to perform the described filtering process, we need the likelihood function of the measurements $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$. This function is a priori unknown in harsh environments, since the distribution of the error term in the measurements model is highly environmental-dependent and varies rapidly with time. In the RSS case, although the error term is usually assumed to be zero-mean Gaussian distributed, this assumption is too naive in realistic scenarios where, for example, only one estimation of the path-loss exponent is available [12,26]. For TOA measurements, this error term has been modeled with several parametric distributions such as Gaussian, Exponential, Gamma or Rayleigh [26-27,43] or by means of specific distributions obtained in each particular propagation environment [28,44]. In the following sections, we propose an adaptive likelihood function for data fusion that dynamically adjusts to the changing propagation conditions from the nature of the measurements collected in real time.

Figure 2.

Density functions involved in filtering process for range and position estimation (darker zones have higher probability): (a) the target with the OBU moves in $t_k$ with respect to its position in $t_{k-1}$; (b) the posterior density in $t_{k-1}$ is known; (c) from the dynamic model we perform the prediction; (d) in $t_k$ the target receives a new set of measurements; (e) from the likelihood we update the prediction to obtain the posterior density in $t_k$.

4.2. Adaptive likelihood for RSS/TOA fusion

The sets of RSS and TOA measurements obtained in each instant consist of samples from the random variable $z_s\left[k\right]$ and $z_{\tau }\left[k\right]$, respectively. As we show below, it is possible to represent the likelihood function in each instant and environment by using the set of samples through a non-parametric representation based on kernels [11,45-46].

A kernel function is a symmetric function (not necessarily positive) whose integral over the entire space is equal to one.

After the reception of $M$ RSS or $M$ TOA measurements $\{z_k^i,i=1,\dots ,M\}$, we can approximate the pdf of $z_s\left[k\right]$ or $z_{\tau }\left[k\right]$ as

$$p\left(z\right)\approx \frac{1}{Mh}\sum _{i=1}^{M}K\left(\frac{z-z_k^i}{h}\right)$$E15

where $K(\bullet )$ is the kernel function and $h$ is a positive number called bandwidth [11,45-46]. Several functions can be chosen for the kernel, where the most common is to use the standard Gaussian kernel [47], i.e.,

$$K\left(x\right)=\frac{1}{\sqrt{2\pi }}{\mathrm{e}}^{-\frac{1}{2}{x}^2}$$E16

By assuming that the distribution of the measurements $z$, has the expression (15) in time $t_k$, we can obtain the likelihood relating distances to measurements in each instant $k$ as the following result shows.

Proposition 1. Let $\mathbf{z}\left[k\right]=\{z_k^i,i=1,\dots ,M\}$ be a set of measurements (samples of $z$) related to the distance $d\left[k\right]$ by a model $\mathrm{ }\mathrm{{\rm E}}\left\{z\right\}=f\left(d\left[k\right]\right)+b$. Then, assuming $z$ follows the distribution given by (15), and calling $\mathrm{ }{\varsigma }_{i,j}=z_k^j-z_k^i+\stackrel{-}{\mathbf{z}\left[\mathrm{k}\right]}$, the likelihood function of the measurements is

$$p\left(\mathbf{z}\left[k\right]|d\left[k\right]\right)=\frac{1}{{\left(2\pi \right)}^{M/2}{\left(Mh\right)}^M}\sum _{\left(i_1,\dots ,i_M\right)\in {\mathcal{M}}^{\mathrm{M}}}{\mathrm{{\rm E}}}_{\mathrm{b}}\{\mathrm{e}\mathrm{x}\mathrm{p}(\frac{-1}{2h^2}\sum _{j=1}^M{({\varsigma }_{i,j}-f(d\left[k\right])-b)}^2\left)\right\}$$E17

where the expectation ${\mathrm{{\rm E}}}_{\mathrm{b}}\{\bullet \}$ is taken with respect to the systematic errors, $b$, in the model.

Proof: see [48].

The Proposition 1 enables to obtain individual likelihoods from a set of measurements. Data fusion from different signal metrics (i.e., RSS and TOA) is carried out by combining these likelihoods. Let ${\mathbf{z}}_s\left[k\right]$ and ${\mathbf{z}}_{\tau }\left[k\right]$ be sets of RSS and TOA measurements, respectively, forming the set of measurements obtained in the instant $k$. Then, assuming that, given the real distance, $d\left[k\right]$, ${\mathbf{z}}_s\left[k\right]$ and ${\mathbf{z}}_{\tau }\left[k\right]$ are independent, we have that,

$$p\left(\mathbf{z}\left[k\right]|d\left[k\right]\right)=p\left({\mathbf{z}}_s\left[k\right]|d\left[k\right]\right)p\left({\mathbf{z}}_{\tau }\left[k\right]|d\left[k\right]\right)$$E18

where the likelihood of each kind of measurement can be dynamically obtained from (17).

In order to describe how the presented adaptive data fusion operates, Figures 3-4 show the histogram of 100 RSS and 100 TOA measurements taken at a fixed distance with the measuring systems described in [12] and [10], respectively. These figures also represent the corresponding Gaussian pdf and the adaptive pdf obtained by means of the kernel-based expression given by (15) and (16).

In Figures 3-4 and in the following, we use a fixed bandwidth of one-half of the resolution of the measuring system [10,12]. This election avoids both undersmoothed curves with too much spurious data artifacts, and oversmoothed densities that obscure the underlying nature of the measurements [46].

From those figures, we can point out that, despite the fact that the true density is unknown, the presented adaptive pdf can express the dynamic behavior of RSS/TOA measurements in harsh environments with better accuracy than histogram and Gaussian density estimates [49-50].

Figure 3.

The adaptive density accurately approximates the complex randomness of RSS measurements in harsh environments.

Figure 4.

The adaptive density accurately approximates the complex randomness of TOA measurements in harsh environments.

In Figure 5 we illustrate the RSS/TOA data fusion process by representing the adaptive likelihood function obtained by means of expressions (17) and (18).

In Figure 5 and in the following sections, we use coarse models for the measurements biases in accordance with previous essays [10,12]. Specifically, the RSS bias is modeled as a Gaussian $\mathcal{N}(0,{\sigma }_s)$ with ${\sigma }_s=3$ dBm, and the TOA bias as a Uniform distribution $\mathcal{U}(0,{\gamma }_{\tau })$ with ${\gamma }_{\tau }=4$ clock cycles.

Figure 5.

The adaptive RSS/TOA data fusion, defined by Proposition 1 and (18), results, in this case, in an improvement of 0.5 meters in ML estimator compared with the Gaussian case, which is equivalent to a reduction of 18% of the error.

From Figure 5, we can point out that the adaptive likelihood function provides more information about the distance than the Gaussian model, by combining the individual adaptive likelihoods obtained with RSS and TOA measurements. Moreover, the height of both functions reflects the more reliable information obtained by adaptive estimation. From that figure, we also observe the improvement achieved by means of data fusion with respect to the individual estimates. This likelihood function leads to the ALPA filter defined in the following section.

4.3. Adaptive likelihood particle filter

Within the framework provided by the HMM, if both dynamic and measurements models are linear-Gaussian, all the posterior distributions are also Gaussian. In this case, all the involved density functions are completely described by their mean vectors and covariance matrices, obtained by a KF [19]. In the case of interest in this chapter, the models in the HMM are neither linear nor Gaussian, and then, the usage of KFs is suboptimal. In order to circumvent this drawback, the classical solution consists of using extended KFs (EKF) [23,25]. However, better performances can be obtained by PFs that let the usage of more general and flexible models [17,19] as the adaptive likelihood described in the previous section.

A PF represents the posterior distribution through a discrete distribution, where the support points and their probabilities are called particles and weights, respectively. To estimate the posterior distribution, we need to iteratively obtain a certain number of samples (particles) and probabilities (weights) capable of representing the posterior distribution. These particles and weights can be obtained by a method known as sequential-importance-sampling (SIS) [19,51], where the weight of the different particles can be determined by evaluating the likelihood function pointwise. Therefore, more realistic models such as the presented adaptive likelihood function for data fusion can be used, leading to the ALPA filtering algorithm describe in Table 1.

i. Initialization: ∙ Initial particles: draw $N$samples $\{{\mathbf{y}}_1^i,i=1,\dots ,N\}$ from the known density function$p\left(\mathbf{y}\right[1\left]\right)$.; ∙ Initial weights:${\omega }_1^i=\frac{1}{N},i=1,\dots ,N$.; ii. Recursive estimation: for $k>1$,; ∙ Particles in instant $k$from particles in instant$k-1$: draw $N$samples $\{{\mathbf{y}}_k^i,i=1,\dots ,N\}$ from the proposal distribution$q\left(\mathbf{y}\right[k\left]\right|{\mathbf{y}}_{k-1}^i,\mathbf{z}\left[k\right])$.; ∙ From RSS measurements and Proposition 1, evaluate the weight of each particle. For $i=1,\dots ,N$

${\tilde{\omega }}_s^i=p\left({\mathbf{z}}_{\mathbf{s}}\right[k\left]\right|{\mathbf{y}}_k^i)$

∙ From TOA measurements and Proposition 1, evaluate the weight of each particle. For $i=1,\dots ,N$

${\tilde{\omega }}_{\tau }^i=p\left({\mathbf{z}}_{\mathbf{\tau }}\right[k\left]\right|{\mathbf{y}}_k^i)$

∙ Evaluate for $i=1,\dots ,N$

${\tilde{\omega }}_k^i={\mathrm{\omega }}_{k-1}^i\frac{{\tilde{\omega }}_s^i{\tilde{\omega }}_{\tau }^{i}p\left({\mathbf{y}}_k^i\right|{\mathbf{y}}_{k-1}^i)}{q\left(\mathbf{y}\right[k\left]\right|{\mathbf{y}}_{k-1}^i,\mathbf{z}\left[k\right])}$

∙ Normalization: for$i=1,\dots ,N$, compute

${\omega }_k^i=\frac{{\tilde{\omega }}_k^i}{\sum _{j=1}^N{\tilde{\omega }}_k^j}$

Table 1.

ALPA filtering.

To implement the algorithm detailed in Table 1, we have to choose a proposal distribution, where the most popular choice is to use the transition prior given by the dynamic model, i.e., $p\left(\mathbf{y}\right[k\left]\right|\mathbf{y}[k-1])$ [19]. This election leads to a rather simple expression for the weights,

$${\tilde{\omega }}_k^i={\mathrm{\omega }}_{k-1}^i{\tilde{\omega }}_s^i{\tilde{\omega }}_{\tau }^i$$E19

Therefore, in order to use this algorithm, we have to obtain samples from the transition prior and evaluate the adaptive likelihood function pointwise. Figure 6 summarizes how this filter works with the proposal distribution chosen. First, we generate particles from the proposal distribution, in this case, the prior distribution, $p\left(\mathbf{y}\left[k\right]|\mathbf{y}\left[k-1\right]\right)$, and then, their weights are updated according to the likelihood function, $p\left(\mathbf{z}\left[k\right]|\mathbf{y}\left[k\right]\right)$. If the support of the proposal distribution does not cover the support of the likelihood function, only few particles will be in the region of importance, thus, the number of particles has to be increased in order to correctly approximate the posterior distribution.

Figure 6.

Transition prior and likelihood functions. Particles are obtained by sampling from the prior and weighting from the likelihood.

In this SIS algorithm, as $k$ increases, the variance of the weights ${\omega }_k^i$ also increases, and therefore, after a certain number of steps, all but one particle will have negligible normalized weights. This problem is known as degeneracy [19]. To overcome this drawback, it is mandatory to perform a resampling step when a severe degeneracy is detected. A measure of degeneracy is the effective sample size $N_{eff}$, estimated as,

$${\widehat{N}}_{eff}=\frac{1}{\sum _{i=1}^N{{(\omega }_k^i)}^2}$$E20

where a small ${\widehat{N}}_{eff}$ indicates a severe degenerancy. Therefore, when degenerancy is detected, $N$ samples with uniform weights are drawn from the discrete representation of the posterior, given by the previous particles and weights, yielding a variant of SIS algorithm called sampling-importance-resampling (SIR) algorithm [19,52].

5. Results

The goal of this section is to quantify the performance of the methods presented in the above sections, leading to the ALPA filter. In order to do that, we obtained experimental data in a real indoor scenario by using the systems described in [10] and [12], and we ran numerous Monte Carlo simulations. In the following, we compare the performance of the introduced techniques with conventional approaches as well as with the CRLB.

We use the dynamic and measurements models above described together with the following state vector and prior information, depending on whether we estimate ranges or positions,

• Range estimation: we use a state vector $\mathbf{y}\left[k\right]=(d\left[k\right],d^{\text{'}}\left[k\right],d^{\text{'}\text{'}}[k\left]\right)$. The standard deviation ${\sigma }_{d^{\left(3\right)}}$ is $1$ m/s³, which is roughly 50% of the maximum [18]. Furthermore, we add prior information about first and second derivatives of the distance, by considering they are distributed as Gaussians $\mathcal{N}(0,{\sigma }_{d^{\text{'}}})$ and $\mathcal{N}(0,{\sigma }_{d^{\text{'}\text{'}}})$ respectively, where ${\sigma }_{d^{\text{'}}}=0.5$ m/s and ${\sigma }_{d^{\text{'}\text{'}}}=0.5$ m/s².

• Position estimation: we use a state vector $\mathbf{y}\left[k\right]=(\mathbf{x}\left[k\right],\mathbf{v}\left[k\right],\mathbf{a}[k\left]\right)$, where $\mathbf{x}\left[k\right]$ consists of the two-dimensional coordinates of the mobile target’s position, and $\mathbf{v}\left[k\right]$ and $\mathbf{a}\left[k\right]$ are the velocity and the acceleration vectors. The same previous values for the deviations of the derivatives of the coordinates are used for dynamic and prior information.

For the experimental data, the target carried a laptop equipped with an IEEE 802.11b/g adapter and the measuring systems described in [10] and [12]. The anchors consisted of IEEE 802.11b/g access points (APs). In the RSS case, the anchors periodically sent beacon frames (at a frequency of $\mathrm{M}$Hz) and the RSS values were obtained based on the RSS indicator at target’s adapter [12]. In the TOA case, the mobile target periodically sent request-to-send frames to each anchor (at a frequency of $\mathrm{M}$Hz), and a counter connected to the WLAN adapter saved the clock-cycles elapsed between the request and the reception of the corresponding clear-to-send frame [10]. For the results presented in this section, we refer as fusion the results of combining RSS and TOA data at every time-step.

5.1. Experimental results

As mentioned above, in a realistic scenario, NLOS propagation together with multipath effects constitute the major drawback of localization in harsh environments. This section illustrates the behavior of the proposed algorithm during a typical path followed by a mobile target in an indoor scenario. We carried out a measurement campaign inside an office building cluttered with clusters of objects and people moving freely in the area of the measurements. The propagation conditions were even harsher than the ones commonly find by an OBU placed within a car. Figure 8 shows the trajectory of 65 meters as well as the position of the 4 APs. It took 100 seconds to complete the whole trajectory, receiving a new set of measurements every second ($∆t=1$s) from all the APs. As reflected in Figure 8, NLOS was always present when measuring with respect to AP3 and AP4, and only in a small percentage of positions there was a LOS between target and anchors AP1 and AP2.

In Table 2, we compare the error achieved with the proposed ALPA range estimation method in the presented scenario to the error obtained with conventional approaches [15,24]. We specify the results for RSS-only and TOA-only cases, and for their fusion. Specifically, we call,

• ML-RSS, ML-TOA, ML-Fusion: the range estimates obtained by means of the ML estimator. We utilize as likelihood function the convolution of the likelihood reported by the measurements (log-normal in the RSS case and Gaussian in the TOA case) and a Gaussian distribution corresponding to the bias.

  In order to guarantee a fair comparison, in Table 2 and in the following experiments, we select the values for the biases in accordance to the ones selected in Section 4. In this way, the RSS bias is modeled as a Gaussian $\mathcal{N}(0,{\sigma }_s)$ with ${\sigma }_s=3$ dBm, and the TOA bias as a Gaussian $\mathcal{N}({\gamma }_{\tau }/2,{\gamma }_{\tau }/4)$ with ${\gamma }_{\tau }=4$ clock cycles.

  The likelihood for the fusion is computed from (18).

• AML-RSS, AML-TOA, AML-Fusion: the ranges that correspond to the result of obtaining the maximum of the adaptive likelihood computed by means of Proposition 1, and (18) in the fusion case.

• EKF-RSS, KF-TOA, EKF-Fusion: the result of applying EKF and KF filters for RSS and TOA measurements, respectively, using the same bias distributions as in the ML case, and the dynamic model given by (3).

• ALPA-RSS, ALPA-TOA, ALPA-Fusion: the range estimates obtained by the ALPA filtering described in Table 1, where $N=10 000$ is the number of particles used.

We summarize for all these methods the quartiles of the absolute error in range estimates as well as the root mean squared error (RMSE), which incorporates both systematic (bias) and random errors. In order to study the influence of the number of measurements, $M$, in the final performance, all these statistics are shown for four different values.

Figure 7 depicts the pdf of the absolute error in range estimation after applying AML-Fusion and ALPA-Fusion methods, taking $10$ RSS and 10 TOA measurements in each one of the positions of the target with respect to the four APs. Figure 7 likewise includes the ML-Fusion and EKF-Fusion methods in order to compare their behavior. Using only $10$ measurements, ML-, AML-, EKF- and ALPA-Fusion obtain an error in range estimation lower than $3$ meters for 55%, 65%, 73%, and 80% of the positions, respectively, which reflects the remarkable performance of the proposed algorithm.

	$M=5$		$M=10$		$M=50$		$M=100$
Quartiles	RMSE	Quartiles	RMSE	Quartiles	RMSE	Quartiles	RMSE
ML-RSS	1.64-3.12-5.45	7.01	1.28-2.94-4.96	5.32	1.36-2.72-4.68	4.34	1.27-2.74-4.74	4.58
ML-TOA	2.09-3.92-7.68	6.42	1.55-3.40-5.64	5.00	1.26-2.69-4.40	3.87	1.12-2.44-4.00	3.55
ML-Fusion	1.52-3.16-5.87	5.25	1.26-2.66-4.73	4.23	1.09-2.24-3.89	3.55	0.87-2.18-3.61	3.26
AML-RSS	1.69-3.25-5.27	5.64	1.44-2.92-5.06	4.71	1.32-2.74-4.64	4.27	1.31-2.70-4.50	4.20
AML-TOA	2.06-3.74-7.38	6.28	1.52-3.31-5.57	4.93	1.18-2.61-4.27	3.81	1.03-2.38-3.86	3.48
AML-Fusion	1.38-2.91-5.19	4.49	1.15-2.32-3.65	3.49	0.86-1.91-3.39	3.06	0.83-1.83-3.26	2.91
EKF-RSS	0.84-2.22-4.26	3.82	1.06-2.59-4.21	3.81	1.21-2.43-4.07	3.76	1.17-2.55-4.04	3.69
KF-TOA	1.11-2.37-3.95	3.60	1.10-2.06-3.63	3.04	0.81-1.76-2.97	2.53	0.86-1.63-2.95	2.36
EKF-Fusion	0.93-1.90-3.24	2.78	0.86-1.82-3.15	2.59	0.82-1.62-2.62	2.25	0.74-1.49-2.55	2.10
ALPA-RSS	0.82-2.33-4.63	3.88	1.17-2.58-4.30	3.79	1.20-2.48-4.18	3.75	1.21-2.64-4.17	3.78
ALPA-TOA	0.94-2.04-3.33	3.11	0.95-1.90-3.06	2.69	0.72-1.48-2.63	2.52	0.76-1.50-2.64	2.32
ALPA-Fusion	0.84-1.72-2.95	2.58	0.80-1.70-2.85	2.35	0.69-1.37-2.36	2.22	0.70-1.45-2.40	2.08

Table 2.

Range estimation error quartiles and RMSE obtained with different algorithms as a function of the number of measurements. All error values are in meters.

Analogously, in Figures 8-9 and Table 3, we summarize the results in position estimation. In this case, we call,

For the results of Figures 8-9 and Table 3, EKF and ALPA filters use a measurements model that directly relates measurements with positions, avoiding the intermediate step of estimating distances and, therefore, removing the trilateration stage.

• ML-RSS, ML-TOA, ML-Fusion: the positions obtained with the ML distances and a trilateration technique based on the radical axis of the circles drawn at each anchor’s position [10,12-13].

• EKF-RSS, EKF-TOA, EKF-Fusion: the positions obtained by means of an EKF whose measurements model relates the measurements to the target’s position.

• PF-RSS, PF-TOA, PF-Fusion: the result of applying the ALPA filter described in Table 1 to the positional-states, with $N=10 000$ particles.

Figure 7.

The height and width of the pdf corresponding to the error achieved by the ALPA filter reflect its better performance in comparison to other conventional range estimation techniques. $10$ RSS and $10$ TOA measurements were taken with respect to each anchor.

	$M=5$		$M=10$		$M=50$		$M=100$
Quartiles	RMSE	Quartiles	RMSE	Quartiles	RMSE	Quartiles	RMSE
ML-RSS	3.83-5.91-8.49	12.99	3.32-5.21-7.49	8.91	3.35-4.94-6.98	6.64	3.26-5.00-6.60	7.43
ML-TOA	3.95-6.14-8.03	7.64	2.80-4.05-6.64	5.70	2.24-3.34-5.16	4.57	1.63-3.20-4.71	4.09
ML-Fusion	3.15-4.95-7.04	6.73	2.40-3.71-6.11	5.10	1.93-3.03-4.93	4.34	1.64-3.03-4.53	3.89
EKF-RSS	2.94-4.46-6.18	5.11	3.47-4.83-6.85	5.83	3.02-4.12-6.33	5.24	3.02-4.24-6.40	5.24
KF-TOA	1.77-2.79-4.25	3.54	2.05-2.80-3.64	3.11	1.54-2.28-3.09	2.61	1.50-2.22-3.14	2.51
EKF-Fusion	2.20-3.24-4.30	3.50	2.08-2.99-3.90	3.25	1.76-2.32-3.00	2.57	1.71-2.13-2.98	2.41
ALPA-RSS	1.93-3.28-5.18	4.36	3.16-3.91-5.18	4.68	2.38-3.09-4.61	4.23	2.72-3.65-4.97	4.37
ALPA-TOA	1.90-2.59-3.76	3.37	1.63-2.54-3.63	2.98	1.08-1.98-3.25	2.66	1.35-2.18-3.05	2.63
ALPA-Fusion	1.77-2.86-3.46	3.14	1.92-2.61-3.34	2.82	1.23-1.85-3.15	2.49	1.28-2.00-2.64	2.40

Table 3.

Position estimation error quartiles and RMSE obtained with several algorithms as a function of the number of measurements. All error values are in meters.

Figure 8.

Trajectory followed by the target and position estimates for different positioning methods. $10$ RSS and $10$ TOA measurements were taken with respect to each anchor.

Figure 9 depicts the pdf of the error in position estimation for the three mentioned RSS/TOA fusion algorithms, taking $10$ RSS and 10 TOA measurements in each one of the positions of the target with respect to the four APs. Using only $10$ measurements, ML-, EKF- and ALPA-Fusion obtain an error in position estimation lower than $3$ meters for 40%, 52%, and 63% of the positions, respectively.

Figure 9.

The proposed ALPA filter obtains the best performance with an error lower than 3 meters for more than 63% of the positions.

Figures 8-9 and Table 3 show the better performance of the proposed ALPA filter for all the analyzed scenarios, resulting, for example, in an RMSE of 2.82 meters for the case of only using $10$ RSS and $10$ TOA measurements, while previous essays obtained RMSEs around 4 meters by using hundreds of measurements [12,28].

5.2. Simulation results

The CRLB provides a lower bound on the minimum achievable mean squared estimation error for any unbiased estimator. In what follows, we use such metric to assess the optimality of the presented ALPA filter against such lower bound.

The Bayesian version of the CRLB is known as the Van Tress CRLB [53], or posterior CRLB, since it is obtained from the posterior distributions of the random state vector [54]. In our case, for each time instant $k$, the CRLB is,

$$\mathrm{{\rm E}}\left\{\right(g\left(\mathbf{Z}\left[k\right]\right)-\mathbf{y}\left[k\right]\left){\left(g\right(\mathbf{Z}\left[k\right])-\mathbf{y}[k\left]\right)}^T\right\}\succcurlyeq {\mathbf{J}}_k^{-1}$$E21

where $g(\mathbf{Z}\left[k\right]$) is an unbiased estimator of $\mathbf{y}\left[k\right]$ and ${\mathbf{J}}_k$ is the Fisher information matrix (FIM) obtained as,

$${\mathbf{J}}_k=-\mathrm{{\rm E}}\left\{{\nabla }_{\mathbf{y}\left[k\right]}{[{\nabla }_{\mathbf{y}\left[k\right]}\log p\left(\mathbf{Z}\left[k\right]|\mathbf{y}\left[k\right]\right)}^{\mathrm{T}}\right]\}$$E22

Tichavský et al. proposed a recursive formula to compute the FIM [55]. For the particular case of the linear-Gaussian dynamic model in (3), being ${\mathbf{Q}}_k$ the covariance matrix in this model, the FIM is given by the recursion [19],

$${\mathbf{J}}_{k+1}={\mathbf{J}}_{k+1}^{\mathbf{z}}+({{\mathbf{Q}}_k+{\mathbf{F}}_k{\mathbf{J}}_k^{-1}{\mathbf{F}}_k^T)}^{-1}$$E23

and

$${\mathbf{J}}_{k+1}^{\mathbf{z}}=-\mathrm{{\rm E}}\left\{{\nabla }_{\mathbf{y}\left[k+1\right]}{[{\nabla }_{\mathbf{y}\left[k+1\right]}\log p\left(\mathbf{z}\left[k+1\right]|\mathbf{y}\left[k+1\right]\right)}^{\mathrm{T}}\right]$$E24

To start this recursion, we assume the initial density as Gaussian, then, the initial FIM coincides with its covariance matrix.

Figure 10 compares the RMSE obtained in range estimation by means of the proposed ALPA-Fusion filter with the RMSE obtained by applying the EKF-Fusion method, and with the square root of the CRLB.

We selected a truncated normal distribution as random error to reflect the limited range of the measuring systems. For the proposed adaptive likelihoods, ${\mathbf{J}}_{k+1}^{\mathbf{z}}$ has no closed-form, then, it was evaluated by Monte Carlo integration.

To obtain such curves, we simulated a trajectory of 85 positions and carried out $1 000$ Monte Carlo experiments. Figure 10 again corroborates the remarkable performance of ALPA filter, since the corresponding curve is much closer to the CRLB than the line corresponding to the EKF error.

Figure 10.

The near-optimal performance of the proposed ALPA filter in harsh environments is corroborated by comparison with the CRLB.

6. Complexity

The key issue in PFs is the exponential growth of computational complexity as a function of the dimension of the state vector, $\mathbf{y}\left[k\right]$, whereas EKF grows as the cube of the dimension [56]. For low dimensional problems, PF remains similar to an EKF, however, for high dimensional problems, PFs suffer from the curse of dimensionality [57]. Then, PFs that track ranges instead of positions can be advantageous from a complexity point of view.

Moreover, from Proposition 1, the complexity of the likelihood grows exponentially with the number of samples. However, this complexity can be reduced by removing redundant components from the RSS and TOA pdfs or from the resulting fusion mixture. To this aim, different criteria such as William‘s criterion [58], Kullback-Leibler distance [59] or clustering [60] can be utilized. Therefore, considering the improvement achieved in range and position estimation with $5$ and $10$ measurements, the proposed ALPA filter could be a good choice for the designing of VANETs that require low consumption. In these cases, in order to save battery, the OBUs transmit only at discrete intervals; therefore, there is more time available for processing a smaller number of samples.

7. Conclusions

In this chapter we have presented an adaptive likelihood function for robust data fusion in localization systems. Based on this likelihood, we have developed the ALPA filter for range and position estimation. This ALPA filter presents several advantages over conventional techniques,

1. it does not assume any parametric statistical model, utilizing the empirical distribution of the measurements at each time by means of Gaussian kernels;

2. it adaptively fuses RSS and TOA data being extensible to any other type of measurement;

3. it takes advantage of the relationship among positions in time by using Bayesian filtering;

4. it addresses the non-linear and non-Gaussian behavior of the measurements by using particle filtering.

These advantages result in a noticeable improvement with respect to other conventional techniques, as corroborated by the experimental and simulation results. Under NLOS and multipath conditions, ALPA filter obtains not only an RMSE in position estimation lower than 3 meters with only 10 RSS and 10 TOA measurements, but also an error remarkably close to the theoretical benchmark provided by the CRLB.

Therefore, ALPA filter is a valuable choice to provide localization in V2I communication systems. Its extension to cooperative localization would make this localization also possible in VANETs based on V2V communication.