Robust Interval-Based Localization Algorithms for Mobile Sensor Networks

Abstract

This paper considers the localization problem in mobile sensor networks. Such a problem is a challenging task, especially when measurements exchanged between sensors may contain outliers, that is, data not matching the observation model. This paper proposes two algorithms robust to outliers. These algorithms perform a set-membership estimation, where only the maximal number of outliers is required to be known. Using these algorithms, estimates consist of sets of boxes whose union surely contains the correct location of the sensor, provided that the considered hypotheses are satisfied. This paper proposes as well a technique for evaluating the number of outliers to be robust to. In order to corroborate the efficiency of both algorithms, a comparison of their performances is performed in simulations using Matlab.

1. Introduction

Mobile Sensor Networks (MSNs) have recently emerged as a challenging research field. An MSN consists of a large number of low-cost smart sensors with limited computational capacities and energy resources [1]. Due to the lack of a fixed infrastructure in MSNs, the sensors are able to move in an uncontrolled manner. For this reason and since sensed data are related to the locations of the sensors in almost all MSN applications, many researchers have focused on the localization problem. A first solution for sensor localization is to equip all sensors with global positioning systems (GPSs) [2]. However, this solution is nonpractical in MSNs, since GPS are expensive, high energy consuming, and having great sizes. The alternative solution consists of equipping a few number of sensors with GPS receivers. These sensors, aware of their locations, are called anchors. The remaining sensors, called nonanchor nodes or simply nodes, have unknown locations, and hence they need to be localized.

Many anchor-based algorithms have been proposed for sensor localization. For instance, Doherty et al. [3] propose a centralized technique for position estimation. Localization in [3] is formalized as a convex optimization problem, having connectivity measurements between sensors as constraints. In [4–6], different approaches requiring few anchors have been proposed. Local maps with relative positions are constructed using measured distances between nodes and their neighbors. A combination of these maps with known positions of anchors lead to absolute positions. Nevertheless, these techniques are not very robust because of errors accumulation while combining the maps. Authors in [7] propose a distributed static algorithm, where each node defines its position as the center of all observed anchors. In a different scenario [8], Galstyan et al. propose an online distributed technique, where nodes use their detection of a moving target to update their position estimates. Blatt and Hero [9] address the problem of source localization using sensors measurements. The problem is formulated as a convex problem that is solved using the aggregated projection onto convex sets (APOCSs) method. In [10–12], dynamic approaches, on the basis of sequential Monte-Carlo [13], are considered to estimate the positions of the nodes. The posterior distribution of the unknown positions is estimated recursively with a set of particles. Terwilliger et al. [14] propose to cover all possible solutions with the smallest enclosing disk. Alternative dynamic algorithms for sensor localization, using interval analysis [15], have been proposed in [16, 17]. Position estimates are boxes covering the possible locations of the sensors.

Existing methods have mainly considered the localization problem with the hypothesis that all measurements are consistent with the considered measurement model. However, in practical situations, outliers, that is, data not matching the measurement model, are encountered. Previously mentioned estimation techniques are not very robust to such outliers. In [18], Jaulin et al. propose a set-membership estimator robust to outliers based on interval analysis. Savarese et al. [19] present a distributed robust algorithm for sensor localization. The method is separated into two phases: the start-up phase, where first estimates of node positions are computed using hop counts to anchors, and the refinement phase, where nodes communicate with their neighbors to update their positions using a least-squares triangulation technique. Nevertheless, a number of factors influence the convergence of the refinement phase, such as the accuracy of first estimates and the magnitude of ranging errors. In [20], Rabbat et al. introduce a robust localization algorithm of an isotropic energy source using kernel averaging techniques. The proposed estimator is more robust than the least-squares estimator under a variety of conditions. Leger and Kieffer [21] present a distributed version of the estimation algorithm [18], assuming that the maximal number of outliers is known. In particular, a static distributed algorithm is proposed for source localization using received signal strength (RSS) measurements. The proposed method adapts the set inversion via interval analysis (SIVIA) algorithm [15] to evaluate a solution set.

In this paper, we propose an original adaptive approach for sensor localization in the presence of outliers. Assuming only that the maximal number of outliers is given, the proposed approach uses connectivity measurements in addition to a mobility model to address the localization problem. The solution is then given using either SIVIA or an alternative combinatorial technique. Another contribution of the paper is that it proposes a technique for evaluating the maximal number of outliers to be robust to. Moreover, using a connectivity-based observation model, the paper compares the performances of both robust localization algorithms.

