Cooperative Localization and Time Synchronization Based on M-VMP Method.

Localization estimation and clock synchronization are important research directions in the application of wireless sensor networks. Aiming at the problems of low positioning accuracy and slow convergence speed in localization estimation methods based on message passing, this paper proposes a low-complexity distributed cooperative joint estimation method suitable for dynamic networks called multi-Gaussian variational message passing (M-VMP). The proposed method constrains the message to be a multi-Gaussian function superposition form to reduce the information loss in the variational message passing algorithm (VMP). Only the mean, covariance and weight of each message need to be transmitted in the network, which reduces the computational complexity while ensuring the information completeness. The simulation results show that the proposed method is superior to the VMP algorithm in terms of position accuracy and convergence speed and is close to the sum-product algorithm over a wireless network (SPAWN) based on non-parametric belief propagation, but the computational complexity and communication load are significantly reduced.

1. Introduction

In recent years, wireless sensor network has been widely used in agriculture, warehousing, production safety, emergency rescue and other fields [1]. The information which sensors collected and transmitted is valuable when combined with the sensors’ location [2,3]. Therefore, location awareness of wireless sensor networks has become one of the most important directions in the development of wireless sensor networks. Global Navigation Satellite System (GNSS) is able to provide location information for sensor nodes, but it is difficult to apply on sensor networks due to high power consumption and high cost [2,3,4]. Furthermore, the poor signal penetration capabilities of GNSS lead to inadequate location information in indoor scenes. Cooperative localization can overcome these problems through ranging and exchanging locations between neighbor nodes [5]. Over the past decade, cooperative localization in wireless sensor networks has drawn considerable attention.

In the recent ten years, cooperative positioning in wireless sensor networks has been widely concerned. The multidimensional scaling (MDS) algorithm [6,7] calculates the shortest path distance between nodes according to the connectivity of the network or distance measurement. Then, the relative coordinate diagram of all nodes is constructed by using the multidimension scale algorithm. Finally, the relative coordinate diagram is transformed into an absolute coordinate graph according to the coordinate of the anchor node. The semi-definite programming (SDP) algorithm [8,9] represents geometric constraints among nodes as a set of linear matrix inequalities, then combines the inequality into a semi-definite programming problem and obtains the location of each node by global optimization. The distance-vector (DV) hop algorithm [10,11] firstly measures the minimum number of hops between the nodes to be located and each anchor node. Then, the average distance of each hop is determined according to the number and coordinate between anchor nodes. After the jump number is converted into distance, the position of the node to be located is obtained according to the trilateral measurement method. The Approximate Perfect Point-In-Triangulation (APIT) algorithm [12,13] firstly selects three anchor nodes connected with the node to be located, and then judges whether they are in the triangle composed of anchor nodes according to the pit test, and then repeats the pit test with a different anchor node combination. Finally, the center of mass position of these triangles is taken as the coordinates of the nodes to be located. The SPAWN algorithm [14] decomposes the posterior probability density function of position variables and represents it as a factor graph. Then, the message is transmitted by particles on the factor graph to calculate the probability of edge distribution of each variable. Finally, position estimation is carried out according to minimum mean square error (MMSE) or map criterion. Although the above cooperative positioning algorithm has high positioning accuracy, it has high computational complexity and high communication cost, which seriously restricts the practical application.

For the above problems, a collaborative location algorithm based on factor graph and VMP is proposed in Reference [15], which has simple message form and small computation. However, Kullback-Leibler (KL) divergence is used to Gauss the non-Gaussian confidence of the nonlinear range measurement model. The first kind of convergence hypergeometric function minimization problem is introduced, which makes the calculation complexity very high. In Reference [16], a hard decision-based cooperative algorithm is proposed to alleviate the effect of the outliers. The geometric relationship between agent position and distance is used to avoid a large deviation coursed by geometric constraints. The authors of Reference [17] propose the distributed particle filtering evolved variational message passing (DPF-E-VMP) algorithm, which improves the convergence speed of positioning estimation by using distributed particle filtering (DPF), but this performance improvement is often accompanied by greater computational consumption. By combining the average consensus method and VMP algorithm, a joint self-localization tracking algorithm called cooperative localization with outlier constraints (CLOC) is proposed in Reference [18], which has better location performance than the separate self-localization algorithm. In Reference [19], KL divergence is minimized by the Newton conjugate gradient method, but the computational complexity is still high. But, the information loss caused by VMP parameterization will have a certain influence on the positioning accuracy. At the same time, the algorithm has been reduced in accuracy and convergence speed due to the clock synchronization between the nodes to be located.

