Greedy Successive Anchorization for Localizing Machine Type Communication Devices

Localization of machine type communication (MTC) devices is essential for various types of location-based applications. In this paper, we investigate a distributed localization problem in noisy networks, where an estimated position of blind MTC machines (BMs) is obtained by using noisy measurements of distance between BM and anchor machines (AMs). We allow positioned BMs also to work as anchors that are referred to as virtual AMs (VAMs) in this paper. VAMs usually have greater position errors than (original) AMs, and, if used as anchors, the error propagates through the whole network. However, VAMs are necessary, especially when many BMs are distributed in a large area with an insufficient number of AMs. To overcome the error propagation, we propose a greedy successive anchorization process (GSAP). A round of GSAP consists of consecutive two steps. In the first step, a greedy selection of anchors among AMs and VAMs is done by which GSAP considers only those three anchors that possibly pertain to the localization accuracy. In the second step, each BM that can select three anchors in its neighbor determines its location with a proposed distributed localization algorithm. Iterative rounds of GSAP terminate when every BM in the network finds its location. To examine the performance of GSAP, a root mean square error (RMSE) metric is used and the corresponding Cramér–Rao lower bound (CRLB) is provided. By numerical investigation, RMSE performance of GSAP is shown to be better than existing localization methods with and without an anchor selection method and mostly close to the CRLB.


Introduction
Machine type communication (MTC) services with location-based information are becoming popular such as combat zone surveillance, health monitoring, fire detection and wild habitat monitoring [1][2][3]. MTC devices (in this paper, we use the term "devices" and "machines" interchangeably) with which sensors are equipped as a basic component are usually distributed randomly, accessed from far-away and configured automatically. Thus, location information in MTC services is crucial. One way to get position information is to fit each MTC device with global positioning system (GPS) that enables it to get its position by interacting with GPS satellites. However, this method may not be cost-effective for large scale implementation and is not even possible for indoor or space applications, where MTC devices cannot receive GPS signals [4].
When GPS is not available, each MTC device will estimate its location based on its relative distances to some of the other MTC devices that already have acquired their positions, which are so-called anchor machines (AMs). Localization algorithms that are used to determine the unknown location of devices are computationally complicated in radar and sonar systems [5,6]. However, Covariance matrix for error inφ v which we assume a white Gaussian process b Bias of an estimator (vector) Covariance matrix of an estimator Ψ Mean square error matrix E Expectation operator Γ Fisher information matrix X D (.) Stress function used in multidimensional scaling L(.) Log-likelihood function g(.) Probability density function (.) −1 Inverse of a matrix (.) T Transpose of a vector or matrix Tr(.) Trace of a matrix

System Model
We consider an MTC network that consists of n BMs and m AMs or VAMs in a two-dimensional space (this work also can be straightforwardly applied for localization in a three-dimensional space). Let N = n + m and φ i = [x i , y i ] T (∈ R 2 ) denote the actual coordinates of machine i, and let S A , S V and S U be index sets of AMs, VAMs and BMs, respectively. We assume that each AM or VAM knows its position as well as its positioning error ϑ p j (j ∈ S A ∪ S V ). It is reasonable to assume that the positioning error of AM is less than or at least equal to that of VAM.
In the MTC environments, each machine is capable of communication with the other machines if they are in the communication range. We use RSS measurement to determine unknown positions of BMs. We assume that each AM and VAM broadcasts reference signals as well as its coordinates with the positioning error (and sometimes with a value of its transmit power). BMs that can hear the information from three or more anchors first select the anchor nodes and then run a localization algorithm to be provided below. It is noted that the localization algorithm is distributed in nature, that is, each BM independently determines its position without cooperation by the other BMs.
We assume that transmit power (denoted by P j , j ∈ S V ∪ S A ) for the broadcast information is known to BMs. Let P ij denote RSS measurement at BM i from anchor j. Then, we can use a path loss model in [18]: where ς ij (φ) = −10β log d ij , β is a distance-power gradient (i.e., a path loss exponent), d ij is the actual distance between machines i and j, and η ij is noisy power due to measurement errors at BM i and unmodeled variability in the fading channel between machines i and j (for example, shadowing). We assume that η ij s are independently distributed zero-mean Gaussian random variables with standard deviation σ ij (dB). We also assume that σ ij is known to BMs. From the RSS measurements, each BM can compute the expected error in measurement denoted by ϑ r ij (i ∈ S U , j ∈ S V ∪ S A ). If BM i receives L RSS samples from anchor j, ϑ r ij = 10 σ ij /10 /L since noisy power is assumed to be independent.