The rest of the paper is organized as follows. Section 2 introduces the localization problem. A description of the SIVIA algorithm and the combinatorial technique is then given in Section 3. Section 4 provides simulation results, whereas Section 5 concludes the paper.

2. Problem Statement

The proposed method is an anchor-based method, where each node exchanges information with anchors to localize itself. Consider a network of $N_a$ anchors and $N_x$ mobile nodes. All sensors are assumed to be in the same plane: their locations at time t are given by ${\mathbf{a}}_i(t)=(a_{i,1}(t),a_{i,2}(t){)}^T$ , $i=1,\dots ,N_a$ , for anchors and ${\mathbf{x}}_j(t)=(x_{j,1}(t),x_{j,2}(t){)}^T$ , $j=1,\dots ,N_x$ , for nodes. In order to reduce the communication costs during the localization process, the proposed method assumes that each mobile node does not exchange information with other nodes. For this reason and without loss of generality, we focus on the localization of one generic mobile node $\mathbf{x}(t)=(x_1(t),x_2(t){)}^T$ , and we thus drop the index j.

2.1. Observation Model

At time t, the mobile node receives signals from a set $J(t)\subseteq \{1,\dots ,N_a\}$ of anchors with specific received signal strengths (RSSs) denoted by ${\rho }_i(t)$ , $i\in J(t)$ . These RSSs are assumed to follow the Okumura-Hata model [22]

${\rho }_i\left(t\right)={\rho }_0-10n_P{\log }_{10}\frac{d_i\left(\mathbf{x}\left(t\right)\right)}{d_0}+{\varepsilon }_i\left(t\right).$

(1)

In (1), ${\rho }_i(t)$ is in $\text{dBm}$ , ${\rho }_0$ is the power measured (in $\text{dBm}$ ) at a reference distance $d_0$ from the anchor ${\mathbf{a}}_i(t)$ , $n_P$ is the path-loss exponent, $d_i(\mathbf{x}(t))=\parallel \mathbf{x}(t){-\mathbf{a}}_i(t)\parallel$ is the Euclidian distance between the anchor ${\mathbf{a}}_i(t)$ and the considered node, and ${\varepsilon }_i(t)$ is the measurement noise, modeled as zero-mean Gaussian with variance ${\sigma }^2$ .

In practice, ${\rho }_0$ and $n_P$ may vary from one anchor to the other, and ${\sigma }^2$ may be quite large. Given the RSS values, the proposed model may lead to inaccurate distance estimates. For this reason, only connectivity information are employed and (1) is only used to determine whether the node is in the vicinity of the ith anchor. Let ${\rho }_r$ be some RSS threshold corresponding to a distance r, which is the sensing range of the sensors. Then, if ${\rho }_i(t)\ge {\rho }_r$ , the distance $d_i(\mathbf{x}(t))$ from the anchor i to the node is deemed less than r. Anchors for which ${\rho }_i(t)\ge {\rho }_r$ are called detected anchors. Only detected anchors are then taken into account for the localization. The observation model is then given by

${\left(x_1\left(t\right)-a_{i,1}\left(t\right)\right)}^2+{\left(x_2\left(t\right)-a_{i,2}\left(t\right)\right)}^2\le r^2,\hspace{1em}i\in I\left(t\right),$

(2)

where $I(t)\subseteq J(t)$ is the set of indices of detected anchors, that is, whose emitted signals have RSS at the node $\mathbf{x}(t)$ larger than ${\rho }_r$ . The observation model is thus given by a set of disk equations centered on the detected anchors and having r as radius.

In real environments, measurements may not follow exactly the observation model. Indeed, due to the additive noise and the inaccuracy of the parameter values, a measured RSS ${\rho }_i(t)$ could be less than ${\rho }_r$ while (2) is satisfied for real and vice versa. In the first case, the anchor is assumed to be out of the vicinity of the node, which is not true, and thus, a correct measurement is omitted, whereas in the second case, an outlier is obtained. The proposed approach takes such outliers into consideration. Using the connectivity-based model, it assumes that the maximal number of outliers is known and denoted by q. In other words, it considers that $|I(t)|-q$ measurements at minimum are correct at each time step.

2.2. Mobility Model