In this paper, a VMP distributed cooperative localization algorithm based on multi-Gaussian is proposed. Based on the non-line of sight (NLOS) environment ranging model, the VMP message passing strategy based on multi-Gaussian is innovated, and the second-order Taylor expansion form of position and time synchronization joint estimation is derived. Its computational complexity is far less than the approximate solution based on KL divergence. The time complexity and communication cost are connected with the traditional VMP algorithm. However, the positioning accuracy and iteration speed have been greatly improved.

The remainder of this paper will be organized as follows: A two-dimensional (2D) wireless network is first established as a system model. Then, a traditional cooperative localization algorithm is introduced and leads to the method proposed in this paper. After that, simulation performance of the proposed method is investigated. Finally, concluding remarks are presented. A list of symbols that are used in the paper is given in Table 1.

Table 1

List of symbols.

Symbol	Meaning	Symbol	Meaning
S_i,n	Set of neighbor anchor nodes of node i at time n	C_i,n	Set of neighbor nodes to be located of node i at time n
x_i,n	Position vector of node i at time n	t˜_i,n	Measurement value of local clock of note i at time n
T_n	Real time value at time n	a_i,n	Slope of local clock of node i at time n
a_ij,n	Relative slope of local clock offset between nodes i,j	N_i,n	Set of neighbor nodes of node i at time n
ρ˜_ij,n	Range measurement between node i,j at time n	ω_ij,n	Measurement noise of ρ˜_ij,n
σd2	Variance of ω_ij,n	ζ_ij,n	NLOS error of ρ˜_ij,n
Ψ_n	Set of all communicable node pairs i,j at time n	Θ_n	Set of all node pairs i,j with NLOS error at time n
x_n	Position vector of all nodes at time n	a_n	Clock offset of all nodes at time n
v_i,n	Average velocity of node i from time n−1 to time n	ω_i,n	Measurement noise of v_i,n
σi,n2	Variance of ω_i,n	𝔃_i,n	Vector to be estimated of note i at time n
𝔃^_i,n	Estimation result of node i at time n	ρ˜_n	Range measurement in all node pairs i,j∈V_n at time n
ζ_n	NLOS error in all node pairs i,j∈Θ_n at time n	𝓧_1:n,𝓐_1:n,𝓟_1:n,𝓑_1:n	Vector sets of x_n,a_n,ρ˜_n,ζ_n from time 1 to time n
μ_f_ij,n→𝔃_i,n𝔃_i,n	Message send from f_ij,n to 𝔃_i,n	b𝔃_i,n	Belief of variable 𝔃_i,n
F·	Confluent Hypergeometric Function of the First Type	E_𝔃_i,n	Mean of belief b𝔃_i,n
V_𝔃_i,n	Covariance of belief b𝔃_i,n	η_M,ij	Weight of the M-th Gaussian distribution
F𝔃_i,n	Fisher Information Matrix (FIM) of 𝔃_i,n	CRLB𝔃_i,n	Cramer-Rao Lower Bound of 𝔃_i,n

2. System Model

The anchor nodes are always deployed at the same height, and high vertical dilution of precision means the system cannot provide reliable vertical positioning results [20], which is usually obtained by other sensors [21]. So, in this paper, a 2D dynamic wireless network is considered, which includes the anchor node with known location and synchronized local time, and the node with inaccurate location information and local time out of sync. The position vector of node at time is , and the measured value of local clock . The slope of local clock at time between node and the external standard clock is . The local clock of all anchor nodes is synchronized with the external reference clock, i.e., . At time , node has neighbor nodes set as , where neighbor anchor node set is node set to be located is , and all communicable node pairs constitute communicable node set .

Considering the local clock drift, node measures time of arrival (TOA) from neighbor node at time as follows: where is the Euclidean distance, between nodes , , is the observed measurement noise, assuming that it obeys the Gaussian distribution, i.e., , and is NLOS error, which is expressed as follows: where [22], is a constant and is the set of all node pairs with NLOS error at time . Where, is the relative slope of local clock offset between nodes , defined as follows:

Define as the position vector of all nodes, as the clock offset slope of all nodes, as the range measurement in all node pairs , and as the NLOS error in all node pairs . Moreover, all the vector sets from time 1 to time are defined as follows: , . The goal of the algorithm is to estimate the accurate node location vector and the clock slope by measuring and the information transmitted between nodes.

3. M-VMP Joint Estimation Algorithm

According to whether the location information is considered as a random variable, the cooperative localization algorithms are divided into two categories: non-Bayesian estimation and Bayesian estimation [5]. In non-Bayesian estimation, the location information is estimated by deterministic methods. The typical algorithms are least square (LS) [23,24] and the maximum likelihood method (ML) [25]. Bayesian estimation is based on the probability model of location information. The typical algorithms include the maximum posterior probability (MAP) estimation [26] and the MMSE method [6]. Location information node estimates by MMSE criteria. The expression is as follows: where is the estimation result and is the vector to be estimated.

3.1. Probability Model

Assume evolve according to a memory-less Gauss–Markov process, then, there are: where is the average velocity of node from time to time , which is measured by the sensor inside nodes. is Gaussian white noise, and its covariance matrix is . Since the motions of all nodes are independent of each other, there are: where is the prior distribution of all nodes’ location information at time 0, which is obtained by GNSS, base station or in other ways [27]. According to the Bayesian rule, the edge probability function in (3) is calculated by the following formula:

The likelihood function in (10) is decomposed into the following formula with the independent observations at different times :

According to the above derivation, (10) is expanded to the following formula:

3.2. Factor Graph Model

Factor graphs intuitively reflect the spatiotemporal relationship between variables [5]. The factor graph is a dichotomous graph, which contains two kinds of nodes: factor node and variable node. Variable node represents the information to be evaluated, and factor node represents the messages passed between variable nodes. The schematic diagram of cooperative localization is shown in Figure 1.

Figure 1

Factor graph model of variational message passing (VMP)-based cooperative localization problem. The circle represents the factor node and the square represents the variable node. represents the localization variable of the node to be located, represents the message that transmits between different times, represents the distance information that the nodes transmit, represents the non-sight distance parameter that affects the information between nodes and represents the non-sight distance error probability function that affects the non-sight distance parameter.

Factor nodes are further assumed as follows:

According to the above definition, edge probability of localization is clearly calculated in the factor graph. Next, a localization algorithm based on the VMP will be introduced.

3.3. VMP-Based Localization Algorithm

VMP uses the exponential model to deliver messages, which greatly reduces communication consumption. According to the VMP message rules based on the factor graph proposed in Reference [28] and the assumptions in Section 3.2, the message delivered from factor node to parameter node at time n is as follows: where the cooperative messages are divided into two categories: one is related to the neighbor anchor node , and the other is related to the neighbor node to be located. The calculations are performed below: where means the mean of with respect to . The belief of is obtained as follows: where is the normalization constant. Obviously, it can be seen that is not a Gaussian function about , and it is difficult to transmit directly between nodes. To reduce communication overhead, Reference [19] approximates by minimizing KL divergence: where represents a Gaussian function, and is the KL divergence between Gaussian function and belief , the calculation formula is as follows:

3.4. M-VMP-Based Joint Estimation Algorithm

In Section 3.3, the clock drift of node or model loss with single Gaussian distribution are not considered, which affects localization accuracy [25]. In this paper, a joint estimation algorithm of time synchronization and localization based on multi-Gaussian distribution VMP is proposed, that the belief of the variable to be estimated is approximately a multi-Gaussian function , and the belief of is . The factor graph between nodes at the time is shown in Figure 2.

Figure 2

FG of a node pair at time n. represent the localization and time synchronization parameter of the node to be estimated, represents the message that transmits between different times, represents the distance information that the nodes transmit, represents the non-sight distance parameter that affects the information between nodes and represents the non-sight distance error probability function that affects the non-sight distance parameter.

In order to approximate to a Gaussian function, Formula (23) is rewritten as follows: where:

Because the existence of two non-linear terms, in and in , results in non-Gaussian of , the two non-linear terms are expanded by a second-order Taylor expansion to obtain the mean and covariance of . For discussion purposes, define and .

Second-order Taylor series expansion for nonlinear term is expanded at by a second-order Taylor to get: where and are the first and second steps of at .

Second-order Taylor series expansion for nonlinear term : is expanded at and by a second-order Taylor to get: where and are the first-order partial derivatives of with respect to and at , and , is the Hessian matrix of at and .