Greedy Selection
BM i selects anchors to be used in the localization process by comparing where j is an index of anchors in its communication range, from which the BM successfully receives the broadcast information. Since three anchors are enough to determine the location of BM in two-space, each BM selects anchors up to smallest ζ ij , which is referred to as greedy selection in this paper.

Removing Collinear Anchors
If anchors selected from Section 3.1.1 are collinear, high localization error is caused mostly due to vertex flipping instead of the position or the measurement error [19,20]. In order to overcome such a problem, we test a sin 2 γ > d min , where a and γ are the shortest side and the smallest angle of the triangle built by the greedy-selected three anchors, and d min is a predetermined constant. This test is provided in [19] as a method of constructing a so-called robust triangle. If the three anchors fail to pass the test, we select another anchor that has the next smallest ζ ij and do the same test for every combination of the selected anchors. If the test fails for all the anchors whose information is available at a BM, the BM detours the localization procedure for the next round in which more anchors are possibly available since new VAMs can be added to the list of anchors.

Other Anchor Selection Methods in the Literature
Though all of the existing anchor selection methods do not consider the position error at an anchor, there are some methods in the literature. In CL [7], BM selects all those anchors that are one-hop neighbors of that BM. In [21], in order to increase the accuracy, a convex-hull anchor selection method is proposed by considering the geometry of anchors, where only those anchors that have the greatest distance from each other are selected. When applying the convex-hull method with RSS measurement, possibly large measurement error is problematic since the increasing distance between anchors usually increases the distance between an anchor and BM, and thus the error [22]. The work of [21] has proposed an advanced convex-hull method where BM selects those AMs that are close to the ordinary convex-hull with the highest virtual location accuracy. The main drawback of this method is that it assumes highly dense anchors. We will numerically compare the performance of proposed greedy selection with that of the existing methods in Section 5.

Iterative Localization Algorithm
Let S S be a set of selected anchors and p i = {p ij = P ij − P j , j ∈ S S } be a column vector of augmented RSS measurement (or observed path loss) at BM i from the selected anchors. BM i (i ∈ S U ) determines (or estimates) its positionφ i bŷ where Furthermore, ς(φ i ) is a column vector, the element of which is A general solution of Equations (3) and (4) can be obtained through a nonlinear optimization method [23][24][25], which iterates to get the optimal value in the feasible region. However, its computational complexity is usually very high [26]. In order to make the computational cost low [27], we linearize ς(φ i ) by using the first-order Taylor series expansion [28] where andφ i is any initially estimated position of BM i (finding a goodφ i will be shortly discussed in Section 3.3). In the above expression and in the sequel, we assume that 1 to j are indexes of the selected anchors for notational brevity. Putting Equation (6) into Equation (4), we have In [29], it is shown that the localization accuracy can be increased by inserting a symmetric weighting matrix R into Equation (8), and the resulting equation is defined as weighted least square (WLS), which has the form Using Equation (9), an approximate solution of Equation (3) is obtained as As a weighting matrix, we use R = diag ζ i1 , ζ i2 , · · · , ζ ij , where ζ ij s are the metrics used in anchor selection in this paper.
The proposed method stops with a solutionφ * i if Θ (WLS) φ * i − Θ (WLS) φ i < δ for a given positive termination threshold δ. Otherwise, we update the initial estimate asφ i =φ * i and find the next solution in Equation (10) again. At the termination, BM i can estimate its position error by ϑ p = 1/Tr((Π T RΠ) −1 ), which will be broadcast by BM i as a VAM after finding its position.