The proposed method takes also advantage of the mobility of the nodes to improve the estimation accuracy. Any available information about the motion of the node could be used to define the mobility model. This paper proposes a very general mobility model, where only the maximal velocity of the node $v_{\max \hspace{0.17em}}$ is assumed to be known. Then, the positions of the generic node at time steps $t-\mathrm{\Delta }t$ and t satisfy

${\left(x_1\left(t\right)-x_1\left(t-\Delta t\right)\right)}^2+{\left(x_2\left(t\right)-x_2\left(t-\Delta t\right)\right)}^2\le \Delta t^2\cdot v_{\max }^2.$

(3)

More generally, the mobility model could be reformulated as follows:

$\mathbf{f}\left(\mathbf{x}\left(t-\Delta t\right),\mathbf{x}\left(t\right),v\right)=0,$

(4)

where v is some parameter only known to belong to some known interval $[v]$ (here, $[v]=[0,v_{\max \hspace{0.17em}}]$ ).

2.3. Description of the Robust Set-Membership Localization

Estimating the location of the sensor at time t consists of finding the set $𝕏(t)$ of all locations consistent with the mobility model (4) and at least $|I(t)|-q$ observation constraints (2). In other words, these locations should be in the vicinity of at least $|I(t)|-q$ detected anchors. A set-membership estimator [23] robust to q outliers [24] at time t is then obtained, since any $|I(t)|-q$ measurements instead of $|I(t)|$ measurements are considered for the estimator.

Assume that $\mathbf{x}(t-\mathrm{\Delta }t)$ belongs to some set $𝕏(t-\mathrm{\Delta }t)$ . According to this approach, and to compute $𝕏(t)$ , a predicted set is first evaluated using the mobility model

${\mathbb{X}}^{*}\left(t\right)=\left\{\mathbf{x}\hspace{0.17em}|\hspace{0.17em}\exists {\mathbf{x}}^{\prime }\in \mathbb{X}\left(t-\Delta t\right),\exists v\in \left[v\right],\mathbf{f}\left({\mathbf{x}}^{\prime },\mathbf{x},v\right)\le 0\right\}\hspace{0.17em}\hspace{0.17em}.$

(5)

Measurements are then taken into account to correct ${𝕏}^{*}(t)$ as follows:

$\mathbb{X}\left(t\right)={\displaystyle \bigcup }_{C{\in \mathbf{C}}_{I\left(t\right)}^{\left|I\left(t\right)\right|-q}}\left(\begin{array}{c}{\displaystyle \bigcap }_{\mathcal{l}\in C}{\mathbb{X}}_{\mathcal{l}}\left(t\right)\end{array}\right),$

(6)

where ${\mathbf{C}}_{I(t)}^{|I(t)|-q}$ is the set of all $(|I(t)|-q)$ -combinations of indices in $I(t)$ , C is a set of indices belonging to ${\mathbf{C}}_{I(t)}^{|I(t)|-q}$ , and

${\mathbb{X}}_{\mathcal{l}}\left(t\right)=\left\{\mathbf{x}\in {\mathbb{X}}^{*}\left(t\right)\hspace{0.17em}|\hspace{0.17em}\Vert \mathbf{x}{-\mathbf{a}}_{\mathcal{l}}\left(t\right)\Vert \le r\right\}.$

(7)

Then, ${\bigcap }_{\ell \in C}{𝕏}_{\ell }(t)$ denotes the set of locations of ${𝕏}^{*}(t)$ that satisfy the specific observation constraints denoted in C, whereas $𝕏(t)$ denotes the set of locations of ${𝕏}^{*}(t)$ that satisfy any $(|I(t)|-q)$ observation constraints. The number of combinations to be considered in (6) is $K(t)=|I(t)|!/(|I(t)|-q)!q!$ , where $q!$ denotes the factorial of q.

An alternative definition of $𝕏(t)$ inspired by [24] is

$\mathbb{X}\left(t\right)=\left\{\mathbf{x}\in {\mathbb{X}}^{*}\left(t\right)\hspace{0.17em}|\hspace{0.17em}{\displaystyle \sum }_{i\in I\left(t\right)}^{}\lambda \left(\mathbf{x},{\mathbf{a}}_i\left(t\right)\right)\ge \left|I\left(t\right)\right|-q\right\},$

(8)

where