After calculation, we can get the mean and covariance of : where , and is the normalization constant. Compared with the approximation method based on KL divergence minimization, the proposed method based on multi-Gaussian approximation with second-order Taylor series expansion greatly reduces the complexity of approximate calculation. When is approximated to a Gaussian function , each node only needs to send its own position vector and covariance matrix to its neighbor nodes, and the communication overhead is much lower than that of the particle message-based method. In addition, since each node has three parameters to be estimated, the network is required to include at least three anchor nodes. The flow of the M-VMP localization algorithm (Algorithm 1) is as follows:

4. Simulation Analysis and Test Results

In this section, Cramer-Rao Lower Bound (CRLB) of is derived first, then the performance of the proposed joint algorithms is analyzed through simulation.

4.1. CRLB Lower Bound of

For node , the Fisher Information Matrix can be derived as: with and CRLB of is .

4.2. Simulation Scenario and Result Analysis

To analyze the performance of the proposed method, we have built a simulation scenario according to the real scene of the zone 1, underground parking lot, Beijing University of Posts and Telecommunications. The real scene is shown in Figure 3a and the top view is shown in Figure 3b.

Figure 3

(a) The real scene of the underground parking lot. (b) Simulation scenario. Zone 1 is the simulation area.

According to the actual size of zone 1, the simulation scene is a rectangular area of 20 × 24 m. By default, the number of anchor nodes is 4, the number of nodes to be located is 20 and the initial position distribution of nodes in the site conforms to the uniform distribution. In order to limit the communication area of nodes to less than 1/2 of the whole simulation area, the maximum communication distance is set to 10 m. The performance of the proposed method is obtained in the line of sight (LOS) environment, except the last simulation, which shows localization accuracy with different NLOS probability. All simulation results are the average of 1000 independent runs.

The initial position measurement error of the anchor node conforms to the Gaussian distribution with the standard deviation of 0.1 m. The initial position measurement error distribution of the node to be located is composed of Gaussian distributions with the standard deviation of 10 m. According to the performance of the crystal oscillator used in the hardware (TG5032CFN), the clock drift of the node to be located conforms to uniform distribution from [−1 ppm, 1 ppm], that is, the maximum distance measurement error caused by the clock drift between adjacent time (1 s) is 300 m. In the process of simulation, anchor nodes remain stationary, and the velocity of the node to be located is consistent with the uniform distribution of the maximum value of 3 m/s, and the direction is random every time step. The position movement measurement of the node conforms to the Gaussian distribution with the standard deviation of 1 m. The default process is 15 s, 10 iterations per second. Root mean square error (RMSE) and cumulative distribution functions (CDF) are used to measure the performance of algorithms.

In Figure 4a, the relationship between localization RMSE and initial location error is shown. The CLOC method [16] treats the uncertainties of nodes’ positions and clock offsets as measurement noise and thus suffers performance degradation. The SPAWN [26] method is implemented by using 4000 particles to represent the messages on FG. The localization accuracy of the VMP method [17] and the proposed M-VMP method are better than the CLOC method, but slightly worse than the SPAWN method.

Figure 4

(a) Root mean square error (RMSE) of position error with different initial location error. (b) Cumulative distribution function (CDF) of position error with different initial location error.

The position error of VMP and the proposed M-VMP method are compared in Figure 4b. The cumulative distribution function of the VMP method and the proposed method under different initial position errors is given in Figure 4b. It is found that the probability of position error less than 1.6 m of the two methods under different initial position errors is basically the same, while the proposed method has a better effect of restraining position error. The result is consistent with the conclusion in Figure 4a.

Figure 5 shows the relationship among the clock drift slope, RMSE of localization accuracy and the number of iterations, respectively. It is seen that three methods converge in the finite iterations, which proves that the algorithms are feasible. Compared with the VMP method [17], the proposed method improves the convergence speed of localization and time synchronization by adding time synchronization parameters into the localization process.

Figure 5

(a) RMSE of the clock drift slope versus iterations and (b) RMSE of position error versus iterations.

In Figure 6, the CDF of the proposed method under different values of the number of Gaussian distribution M are compared. It is seen that the localization accuracy (3) of the proposed method is better than that of the VMP method when M is greater than 2. With the increase of M, the localization accuracy is also improved.

Figure 6

CDF of position error with different values of the number of Gaussian distribution M in line of sight (LOS) environment.