Initial Location Estimation
If BM has received information from a sufficient number of AMs or VAMs and starts to determine its location, it needs a certain initial estimate of its location when the linearization provided in Equation (6) is evoked. Though the proposed algorithm in Section 3.2 works well with any arbitrary initial point, we describe a method of finding a good initial point in the following and will show that it improves the convergence speed significantly in Section 5. The following method is based on MDS in [17] and implemented in a distributed way.
Let N − 1 be the number of AMs or VAMs, the signals from which BM can hear, and let index 1 denote the current BM that wants to have its initial location. In addition, 2, 3, · · · , N indicate the AMs and the VAMs for notational simplicity. Letd ij be the distance estimate between machine i and j. Since BM 1 can measure RSS from the anchors, we can have, from Equation (1), Since the BM can see the position of anchors in the broadcast information, letting (x j ,ỹ j ) be the broadcast coordinate from anchor j, it can estimate the distance between those anchors bỹ If we further letd jj = 0, then we can construct a so-called proximity information matrix (PIM) D=[D] ij =d 2 ij . It is noted that D is an N × N square symmetric matrix. Using PIM D, an objective of MDS is to find relative coordinatesΦ of the machines in the hearing range of BM 1, such that d ij (Φ) = φ i −φ j matchesd ij as well as possible. For this purpose, an optimization problem to be solved is where is a so-called stress function provided in [30]. A solution is obtained by decomposing a double centered matrix Ω so as to have Ω ∆ =Φ TΦ [31], where each element of Ω is constructed by It is noted that the rank of Ω is 2. Let Λ = diag(λ 1 , λ 2 , · · · , λ N ) be a diagonal matrix of the eigenvalue λ j s of Ω. In addition, let C = [c 1 , c 2 , ..., c N ] be a matrix of the eigenvectors c j corresponding to λ j . Then, taking eigendecomposition of Ω gives Ω = CΛC T , and we have a solutioñ (16) and can useφ 1 as an initial position of BM.

Cramér-Rao Lower Bound
In this section, we derive a CRLB of the localization error for multiple n b BMs that commonly use m AMs and n v VAMs. Let φ b denote the position vector of BMs. We assume that AMs have their locations without position error and without loss of generality. Let us denote actual and estimated positions of VAMs by φ v andφ v , respectively. The errors in measurement p as well as position φ v are Gaussian-distributed and mutually independent. The log-likelihood function of data vector T represents the unknown vector, K is a constant that is independent of the unknowns, and Q φ v is a covariance matrix of position errors inφ v . Letˆ =ˆ (x) denote a position estimator. Then, the error performance can be determined by its bias ) T and its mean square error Ψ(ˆ ) = (ˆ ) + b(ˆ )b T (ˆ ) [28,32]. When RSS noise is sufficiently small, an estimator could be unbiased b(ˆ ) = 0 and Ψ(ˆ ) = (ˆ ). In this case, the following inequality also holds where Γ is the fisher information matrix (FIM) [32]. For multiple BMs, equivalent FIM (EFIM) is developed in [33] and has a form of where By using a formula of partitioned matrix inversion similarly in [28], we can obtain Υ −1 is similar to the CRLB of φ b given in [34] when an exact value of φ v is known. Consequently, the second term in Equation (21) is the cost in CRLB when φ v is not known or known with errors. In [35], CRLB when n b = 1 is provided, which is a special case of Equation (21).