$\lambda \left(\mathbf{x},{\mathbf{a}}_i\left(t\right)\right)=\{\begin{array}{cc}1,\hfill & \text{if}\hspace{0.17em}\hspace{0.17em}\Vert \mathbf{x}-{\mathbf{a}}_i\left(t\right)\Vert \le r,\hfill \\ 0,\hfill & \text{otherwise}.\hfill \end{array}$

(9)

This definition does not involve any combinatorial. One may easily prove that (6) and (8) are equivalent. This technique would be used in the following to solve the localization problem.

3. Localization Algorithms

Solving the localization problem in a guaranteed way consists of finding the set of all node locations that satisfy the problem constraints while being robust to outliers. In this paper, interval analysis [15] is employed to achieve this goal. At each time step, the proposed method computes a set of nonoverlapping boxes, called subpaving [15], whose union covers the solution set $𝕏(t)$ . Assume that $[𝕏](t)$ is the solution subpaving containing the actual position of the generic node at time t. One has

$\left[\mathbb{X}\right]\left(t\right)={\displaystyle \bigcup }_{1\le j\le n\left(t\right)}\left[{\mathbf{x}}_j\right]\left(t\right),$

(10)

where $[{\mathbf{x}}_j](t)=[x_{j,1}](t)\times [x_{j,2}](t)$ is a two-dimensional box and $n(t)$ is the number of boxes in $[𝕏](t)$ . As shown in Section 2.3, finding $[𝕏](t)$ involves a prediction phase followed by a correction phase.

3.1. Prediction Phase

Assume that $[𝕏](t-\mathrm{\Delta }t)$ is the subpaving obtained at time $t-\mathrm{\Delta }t$ . Computing the predicted set $[𝕏{]}^{*}(t)$ may be done by evaluating (5), where $𝕏(t-\mathrm{\Delta }t)$ is replaced by $[𝕏](t-\mathrm{\Delta }t)$ , which is quite difficult. Nevertheless, for each box ${[\mathbf{x}}_j](t-\mathrm{\Delta }t)=[x_{j,1}](t-\mathrm{\Delta }t)\times [x_{j,2}](t-\mathrm{\Delta }t)\in [𝕏](t-\mathrm{\Delta }t)$ , the corresponding box ${[\mathbf{x}}_j{]}^{*}(t)=[x_{j,1}{]}^{*}(t)\times [x_{j,2}{]}^{*}(t)\in [𝕏{]}^{*}(t)$ has to be compliant with the mobility model (3), leading to the following constraint:

$\begin{array}{cc}\hfill & {\left({\left[x_{j,1}\right]}^{*}\left(t\right)-\left[x_{j,1}\right]\left(t-\Delta t\right)\right)}^2\hfill \\ \hfill & \hspace{1em}\hspace{1em}+{\left({\left[x_{j,2}\right]}^{*}\left(t\right)-\left[x_{j,2}\right]\left(t-\Delta t\right)\right)}^2\subseteq \left[0,\Delta t^2\cdot v_{\max }^2\right].\hfill \end{array}$

(11)

Relaxing (11) yields the following expressions of $[x_{j,1}{]}^{*}(t)$ and $[x_{j,2}{]}^{*}(t)$ :

$\begin{array}{cc}\hfill {\left[x_{j,i}\right]}^{*}\left(t\right)& =\left[x_{j,i}\right]\left(t-\Delta t\right)\hfill \\ \hfill & \hspace{1em}+\left[-\Delta t\cdot v_{\max },\Delta t\cdot v_{\max }\right],\hspace{1em}i=\mathrm{1,2}.\hfill \end{array}$

(12)

Let $[\widetilde{𝕏}{]}^{*}(t)={\bigcup }_j[{\mathbf{x}}_j{]}^{*}(t)$ be the set of all boxes evaluated with (12). One may prove that $[𝕏{]}^{*}(t)\subset [\widetilde{𝕏}{]}^{*}(t)$ . The convex hull $[\mathbf{x}{]}^{*}(t)=[x_1{]}^{*}(t)\times [x_2{]}^{*}(t)$ [15] of $[\widetilde{𝕏}{]}^{*}(t)$ is the smallest box containing $[\widetilde{𝕏}{]}^{*}(t)$ . Its components are defined as

${\left[x_i\right]}^{*}\left(t\right)=\left[\underset{j\le n\left(t-1\right)}{\min }\left({\underline{\mathbf{x}}}_{j,i}^{*}\left(t\right)\right),\underset{j\le n\left(t-1\right)}{\max }\left({\bar{\mathbf{x}}}_{j,i}^{*}\left(t\right)\right)\right],\hspace{1em}i=\mathrm{1,2},$

(13)

where ${\underline{\mathbf{x}}}_{j,1}^{*}(t)$ and ${\bar{\mathbf{x}}}_{j,1}^{*}(t)$ are the low and high endpoints of $[{\mathbf{x}}_{j,1}{]}^{*}(t)$ , respectively. ${[\mathbf{x}]}^{*}(t)$ is a rectangular area that covers all possible locations that could be taken by the node at time t according to its mobility model. The convex hull is used in the correction phase instead of $[\widetilde{𝕏}{]}^{*}(t)$ to reduce the computational complexity.

3.2. Correction Using the SIVIA Algorithm

The SIVIA algorithm [15] performs a succession of bisections and selections of boxes compliant with the localization constraints. Let $[\mathbf{x}]$ be a box, set initially to $[\mathbf{x}{]}^{*}(t)$ . The following test function is evaluated:

$\gamma \left(\left[\mathbf{x}\right]\right)=\{\begin{array}{cc}1\hfill & \text{if}\hspace{0.17em}\hspace{0.17em}{\displaystyle \sum }_{i\in I\left(t\right)}\lambda (\left[\mathbf{x}\right],{\mathbf{a}}_i)\ge \left|I\left(t\right)\right|-q,\hfill \\ -1\hfill & \text{if}\hspace{0.17em}\hspace{0.17em}{\displaystyle \sum }_{i\in I\left(t\right)}\tilde{\lambda }(\left[\mathbf{x}\right],{\mathbf{a}}_i)>q,\hfill \\ 0\hfill & \text{otherwise},\hfill \end{array}$

(14)

where $\lambda ([\mathbf{x}],\mathbf{a})$ is equal to 1 if $\sup \hspace{0.17em}(\parallel [\mathbf{x}]-\mathbf{a}\parallel )\le r$ and 0 otherwise, while $\tilde{\lambda }([\mathbf{x}],\mathbf{a})$ is equal to 1 if $\inf \hspace{0.17em}(\parallel [\mathbf{x}]-\mathbf{a}\parallel )>r$ and 0 otherwise. Graphically, $\lambda ([\mathbf{x}],\mathbf{a})=1$ means that the box $[\mathbf{x}]$ is entirely included in the connectivity disk centered on the anchor a, whereas $\tilde{\lambda }([\mathbf{x}],\mathbf{a})=1$ means that the box $[\mathbf{x}]$ is entirely outside the connectivity disk centered on the anchor a.

The box $[\mathbf{x}]$ is added to the solution $[𝕏](t)$ if $\gamma ([\mathbf{x}])=1$ , meaning that all $\mathbf{x}\in [\mathbf{x}]$ satisfy at least $|I(t)|-q$ observation constraints. Boxes inconsistent with more than q observation constraints ( $\gamma ([\mathbf{x}])=-1$ ) are withdrawn, whereas others having a nonempty intersection with the solution set ( $\gamma ([\mathbf{x}])=0$ ) are bisected. The box $[\mathbf{x}]$ is bisected into two subboxes of equal area $[\mathbf{x}{]}_1$ and $[\mathbf{x}{]}_2$ along the dimension having the largest width. The subboxes are then tested, kept in the solution, withdrawn, or bisected until the maximal width of the resulting subboxes is less than a given threshold δ. An illustration of the proposed method is given in Figure 1. It shows four detected anchors, one of them being an outlier ( $q=1$ ). The exact solution of the problem is given in light gray, whereas the subpaving provided by SIVIA is given in both light and dark gray.

Figure 1

Robust localization with the SIVIA algorithm.

3.3. Correction Using the Combinatorial Technique

In this algorithm, both the combinatorial formulation (6) and the convex hull $[\mathbf{x}{]}^{*}(t)$ including all propagated boxes using (13) are considered. Based on (6), the proposed algorithm consists of contracting the initial domain $[\mathbf{x}{]}^{*}(t)$ with each combination $C{\in \mathbf{C}}_{I(t)}^{|I(t)|-q}$ of observations. For this reason, all observation equations indicated in C are iterated in the forward-backward contractor [15]. This contractor iterates all constraints without any prior order until no contraction is possible. The resulting region is the smallest box covering the intersection of $[\mathbf{x}{]}^{*}(t)$ with all the observation disks of C. In order to use each constraint of (2) in the forward-backward contractor, one should express $x_1(t)$ as a function of $x_2(t)$ , and vice versa as follows:

$\begin{array}{c}{a}_{i,1}\left(t\right)-b_{i,1}\left(t\right)\le x_1\left(t\right)\le a_{i,1}\left(t\right)+b_{i,1}\left(t\right),\\ a_{i,2}\left(t\right)-b_{i,2}\left(t\right)\le x_2\left(t\right)\le a_{i,2}\left(t\right)+b_{i,2}\left(t\right),\end{array}$

(15)

for each detected anchor i where $b_{i,1}(t)=\sqrt{r^2-(x_2(t)-a_{i,2}(t){)}^2}$ and $b_{i,2}(t)=\sqrt{r^2-(x_1(t)-a_{i,1}(t){)}^2}$ . Using intervals and having an initial box $[\mathbf{x}]$ , these inequalities would lead to the contracted box $[{\mathbf{x}}^{\prime }]$ defined as follows:

$\begin{array}{c}\hspace{0.17em}\hspace{0.17em}\left[x_1^{\prime }\right]=\left[a_{i,1}\left(t\right)-{\bar{b}}_{i,1}\left(t\right),a_{i,1}\left(t\right)+{\bar{b}}_{i,1}\left(t\right)\right],\\ \hspace{0.17em}\hspace{0.17em}\left[x_2^{\prime }\right]=\left[a_{i,2}\left(t\right)-{\bar{b}}_{i,2}\left(t\right),a_{i,2}\left(t\right)+{\bar{b}}_{i,2}\left(t\right)\right],\end{array}$

(16)

where ${\bar{b}}_{i,1}(t)=\sup \hspace{0.17em}([\sqrt{r^2-([x_2]-a_{i,2}(t){)}^2}])$ and ${\bar{b}}_{i,2}(t)=$ $\sup \hspace{0.17em}([\sqrt{r^2-([x_1]-a_{i,1}(t){)}^2}])$ . Then, considering the combination C of constraints, and starting with the predicted domain (initially ${[\mathbf{x}]}_C={[\mathbf{x}]}^{*}(t)$ ), the contracted box $[\mathbf{x}{]}_C$ would be obtained by performing the following steps (Algorithm 1).

Algorithm 1: Computation of the contracted box using the forward-backward contractor.

while ${[\mathbf{x}]}_C$ is contracted do

  for $i\in C$   do

     ${\bar{b}}_{i,1}(t)=\sup \hspace{0.17em}([\sqrt{r^2-([x_2{]}_C-a_{i,2}(t){)}^2}])$ ;

     $[x_1{]}_C=[x_1{]}_C\cap [a_{i,1}(t)-{\bar{b}}_{i,1}(t),a_{i,1}(t)+{\bar{b}}_{i,1}(t)]$ ;

     ${\bar{b}}_{i,2}(t)=\sup \hspace{0.17em}([\sqrt{r^2-([x_1{]}_C-a_{i,1}(t){)}^2}])$ ;

     $[x_2{]}_C=[x_2{]}_C\cap [a_{i,2}(t)-{\bar{b}}_{i,2}(t),a_{i,2}(t)+{\bar{b}}_{i,2}(t)]$ ;

  end

end

Each nonempty box $[\mathbf{x}{]}_C$ , $C\in {\mathbf{C}}_{I(t)}^{|I(t)|-q}$ , is added to the solution $[𝕏](t)$ at time t. The boxes in $[𝕏](t)$ may have nonempty intersections. Consequently, in order to obtain nonoverlapping boxes, one could apply the following procedure.

• (i)

Consider an empty final set.

(ii)

Sort all boxes of $[𝕏](t)$ according to their decreasing areas.

(iii)

Add the largest box to the final set.

(iv)

Select the following box in the sorted list.

(v)

Deprive it from all boxes already added to the final set.

(vi)

Add the result to the set.

Steps (iv) to (vi) are repeated until all sorted boxes are considered. Recall that depriving a box $[\mathbf{x}]$ from a box $[\mathbf{y}]$ yields a set of nonoverlapping boxes covering all the points x of $[\mathbf{x}]$ not included in $[\mathbf{y}]$ . An illustration of the proposed method is given in Figure 2. It shows four detected anchors, one of them yielding an erroneous observation. The first solution leads to two boxes, the box    ①    and the one in dashed line. Using the depriving technique, three nonoverlapping boxes    ①,    ②, and    ③    (in bold line) are then selected covering the exact solution (in light gray).

Figure 2

Robust localization with the combinatorial technique.

3.4. Evaluation of the Number of Outliers to Tolerate

Consider the vector $\mathit{\rho }$ of RSS measurements. Determining from $\mathit{\rho }$ the maximal number of outliers q, which have to be tolerated, may be done by choosing q such that $\text{Pr}(Q\le q\hspace{0.17em}|\hspace{0.17em}\mathit{\rho })>1-\nu$ , where $\text{Pr}(Q\le q\hspace{0.17em}|\hspace{0.17em}\mathit{\rho })$ is the probability that q or less outliers have occurred knowing $\mathit{\rho }$ and $\nu \in [\mathrm{0,1}]$ is some tuning parameter. One has

$\begin{array}{cc}\hfill \text{Pr}\left(Q\le q\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right)& ={\displaystyle \sum }_{k=0}^q\text{Pr}\left(Q=k\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right)\hfill \\ \hfill & ={\displaystyle \sum }_{k=0}^q{\displaystyle \sum }_{\left(\mathbf{s}\in {\left\{\mathrm{0,1}\right\}}^{N_a},{\sum }_i^{}{s}_i=k\right)}\text{Pr}\left(Q=k,\mathbf{s}\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right),\hfill \end{array}$

(17)

where s is some pattern indicating whether anchor $i=1,\dots ,N_a$ is providing an outlier ( $s_i=1$ ) or a reliable measurement ( $s_i=0$ ). Now,

$\text{Pr}\left(Q=k,\mathbf{s}\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right)=\text{Pr}\left(Q=k\hspace{0.17em}|\hspace{0.17em}\mathbf{s},\mathbf{\rho }\right)\cdot \text{Pr}\left(\mathbf{s}\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right).$

(18)

Let ${\rho }_i$ and ${\rho }_i^{*}$ be the noisy and noiseless RSS measurements provided by the ith anchor. The probability $p_i$ that the ith anchor provides an outlier is null ( $p_i=0$ ) if ${\rho }_i<{\rho }_r$ , since in this case, the anchor i is not detected and will thus not provide any outlier. Otherwise, $p_i$ is given as follows:

$\begin{array}{cc}\hfill p_i=\text{Pr}\left({\rho }_i^{*}<{\rho }_{\text{r}}\hspace{0.17em}|\hspace{0.17em}{\rho }_i\ge {\rho }_{\text{r}}\right)& =\text{Pr}\left({\rho }_i-{\rho }_i^{*}>{\rho }_i-{\rho }_{\text{r}}\right)\hfill \\ \hfill \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}& =\text{Pr}\left({\varepsilon }_i>{\rho }_i-{\rho }_{\text{r}}\right)\hfill \\ \hfill & =\frac{1}{2}\left(1-\mathrm{erf}\left(\frac{{\rho }_i-{\rho }_{\text{r}}}{\sqrt{2{\sigma }^2}}\right)\right),\hfill \end{array}$

(19)

where ${\varepsilon }_i$ is the ith measurement noise. Then, assuming that all measurement noise samples are independent,

$\text{Pr}\left(\mathbf{s}\hspace{0.17em}|\hspace{0.17em}\mathbf{\rho }\right)={\displaystyle \prod }_{i=1}^{N_a}\left(p_i{s}_i+\left(1-p_i\right)\left(1-s_i\right)\right).$

(20)

Now, since

$\begin{array}{c}\text{Pr}\left(Q=k\hspace{0.17em}|\hspace{0.17em}\mathbf{s},\mathbf{\rho }\right)=\text{Pr}\left(Q=k\hspace{0.17em}|\hspace{0.17em}\mathbf{s}\right)=\{\begin{array}{cc}1,\hfill & \text{if}\hspace{0.17em}\hspace{0.17em}{\displaystyle \sum }_i{s}_i=k,\hfill \\ 0,\hfill & \text{otherwise},\hfill \end{array}\end{array}$

(21)

one is able to evaluate $\text{Pr}(Q\le q\hspace{0.17em}|\hspace{0.17em}\mathit{\rho })$ and to choose q.

4. Simulations

In this section, we compare the SIVIA-based method (SBL) to the combinatorial-based one (CBL). We consider a group trajectory model, where sensors are moving along similar trajectories over $100\hspace{0.17em}\text{s}$ . We deploy 31 sensors in a $100\hspace{0.17em}\text{m}\times 100\hspace{0.17em}\text{m}$ area, 30 of them being anchors. Since nodes use only anchors information to localize themselves, we consider the localization of a single mobile node. The simulated trajectory of the node is with a maximal velocity of $2.3\hspace{0.17em}\text{m}·{\text{s}}^{-1}$ and a localization step $\mathrm{\Delta }t$ of $1\hspace{0.17em}\text{s}$ . RSS measurements are generated using the distances between the considered node and all anchors and model (1) with ${\rho }_0=100\hspace{0.17em}\text{dBm}$ , $d_0=1\hspace{0.17em}\text{m}$ , and $n_P=4$ . Moreover, r is set to $10\hspace{0.17em}\text{m}$ and ${\rho }_r$ to $60\hspace{0.17em}\text{dBm}$ . In order to compare the SBL algorithm to the CBL algorithm, we take different values of the variance ${\sigma }^2$ of the measurement noise, leading to different numbers of outliers. In fact, for each value of ${\sigma }^2$ , q is evaluated using the results in Section 3.4, with $\nu =0.1$ .

Note that the initial position of the node is supposed to be known. One may also use the whole deployment area as initial domain. All simulations are performed on an Intel(R) Core(TM)2 CPU at $2.40\hspace{0.17em}\text{GHz}$ and $1\hspace{0.17em}\text{GB}$ RAM, using MATLAB 6.1.

We first set $\sigma =3\hspace{0.17em}\text{dBm}$ , yielding either none or only one outlier per time step. Applying the results of Section 3.4 one obtains $\text{Pr}(q=0)=0.4161$ , $\text{Pr}(q\le 1)=0.9109$ , and $\text{Pr}(q\le 2)=0.9958$ . With $\nu =0.1$ , one gets $q=1$ . With $\nu =0.01$ , one would take $q=2$ , which leads to less accurate estimates (less measurements are taken into account) but still containing the solution set. Note that if too less outliers are tolerated, one may obtain empty solution sets or sets not containing the actual location of the node.

With $q=1$ , Figure 3 shows the subpavings obtained with both SBL and CBL methods. Note that the threshold δ of SIVIA is set to $1\hspace{0.17em}\text{m}$ . The plot shows that both results cover the actual position of the node. The average ratio of the areas of subpavings obtained with SBL over those obtained with CBL is equal to 0.827. SBL leads to a more accurate estimate. However, the average time required for the localization process is equal to $0.471\hspace{0.17em}\text{s}$ per time step with SBL and to $0.051\hspace{0.17em}\text{s}$ with CBL. This difference is due to the limited number of combinations considered in CBL at each time step ( $6\le |I|\le 10$ with $q=1$ ). The average number of boxes per subpaving is equal to 120 with SBL, whereas it is equal to 3 with CBL. Here, CBL is less memory consuming. Note that with a precision parameter higher than $1\hspace{0.17em}\text{m}$ in SIVIA, SBL needs less computing time but provides larger subpavings.

Figure 3

An illustration of the subpavings obtained with the SBL and the CBL algorithms.

In a second set of experiments, σ varies from $1\hspace{0.17em}\text{dBm}$ to $15\hspace{0.17em}\text{dBm}$ . Figure 4 shows the maximal number of outliers q, the total number of considered anchors $|I|$ , the average computing time per step, and the ratio of the average subpaving areas obtained with SBL over CBL as a function of σ. The simulated data are generated ten times for each σ, and the results are thus average values over the set of simulations. The plot shows that CBL is faster than SBL when the standard deviation σ is less than $8\hspace{0.17em}\text{dBm}$ . In these cases, q is less than 4, and the maximal number of considered anchors is less than 13. When the noise variance increases, the computation time of the CBL method becomes quite large compared to SBL. Choosing one algorithm or the other depends on the anchor density and on the proportion of outliers.

Figure 4

An illustration of q and $|I|$ at (a), the computing times are in (b) and the ratio of the subpavings areas (SBL over CBL) in (c).

5. Conclusion

This paper proposes and compares two techniques for mobile sensor localization that are robust to any fixed number of erroneous measurements. Using interval analysis, the estimates are sets of nonoverlapping boxes containing the actual location. The SIVIA-based algorithm (SBL) bisects the search region leading to many boxes describing efficiently the solution set; the combinatorial method (CBL) leads to larger boxes including the solution as well. In terms of computing time, CBL is more efficient than SBL for a small number of outliers, whereas the complexity of SBL is almost constant whatever the number of tolerated outliers.