Figure 7 reflects the influence of communication distance on the localization accuracy of the proposed method under the same node density. Obviously, localization accuracy increases with the increase of communication distance, which increases the connectivity of nodes and the redundancy of information. At the same time, it can be seen that the improvement of localization accuracy is not obvious when the communication distance is 20 and 30 m. The difference of node connectivity between those two conditions is little, so the localization accuracy cannot be greatly improved. Because the communication distance is proportional to the transmission power of the node, it is necessary to determine the appropriate communication distance in practical application by considering the localization accuracy, node power consumption and other indicators.

Figure 7

CDF of position error with different communication distance.

In the simulation of Figure 8, in order to amplify the influence of node density on positioning accuracy, we reset some of the simulation parameters, and set the communication distance to 15 m and the variance of the distance observation and the initial position error of the node to be located to 20 m2. It can be seen from Figure 8 that as the density of nodes increases, the overall positioning accuracy is decreasing. However, when the density of nodes to be located is too large, the positioning error will increase. At the same time, when there are more nodes to be located, the speed at which the positioning error decreases with more anchor nodes also decreases. This is due to the increase in the density of nodes to be located, and the overall weight of the information contained in neighboring nodes to be located in the proposed method becomes larger. When the number of anchor nodes is insufficient, the increase of neighboring nodes to be located helps to improve the positioning accuracy, but the effect is limited.

Figure 8

Relationship between position error and nodes’ density.

In Figure 9, the relationship between the RMSE of localization accuracy and NLOS occurrence probability of the VMP method and the proposed method is compared. Furthermore, it compares with the method that only deals with LOS. It is seen that both the VMP and the proposed method effectively suppress NLOS error. Because the two-way NLOS parameters between nodes are added in the proposed method, the NLOS error is better suppressed.

Figure 9

(a) RMSE of position accuracy with different NLOS probability. (b) CDF of position error with different NLOS probability.

4.3. Computational Complexity and Communication Overhead Analysis

Because the proposed method in this paper and the comparison methods are distributed algorithms, the computation and communication are done by the nodes independently, so only one node’s time complexity and communication overhead need to be considered in the analysis.

The time complexity is evaluated by the number of operations in local computing, and the communication cost is evaluated by the number of information parameters broadcast by nodes. The time complexity of the CLOC method is related to its neighbor nodes number , and its communication cost is . The SPAWN method is based on particles. When the number of particles in each message is , its time complexity is , and its communication cost is . Compared with the VMP method, the number of Gaussian messages in the proposed method is M times, and every message from neighbor nodes needs to be processed once. The time complexity of dimension reduction depends on , so its time complexity is . In the proposed method, nodes transmit a localization vector including and a covariance matrix, so the communication cost is . The computational complexity, run-time and communication Overhead of the three algorithms is shown in Table 2.

Table 2

Comparisons of different methods for each node at each iteration.

Method	Computational Complexity	Run-Time	Communication Overhead
M-VMP (proposed method)	OM·N_inb+ON_inb	84.449 ms	3·OM+O1
CLOC	ONinb2	288.500 ms	3·O1
SPAWN	ON_p+Np2·N_inb	5469.089 ms	ON_p+2·O1

4.4. Future Research Directions

In the next stage of research, we will mainly focus on the following aspects: First, analyze the performance of the positioning and simultaneous joint estimation problem and its influencing factors in principle. Secondly, expand the positioning scene from 2 dimensions to 3 dimensions, and reduce the problem of positioning accuracy degradation in dense scenes with nodes to be located through methods such as signal quality screening. Finally, the proposed algorithm will be implemented based on the hardware platform, and the measured results will be compared with the simulation results to further improve the performance of the proposed method.

5. Conclusions

This paper presented a M-VMP-based TOA localization and time synchronization joint estimation algorithm for mixed LOS and NLOS environments. Firstly, according to the VMP method, a message propagation model based on a factor graph was constructed for localization and time synchronization. Owing to the existence of nonlinear terms, it is difficult to represent the message in a closed form, so Taylor expansion was used for the linearization of nonlinear terms. Moreover, to reduce the communication cost, all messages were expressed in the form of multi-Gaussian distribution, and only the mean value, covariance and weight of each Gaussian distribution need to be transferred in the message transmission. The simulation results showed that the accuracy and convergence speed of the proposed method were close to that of SPAWN method, and the time complexity and communication cost were greatly reduced: the run-time was only 1.5% that of SPAWN.