Numerical Results
For numerical investigation, we assume 17 AMs, and 120 BMs are randomly distributed inside a square-region: {(x, y) : 0 ≤ x ≤ S, 0 ≤ y ≤ S}, where S = 15 unit of length. We use two different configurations of AM distribution: (a) uniformly distributed over the whole area (as shown in Figure 1a) and (b) uniformly distributed in the lower left corner (in Figure 1b). We assume a distance-power gradient of β = 2, a standard deviation of measurement of σ ij = 2 dB and transmit power of P j = 40 dBm for all of the transmitting machines. We also assume that BMs can hear the signals from those AMs or VAMs, the actual distance to which is not greater than R, which is referred to as a communication range. As a performance metric, we use a normalized RMSE where Q is the number of Monte Carlo runs,φ (i) j is the estimated coordinate of BM j in the i-th run whose actual coordinate is φ (i) j and n is the total number of BMs. For each result in the figures below, Q = 10, 000 runs were done. The performance of proposed GSAP is compared with that of existing methods CL [7], WCL [8], LLS [10], subspace [11], and the CRLB provided in Equation (21). For a termination criterion of the localization algorithm, we have used δ = 10 −4 .
We have suggested an algorithm of determining a good initial position of BMs using MDS. Figure 2 compares convergence speeds of the localization algorithm provided in this paper when a random initial position or a good initial position is used. For the test, an actual location of BM φ = [4.0211, 6.5256] T is investigated and R is assumed to be 2. A randomly generated initial point is φ = [0.40171, 2.3843] T and an initial point obtained from MDS isφ = [3.3503, 5.7309] T . From Figure 2, we can see that the convergence is faster if the initial point from MDS is used. When we take δ = 10 −4 for a stopping criterion of the algorithm, a good initial point terminates the localization at nine iterations compared to 20 iterations with the random initial point. Both are convergent at the same position [4.0178, 6.5206] T . In Figure 3a,b, the effect of AM distributions on the localization speed is investigated. They plot the number of machines marked as anchors (AMs or VAMs) and BMs at the start of each round in GSAP. In Figure 3a, at round 5, all of the 120 BMs can determine their locations when R = 2. If R = 4, only three rounds are enough. If AMs are located in the lower left corner, seven and four rounds are needed to determine the location of 120 BMs when R = 2 and R = 4, respectively. It is seen that a biased-distribution of AMs can severely degrade the localization speed, especially when the communication range is relatively small.  Figure 4 shows the effect of σ dB = σ ij on the performance of localization algorithms. R = 2 is assumed. In the simulation, identical σ dB is applied for each transmitter-receiver pair. For CL and WCL, anchor selection methods (ASM) provided in [20] are used, and the greedy anchor selection (GAS) proposed in this paper is used for LLS and Subspace techniques. The RMSE degrades (gets greater) as σ dB increases for every algorithm tested. The proposed GSAP, however, provides the best performance regardless of the assumptions on the AM distribution. Comparing the results from the two AM distributions, we can see that if AMs are uniformly distributed over the whole area, RMSE gain is about 4.165 dB compared with the locally distributed case at σ dB = 3. In both of the cases, the performance gap between GSAP and CRLB increases according to the increasing uncertainty.   Figure 5 shows the effect of communication range R on the RMSE performance. When the communication range increases, each BM has the reference information from more anchors. However, since the location information from VAMs has some errors, the larger number of anchors is not always helpful. To make the comparison fair, we also have applied similar weight matrices, which are used in WCL, to the LLS and subspace techniques, referred to as weighted LLS and weighted subspace in Figure 5, respectively. It can be seen that, without anchor selection, the RMSE performance degrades though R increases in both Figure 5a,b. However, with anchor selection, it is improved with increasing R. If R is large enough (for example, when R ≥ 10), most of the BMs can communicate with more than two anchors at the first round, and the RMSE can be reduced though anchor selection is not used. Comparing the performances with and without anchor selection, the localization algorithms with certain anchor selection certainly achieve better RMSE than without the anchor selection if the communication range is limited.

Conclusions
In this paper, we have proposed and investigated GSAP. In an MTC network that consists of a certain number of AMs and many BMs, GSAP eventually estimates the location of all the BMs. Though some BMs cannot find enough anchors in their respective neighbors in an initial round of GSAP, they can find anchors in the following rounds by regarding positioned BMs as VAMs. In this procedure, RMSE is, of course, increasing by the positioning error in VAMs, but it is numerically shown that GSAP provides RMSE better than the existing methods (CL, WCL, LLS and subspace techniques).
Specifically, the proposed method achieves about 52%-97% improvement in terms of RMSE when AMs are uniformly distributed over the area. If AMs are located in a corner, the improvement increases to about 57%-98%. To reduce the effect of positioning error in VAMs, GSAP selects and uses only the three best anchors in terms of the errors. The RMSE performance of GSAP is also shown to be close to the CRLB. GSAP is a generic procedure in finding the location of BMs in a relatively large area with a limited number of anchors and is promising in determining the location of MTC devices in future Internet of Things applications.