# A Mobile Localization Method in Smart Indoor Environment Using Polynomial Fitting for Wireless Sensor Network

Authors: Yan Wang; Yang Yan; Zhengjian Li; Long Cheng
Journal: Journal of Sensors
 (2020)

---

## Abstract

The main factor affecting the localization accuracy is nonline of sight (NLOS) error which is caused by the complicated indoor environment such as obstacles and walls. To obviously alleviate NLOS effects, a polynomial fitting-based adjusted Kalman filter (PF-AKF) method in a wireless sensor network (WSN) framework is proposed in this paper. The method employs polynomial fitting to accomplish both NLOS identification and distance prediction. Rather than employing standard deviation of all historical data as NLOS detection threshold, the proposed method identifies NLOS via deviation between fitted curve and measurements. Then, it processes the measurements with adjusted Kalman filter (AKF), conducting weighting filter in the case of NLOS condition. Simulations compare the proposed method with Kalman filter (KF), adjusted Kalman filter (AKF), and Kalman-based interacting multiple model (K-IMM) algorithms, and the results demonstrate the superior performance of the proposed method. Moreover, experimental results obtained from a real indoor environment validate the simulation results.


---

## Body

## 1. Introduction

A wireless sensor network (WSN) is a distributed sensor network which consists of numerous tiny sensor nodes. Numerous research studies have been done on WSN, particularly on localization, since WSN is a significant way to approach the localization problem. It can be used to estimate the position of a mobile target combining distance measurements to multiple nodes of the WSN. In general, nodes in a WSN localization context are mainly classified as beacon nodes (the location is known) and unknown nodes (the location is unknown). The main task of WSN localization is to calculate the coordinates of the unknown node via biased measurements obtained from the beacon nodes. Localization techniques of WSN can be divided into two types: range-based localization algorithms [1] and nonranging (range-free) based localization algorithms [2]. Some classic methods of range-based indoor localization technology are time of arrival (TOA) [3, 4], time difference of arrival (TDOA) [5, 6], angle of arrival (AOA) [7], and received signal strength (RSS) [8, 9]. Some advanced approaches which combine the methods above (data fusion) are also applied practically [10–12]. The proposed method in this paper employs TOA in consideration of its low complexity and high accuracy [13].Distance measurements are prone to be affected by NLOS. If the propagation condition between two types of nodes is line of sight (LOS), some classic localization approaches such as Kalman filtering (KF), which assume estimates are Gaussian, unbiased, and homoscedastic, could achieve high accuracy. However, in indoor environments, the signal between two types of nodes may be reflected or diffracted by obstacles. This type of scenario is classified as “nonline of sight” (NLOS). A portion of the measurements produced in the presence of NLOS are incorrect (due to positive deviation) [14–16], which presents a problem when designing localization algorithms. Moreover, NLOS obeys different statistical distributions for different scenes.Some methods have been applied to improve the applicability of KF, which effectively alleviated fluctuation of RSS caused by NLOS [17, 18]. Apart from this, some methods have appeared to further mitigate the NLOS error. Those methods can be divided into two categories, one needs to conduct NLOS identification, while the other does not.Recently, algorithms for computing the position of mobile nodes without NLOS identification are studied [19–21]. Chen proposed the residual weighting algorithm (RWGH) [19]. The reciprocal of the residuals is employed as the weight to perform weighted summation, and NLOS error is significantly attenuated. Considerable works deal with nonlinear and non-Gaussian problems in localization framework. One of them is the hierarchical voting-based mixed filter (HVMF) which mixes the Square Root Unscented Kalman Filter (SRUKF) and the particle filter (PF) [20]. HVMF can efficiently reduce the effect of NLOS errors no matter the amount of measurement noise. However, since the RWGH algorithm exhaustively combines measured values and HVMF is a combination of several algorithms, they consume a lot of time and energy. Methods based on the interacting multiple model (IMM) achieves good performance too, which are also some typical cases without NLOS identification. A Kalman-based IMM filtering algorithm has been proposed [21]. IMM is employed as a switch between the line of sight (LOS) and non-LOS (NLOS) states with two interactive modes. The LOS/NLOS smoothers, which correspond to two modes, are used simultaneously. Then, results obtained from both smoothers are combined. Compared to HVMF, K-IMM takes less computation time; however, its computational complexity is higher than methods with NLOS identification and switching to a smoother channel.The methods for NLOS identification to compute the position of mobile nodes have been studied [14, 22–26]. It is necessary to determine whether the current measurement is obtained under NLOS conditions. If the identification of the propagation condition is correct, the localization accuracy can be achieved to a high level. A least square Support Vector Machine-based NLOS identification algorithm [14] shows nice performance in mitigating NLOS-affected UWB data. Moreover, the Wylie algorithm [22] is one of the aforementioned NLOS identification methods [22–24]. It fits several historical measurements and compares measurement noise with the root mean square residual of the fitted value and measurement. Another method conducts NLOS identification via comparing measurement noise with the residual of prediction and measurement [25]. These two methods above have achieved good positioning accuracy. They tried to correct the NLOS error through two different NLOS identification methods. The Wylie algorithm does not require prior information of NLOS error; however, it is not sensitive enough when LOS/NLOS switching. The method in [25] requires prior NLOS information because it has to obtain the parameters of RSS indication. Therefore, in order to make a compromise between sensitivity and robustness, this paper states a TOA-based fitting NLOS identification via distance prediction.Some methods employ fitting for NLOS mitigation. The modified Kalman filter algorithm is utilized to reduce NLOS error, and on the basis of it, the least squares method (LSM) is employed to reconstruct the measured value [26]. In the method stated in this paper, a curve is fitted through LSM based on the historical measurements, and this curve is used not only for NLOS identification but also for alleviation of NLOS errors, improving localization accuracy without consuming too much computational time. Moreover, the algorithm in [26] is only used for the stationary localization of WSN nodes, while the proposed method is designed for mobile targets.In this paper, a method is proposed to deal with NLOS noise in an indoor localization framework. The method achieves a balance between low complexity and high accuracy, without requiring prior information on NLOS errors. Simulation on various node dispositions and noise distributions demonstrates good performance of the method, particularly under large NLOS noise cases. Besides, the proposed method is suitable for 2-dimensional situations, as extension to 3-dimensional situations is rather straightforward.The rest of the paper is framed as follows. In Section2, the system model is illustrated. In Section 3, the proposed method is discussed in detail. Simulation analysis and experiment results are shown in Section 4. Section 5 draws conclusions.


## 2. System Model

As shown in Figure1, it is assumed that M beacon nodes are randomly deployed in the radio propagation environment to detect the signal from the mobile node. The coordinate of the m-th beacon node BSm is given as follows, ξm=xm,ymT for m=1,2,⋯,M, and the location of the obstacles is unknown. The mobile node MS moves randomly in the measuring site, whose coordinate at time k is Zk=xk,ykT for k=1,2,⋯,K. This paper considers a 2-dimensional scenario.Figure 1
Data exchanges between beacon node and mobile node.The beacon nodes receive the signal sent from the mobile node, and the distance information is converted from it. Synchronization amongBSs and MS is assumed for analytical convenience [27]. The estimation model based on the TOA method in LOS propagation conditions is modeled as follows [28]:
(1)tTOA=tTOAo+nTOA,where tTOAo is the real TOA between the BSm and the MS and nTOA is the measurement noise modeled as zero-mean white Gaussian with variance σ2.The true distance between them-th beacon node BSm and the mobile node MS at time k is written as
(2)dmok=xm−xk2+ym−yk2.The measured distancedmk corresponding to the TOA data between the BSm and the MS at time k in LOS propagation conditions is modeled as [29]
(3)dmk=c⋅tTOA=dmok+nm,where c is the propagation velocity of the wireless signal and nm is the measurement noise modeled as zero-mean white Gaussian with variance σm2.In NLOS propagation conditions, since the obstacles are located between theBSm and the MS, the signal arrives at BSm after reflection and diffraction effects. Thus, the distance measurement dmk corresponding to the TOA data between the BSm and the MS at time k in NLOS propagation conditions is modeled as [30, 31]
(4)dmk=dmok+nm+bNLOS,where bNLOS is the NLOS error which is assumed independent with the measurement noise nm, and bNLOS is positive on account of longer path of indirect propagation. The NLOS error bNLOS is assumed to obey a Gaussian, uniform, or exponential distribution, and its distribution and parameters vary in different indoor radio environments and measurement methods.The probability density function (PDF) ofnm can be expressed by
(5)fG0nm=12πσm2exp−nm22σm2.PDF ofbNLOS when it obeys the Gaussian distribution bNLOS~NμNLOS,σNLOS2 is given by
(6)fGNbNLOS=12πσNLOS2exp−bNLOS−μNLOS2σNLOS2.PDF ofbNLOS when it obeys the uniform distribution bNLOS~Uumin,umax is given by
(7)fUbNLOS=1umax−umin,umin≤bNLOS≤umax,0,else.PDF ofbNLOS when it obeys the exponential distribution bNLOS~Eλ is given by
(8)fEbNLOS=λ−1exp−bNLOSλ,bNLOS≥0,0,else.


## 3. Proposed Method

As shown in Figure1, the proposed method is performed with M beacon nodes. Its input is the distance measurement dmk, and output is the coordinate estimation Z^k=x^ky^kT of mobile node MS.


### 3.1. General Concept

From the distance measurement models, a state vector can be obtained according to theBSm which is defined by θmk=dmokvmokT for k=1,2,⋯,K, where vmok is the velocity of the MS according to the BSm. The state space model of the MS signal under LOS/NLOS environment can be described as
(9)θmk+1=Fθmk+BaMSk,fork=1,2,⋯,K,where
(10)F=1T01,B=T22T,T=tk−tk−1 is the sampling period, and aMSk is the random process noise due to the acceleration of the MS that is modeled as zero-mean white Gaussian distribution with diagonal covariance matrix Q=σm2I. The error covariance of the system is defined as
(11)Pmk=Eθmk−θ^mk∣kθmk−θ^mkkT,where θ^mkk is the estimation of the state vector θmk and E represents the mathematical expectation of the expression contained within.The measurement equation of theBSm under LOS/NLOS environment can be described as
(12)zmk=Hθmk+Vm,where zmk is the observation vector according to the BSm, H is the observation matrix, Vm is the observation noise from the BSm that is modeled as
(13)Vm=Vm~N0,σm2,LOS,Vm~NμNLOS,σm2+σNLOS2,NLOS.As shown in Figure2, the proposed method consists of four main stages:
(1)
NLOS Detection through Polynomial Fitting. After obtaining the distance measurement dmk at measurement time k, historical measurement data at time k−1 to time k−L are selected. Thus, a total of L measurement data is used to execute polynomial fitting. Then, a curve s^mLk−1 is obtained from it. The fitted curve gives the distance prediction value dmpk of the mobile node MS at time k. dmk and dmpk two values are subtracted to obtain a residual. Then, the residual is employed to identify whether the signal propagation condition between the BSm and the MS is LOS or NLOS.
(2)
Adjusted Kalman Filter. Part of the state space model is adjusted with the primary and the quadratic coefficient of the fitted curve. Then, the Kalman filter is performed based on the adjusted state space model and obtains the state estimation θ^mkk. Finally, the distance estimation d^mk is obtained from the state estimation θ^mkk.
(3)
Combination through Weighting Filter. In this stage, if the propagation condition is judged as NLOS in the first stage, the state estimation θ^mkk obtained from the previous stage and the fitted curve s^mLk are weighted to obtain a weighted state estimation θ^mWkk. Finally, the method yields total distance estimation d^mk in the case of LOS/NLOS.
(4)
Location Estimation. With the distance estimation d^mk obtained from the previous stage, the coordinate estimate Z^k=x^ky^kT of the MS is calculated in this stage by the maximum likelihood estimation localization method.Figure 2
Flowchart of proposed method.


### 3.2. NLOS Detection through Polynomial Fitting

NLOS detection is implemented by quadratic polynomial fitting which is performed on the lastL historical measurement data based on the current measurement time. Then, the distance prediction of the mobile node MS is obtained by the curve fitted from previous polynomial fitting, which is illustrated in Figure 3. A residual is obtained from the prediction value and the distance measurement at the next measurement time. Finally, the residual is used as the threshold to identify whether the signal propagation condition between the beacon node BSm and the MS is LOS or NLOS.Step 1.
A curve is fitted using the lastL historical measurement data according to the BSm, and the curve can be described as s^mLk=∑n=0N−1a^mLnkn, where a^mLn is the nth coefficient of the fitted polynomial, m represents the BSm, and k is the time step. The unknown coefficients are obtained via least squares method:
(14)I=min∑k=0L−1s^mLk−dmk2=min∑k=0L−1∑n=0N−1a^mLnkn−dmk2.
I=min∑k=0L−1∑n=0N−1a^mLnkn−dmk2 is the multivariate function of a^mL0,a^mL1,⋯,a^mLN−1. Since a function has a partial derivative at a point and has an extreme value at that point, ∂I/∂a^mLα=0, for α=0,1,⋯,N−1; then, a^mL0,a^mL1,⋯,a^mLN−1 can be solved. Thus, we obtain the fitted curve:
(15)s^mLk=∑n=0N‐1a^mLnkn.Step 2.
After obtaining distance measurementdmk at time k, L measurements dmk−L⋯dmk−1 are selected to fit a curve s^mLk−1. Then, the distance between the BSm and the MS at time k can be predicted by the fitted curve. The one step distance prediction obtained by the BSm at time k is
(16)dmpk=s^mLkk−1,where s^mLkk−1 represents the value of the curve s^mLk−1 when horizontal coordinate is k. Then, we calculate the residual by
(17)empk=dmpk−dmk.
Based on the premise that the standard deviation of the distance measurement samples in NLOS is larger than that in LOS, it can be predicted whether there exists NLOS propagation between theBSm and the MS at time k:
(18)H0:empk<n⋅σm,H1:empk≥n⋅σm,2≤n≤5.
The hypothesisH0 holds true, if the propagation condition between the BSm and the MS is LOS at time k. The alternative H1 holds when the condition is NLOS. n is a positive integer.
In the simulation,n=3 is assumed. When H0 holds true, the probability of dmk in the interval dmok−3σm,dmok+3σm is Pdmk−dmok<3σm. According to (3) and (5), Pdmk−dmok<3σm is P=∫−33fxdx=1−2Φ3=0.9974. In this case, it is a small probability event when dmk is outside the interval dmok−3σm,dmok+3σm.Figure 3
Schematic diagram of polynomial fitting.


### 3.3. Adjusted Kalman Filter

In the Kalman filter section, part of the state space model is adjusted as follows.The linear term coefficient of the expression of the fitted curve can be understood as the velocity of the mobile nodeMS; that is, v^mk=a^mL1. Thus, the state vector can be rewritten as
(19)θmk=dmoka^mL1T,fork=1,2,⋯,K.Similarly, the quadratic coefficient of the fitted curve can be assumed as the acceleration of theMS; that is, a^MS=a^mL2. Thus, the state space model can be rewritten as
(20)θmk=Fθmk−1+Ba^mL2.Then, the Kalman filter is performed based on this model.Step 1.
Kalman Prediction. The state prediction and the prediction error covariance matrix can be obtained by the Kalman filter based on the rewritten state space model:(21)θ^mkk−1=Fθ^mk−1k−1,(22)Pmkk−1=FPmk−1k−1FT+σm2BBT,where θ^mkk−1 is the prediction of the state vector, θ^mk−1k−1 is the estimation of the state vector at time k−1, Pmkk−1 is the prediction error covariance matrix, and Pmk−1k−1 is the estimation error covariance matrix at time k−1. The Kalman residual is
(23)emk=zmk−z^mkk−1,(24)z^mkk−1=Hθ^mkk−1,where H=10. The Kalman gain is
(25)Kmk=Pmkk−1HTHPmkk−1HT+V−1.
It is employed to minimize the covariance of the distance prediction errorEemkemTk. V is the observation error covariance matrix.Step 2.
Kalman Update. The state estimation and the estimation error covariance matrix can be obtained by the following formula:(26)θ^mkk=θ^mkk−1+Kmkemk,(27)Pmkk=Pmkk−1−KmkHPmkk−1HT+VKmkT.
The distance estimated^mk between the beacon node BSm and the MS can be obtained by the state estimate θ^mkk:
(28)d^mk=Dθ^mkk,LOS,Dθ^mWkk,NLOS,where D=10, θ^mWkk is the weighted value of the state estimate under NLOS, and the weighting filter method will be introduced in the next section.


### 3.4. Combination through Weighting Filter

In the NLOS propagation condition, in order to alleviate the NLOS error, a weighted summation can be employed with the fitted curves^mLk and the state estimate θ^mkk. s^mLk is obtained from dmk−L+1⋯dmk, which is different from the curve s^mLk−1 calculated in NLOS detection stage.First, for further trust on the fitted value when the filtered value badly deviates from the real distance,wmk=exp−θ^mkk/s^mLkk is defined as the weight, where s^mLkk represents the value of curve s^mLk when horizontal coordinate is k. Then, the weighted sum can be expressed as
(29)θ^mWk∣k=s^mLkkv^mk1−wmk+θ^mkkwmk,where v^mk=Eθ^mkk,E=01. Thus, the distance measurement dmk, according to the beacon node BSm at time k, converts into the distance estimate d^mk which can be described as
(30)d^mk=Dθ^mkk,LOS,Ds^mLk∣kv^mkTwmk+θ^mkk1−wmk,NLOS,Algorithm 1: PF-AKF.
1:fork=1,2,⋯,Kdo
2: obtain distance measurementdmk.
3: Calculate curves^mLk−1 and s^mLk by polynomial fitting using Eq. (14) and Eq. (15).
4: Obtain distance predictiondmpk+1 using Eq. (16).
5: NLOS detection by hypothesis test using Eq. (18).
6: Adjusted Kalman Filter using Eq. (19) and Eq. (20).
7: Kalman prediction using Eq. (21) and Eq. (22).
8: Kalman update using Eq. (26) and Eq. (27).
9:if NLOS detected then
10:  calculate weighted state estimationθ^mWkk by weighting filter using Eq. (29).
11:end if
12: obtain distance estimationd^mk by combination using Eq. (28).
13: calculate coordinate estimationZ^k by maximum likelihood localization using Eq. (34).
14:end for


### 3.5. Location Estimation

In this section, maximum likelihood localization method [32, 33] is briefly introduced. As mentioned above, the coordinates of the beacon nodes are x1,y1,x2,y2,⋯,xM,yM. The coordinate of the mobile node MS at time k is Z=xk,ykT for k=1,2,⋯,K. The output of the proposed method is d^mk. Then, the following equation set is obtained:
(31)x1−xk2+y1−yk2=d^1k2,⋮xM−xk2+yM−yk2=d^Mk2.The equation set above can be rewritten as follows:(32)2xkx1−x2+2yky1−y2=d^2k2−d^1k2+x12+y12−x22+y22,⋮2xkx1−xM+2yky1−yM=d^Mk2−d^1k2+x12+y12−xM2+yM2.It can be represented by a linear equationAZ^k=B, where A and B are given by
(33)A=2x1−x2y1−y2x1−x3y1−y3⋮⋮x1−xMy1−yM,B=d^2k2−d^1k2−x22+y22+x12+y12d^3k2−d^1k2−x32+y32+x12+y12⋮d^Mk2−d^1k2−xM2+yM2+x12+y12.The coordinate matrix of theMS is calculated as follows
(34)Z^k=x^ky^kT=ATA−1ATBThe overall PF-AKF procedure is summarized in the pseudocode in Algorithm1.


## 3.1. General Concept

From the distance measurement models, a state vector can be obtained according to theBSm which is defined by θmk=dmokvmokT for k=1,2,⋯,K, where vmok is the velocity of the MS according to the BSm. The state space model of the MS signal under LOS/NLOS environment can be described as
(9)θmk+1=Fθmk+BaMSk,fork=1,2,⋯,K,where
(10)F=1T01,B=T22T,T=tk−tk−1 is the sampling period, and aMSk is the random process noise due to the acceleration of the MS that is modeled as zero-mean white Gaussian distribution with diagonal covariance matrix Q=σm2I. The error covariance of the system is defined as
(11)Pmk=Eθmk−θ^mk∣kθmk−θ^mkkT,where θ^mkk is the estimation of the state vector θmk and E represents the mathematical expectation of the expression contained within.The measurement equation of theBSm under LOS/NLOS environment can be described as
(12)zmk=Hθmk+Vm,where zmk is the observation vector according to the BSm, H is the observation matrix, Vm is the observation noise from the BSm that is modeled as
(13)Vm=Vm~N0,σm2,LOS,Vm~NμNLOS,σm2+σNLOS2,NLOS.As shown in Figure2, the proposed method consists of four main stages:
(1)
NLOS Detection through Polynomial Fitting. After obtaining the distance measurement dmk at measurement time k, historical measurement data at time k−1 to time k−L are selected. Thus, a total of L measurement data is used to execute polynomial fitting. Then, a curve s^mLk−1 is obtained from it. The fitted curve gives the distance prediction value dmpk of the mobile node MS at time k. dmk and dmpk two values are subtracted to obtain a residual. Then, the residual is employed to identify whether the signal propagation condition between the BSm and the MS is LOS or NLOS.
(2)
Adjusted Kalman Filter. Part of the state space model is adjusted with the primary and the quadratic coefficient of the fitted curve. Then, the Kalman filter is performed based on the adjusted state space model and obtains the state estimation θ^mkk. Finally, the distance estimation d^mk is obtained from the state estimation θ^mkk.
(3)
Combination through Weighting Filter. In this stage, if the propagation condition is judged as NLOS in the first stage, the state estimation θ^mkk obtained from the previous stage and the fitted curve s^mLk are weighted to obtain a weighted state estimation θ^mWkk. Finally, the method yields total distance estimation d^mk in the case of LOS/NLOS.
(4)
Location Estimation. With the distance estimation d^mk obtained from the previous stage, the coordinate estimate Z^k=x^ky^kT of the MS is calculated in this stage by the maximum likelihood estimation localization method.Figure 2
Flowchart of proposed method.


## 3.2. NLOS Detection through Polynomial Fitting

NLOS detection is implemented by quadratic polynomial fitting which is performed on the lastL historical measurement data based on the current measurement time. Then, the distance prediction of the mobile node MS is obtained by the curve fitted from previous polynomial fitting, which is illustrated in Figure 3. A residual is obtained from the prediction value and the distance measurement at the next measurement time. Finally, the residual is used as the threshold to identify whether the signal propagation condition between the beacon node BSm and the MS is LOS or NLOS.Step 1.
A curve is fitted using the lastL historical measurement data according to the BSm, and the curve can be described as s^mLk=∑n=0N−1a^mLnkn, where a^mLn is the nth coefficient of the fitted polynomial, m represents the BSm, and k is the time step. The unknown coefficients are obtained via least squares method:
(14)I=min∑k=0L−1s^mLk−dmk2=min∑k=0L−1∑n=0N−1a^mLnkn−dmk2.
I=min∑k=0L−1∑n=0N−1a^mLnkn−dmk2 is the multivariate function of a^mL0,a^mL1,⋯,a^mLN−1. Since a function has a partial derivative at a point and has an extreme value at that point, ∂I/∂a^mLα=0, for α=0,1,⋯,N−1; then, a^mL0,a^mL1,⋯,a^mLN−1 can be solved. Thus, we obtain the fitted curve:
(15)s^mLk=∑n=0N‐1a^mLnkn.Step 2.
After obtaining distance measurementdmk at time k, L measurements dmk−L⋯dmk−1 are selected to fit a curve s^mLk−1. Then, the distance between the BSm and the MS at time k can be predicted by the fitted curve. The one step distance prediction obtained by the BSm at time k is
(16)dmpk=s^mLkk−1,where s^mLkk−1 represents the value of the curve s^mLk−1 when horizontal coordinate is k. Then, we calculate the residual by
(17)empk=dmpk−dmk.
Based on the premise that the standard deviation of the distance measurement samples in NLOS is larger than that in LOS, it can be predicted whether there exists NLOS propagation between theBSm and the MS at time k:
(18)H0:empk<n⋅σm,H1:empk≥n⋅σm,2≤n≤5.
The hypothesisH0 holds true, if the propagation condition between the BSm and the MS is LOS at time k. The alternative H1 holds when the condition is NLOS. n is a positive integer.
In the simulation,n=3 is assumed. When H0 holds true, the probability of dmk in the interval dmok−3σm,dmok+3σm is Pdmk−dmok<3σm. According to (3) and (5), Pdmk−dmok<3σm is P=∫−33fxdx=1−2Φ3=0.9974. In this case, it is a small probability event when dmk is outside the interval dmok−3σm,dmok+3σm.Figure 3
Schematic diagram of polynomial fitting.


## 3.3. Adjusted Kalman Filter

In the Kalman filter section, part of the state space model is adjusted as follows.The linear term coefficient of the expression of the fitted curve can be understood as the velocity of the mobile nodeMS; that is, v^mk=a^mL1. Thus, the state vector can be rewritten as
(19)θmk=dmoka^mL1T,fork=1,2,⋯,K.Similarly, the quadratic coefficient of the fitted curve can be assumed as the acceleration of theMS; that is, a^MS=a^mL2. Thus, the state space model can be rewritten as
(20)θmk=Fθmk−1+Ba^mL2.Then, the Kalman filter is performed based on this model.Step 1.
Kalman Prediction. The state prediction and the prediction error covariance matrix can be obtained by the Kalman filter based on the rewritten state space model:(21)θ^mkk−1=Fθ^mk−1k−1,(22)Pmkk−1=FPmk−1k−1FT+σm2BBT,where θ^mkk−1 is the prediction of the state vector, θ^mk−1k−1 is the estimation of the state vector at time k−1, Pmkk−1 is the prediction error covariance matrix, and Pmk−1k−1 is the estimation error covariance matrix at time k−1. The Kalman residual is
(23)emk=zmk−z^mkk−1,(24)z^mkk−1=Hθ^mkk−1,where H=10. The Kalman gain is
(25)Kmk=Pmkk−1HTHPmkk−1HT+V−1.
It is employed to minimize the covariance of the distance prediction errorEemkemTk. V is the observation error covariance matrix.Step 2.
Kalman Update. The state estimation and the estimation error covariance matrix can be obtained by the following formula:(26)θ^mkk=θ^mkk−1+Kmkemk,(27)Pmkk=Pmkk−1−KmkHPmkk−1HT+VKmkT.
The distance estimated^mk between the beacon node BSm and the MS can be obtained by the state estimate θ^mkk:
(28)d^mk=Dθ^mkk,LOS,Dθ^mWkk,NLOS,where D=10, θ^mWkk is the weighted value of the state estimate under NLOS, and the weighting filter method will be introduced in the next section.


## 3.4. Combination through Weighting Filter

In the NLOS propagation condition, in order to alleviate the NLOS error, a weighted summation can be employed with the fitted curves^mLk and the state estimate θ^mkk. s^mLk is obtained from dmk−L+1⋯dmk, which is different from the curve s^mLk−1 calculated in NLOS detection stage.First, for further trust on the fitted value when the filtered value badly deviates from the real distance,wmk=exp−θ^mkk/s^mLkk is defined as the weight, where s^mLkk represents the value of curve s^mLk when horizontal coordinate is k. Then, the weighted sum can be expressed as
(29)θ^mWk∣k=s^mLkkv^mk1−wmk+θ^mkkwmk,where v^mk=Eθ^mkk,E=01. Thus, the distance measurement dmk, according to the beacon node BSm at time k, converts into the distance estimate d^mk which can be described as
(30)d^mk=Dθ^mkk,LOS,Ds^mLk∣kv^mkTwmk+θ^mkk1−wmk,NLOS,Algorithm 1: PF-AKF.
1:fork=1,2,⋯,Kdo
2: obtain distance measurementdmk.
3: Calculate curves^mLk−1 and s^mLk by polynomial fitting using Eq. (14) and Eq. (15).
4: Obtain distance predictiondmpk+1 using Eq. (16).
5: NLOS detection by hypothesis test using Eq. (18).
6: Adjusted Kalman Filter using Eq. (19) and Eq. (20).
7: Kalman prediction using Eq. (21) and Eq. (22).
8: Kalman update using Eq. (26) and Eq. (27).
9:if NLOS detected then
10:  calculate weighted state estimationθ^mWkk by weighting filter using Eq. (29).
11:end if
12: obtain distance estimationd^mk by combination using Eq. (28).
13: calculate coordinate estimationZ^k by maximum likelihood localization using Eq. (34).
14:end for


## 3.5. Location Estimation

In this section, maximum likelihood localization method [32, 33] is briefly introduced. As mentioned above, the coordinates of the beacon nodes are x1,y1,x2,y2,⋯,xM,yM. The coordinate of the mobile node MS at time k is Z=xk,ykT for k=1,2,⋯,K. The output of the proposed method is d^mk. Then, the following equation set is obtained:
(31)x1−xk2+y1−yk2=d^1k2,⋮xM−xk2+yM−yk2=d^Mk2.The equation set above can be rewritten as follows:(32)2xkx1−x2+2yky1−y2=d^2k2−d^1k2+x12+y12−x22+y22,⋮2xkx1−xM+2yky1−yM=d^Mk2−d^1k2+x12+y12−xM2+yM2.It can be represented by a linear equationAZ^k=B, where A and B are given by
(33)A=2x1−x2y1−y2x1−x3y1−y3⋮⋮x1−xMy1−yM,B=d^2k2−d^1k2−x22+y22+x12+y12d^3k2−d^1k2−x32+y32+x12+y12⋮d^Mk2−d^1k2−xM2+yM2+x12+y12.The coordinate matrix of theMS is calculated as follows
(34)Z^k=x^ky^kT=ATA−1ATBThe overall PF-AKF procedure is summarized in the pseudocode in Algorithm1.


## 4. Performance Evaluation and Analysis

In this section, extensive simulation experiments are carried out with the help of MATLAB2017a installed on a PC with Intel(R) Core(TM) i5-6300HQ CPU @ 2.30 GHz and 8.00 GB RAM. To prove the efficacy of the proposed PF-AKF algorithm, the simulation results of the proposed method are compared with those of the Kalman filter (KF), adjusted Kalman filter (AKF), and Kalman-based interacting multiple model (K-IMM) [15] techniques. KF is a classic method on which the proposed method is established. K-IMM, whose simulation preset is similar to the proposed method, is an algorithm which can reach high localization accuracy. The beacon nodes are randomly deployed in a square space of 120m×100m, and the propagation condition between the mobile node and the beacon nodes is randomly generated with LOS probability α. In the simulation, just one parameter changes in fixed increments every time 1000 Monte Carlo operations are conducted, and the others remain default. Root mean square error (RMSE) is considered a performance indicator, which is given by
(35)RMSE=1K⋅tn∑i=1tn∑k=1Kxk−x^ik2+yk−y^ik2,where tn=1000, K=91, xk,yk is the real location of the mobile node at time k, and x^ik,y^ik represents the location estimation for i-th trial at time k. Additionally, error bars are imparted, showing standard deviation, to evaluate the performance in a comprehensive approach.Figure4 shows the position of beacon nodes and also shows that the mobile node moves along a fixed single trajectory. In each simulation, the trajectory of the mobile node is fixed and it is assumed that the mobile node has a constant velocity of 1 m/s.Figure 4
Position of the mobile node simulated at each time step.


### 4.1. The NLOS Errors Obey a Gaussian Distribution

The default parameter values in the simulation are shown in Table1.Table 1
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
NμNLOS,σNLOS2
N2,42
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure5 shows the effect of changing the number of beacon nodes on the localization accuracy of PF-AKF, K-IMM, KF, and AKF. As the number of beacon nodes increases, the localization accuracy of the four algorithms increases. The localization accuracy of PF-AKF is always better than that of the other three, and AKF always performs better than KF. In the case of M=5, the localization accuracy of PF-AKF is 12.190%, 17.780%, and 22.300% higher than K-IMM, AKF, and KF, respectively.Figure 5
RMSE versus number of beacon nodes.Figure6 shows the effect of changing the probability of LOS propagation α on the localization accuracy. As α increases, RMSE decreases rapidly. The localization accuracy of PF-AKF is higher than K-IMM and KF. When α=0.6, the localization accuracy of PF-AKF is 11.500%, 18.020%, and 21.18% higher than that of K-IMM, AKF, and KF, respectively.Figure 6
RMSE versus probability of LOS propagationα.Figure7 shows the effect of changing the deviation of measurement noise σm on the localization accuracy. With the increase of σm, the localization accuracy of the four algorithms gradually degrades. In the case of small σm, PF-AKF performs better than the other three algorithms; however, under a large σm condition specifically σm greater than 3.5 m, the localization accuracy of K-IMM is higher than PF-AKF. When σm=1, the PF-AKF localization accuracy is 9.570%, 15.560%, and 19.960% higher than K-IMM, AKF, and KF, respectively. Meanwhile, when σm=4, the PF-AKF localization accuracy is 3.510% lower than K-IMM and 10.770% and 15.42% higher than AKF and KF, respectively. When σm reaches a certain value, performance of PF-AKF becomes poorer than K-IMM, which results from relatively low accuracy of NLOS detection.Figure 7
RMSE versus standard deviation of measurement noise.The effect of varying NLOS errors’ meanμNLOS and standard deviation σNLOS on the four algorithms is presented in Figures 8 and 9, respectively. As shown in Figure 8, as the mean of NLOS errors increases, the RMSE of the four algorithms increases. The average RMSE of PF-AKF is 2.7335 meters, and those of KF, AKF, and K-IMM algorithms are 3.0411 meters, 3.2599 meters, and 3.4464 meters, respectively. As shown in Figure 9, with the increase of the standard deviation of NLOS error, the localization accuracy of the four algorithms degrades. When σNLOS=2, the localization accuracy of PF-AKF improves about 13.480%, 17.180%, and 19.190% in comparison with K-IMM, AKF, and KF algorithms, respectively.Figure 8
RMSE versus mean of NLOS error.Figure 9
RMSE versus standard deviation of NLOS error.


### 4.2. The NLOS Errors Obey a Uniform Distribution

The default parameter values in the simulation are shown in Table2.Table 2
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
Uumin,umax
U2,5
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure10 shows the cumulative distribution function of the localization error. 97.5% of the localization error of PF-AKF, K-IMM, AKF, and KF goes below 1.960 meters, 2.067 meters, 2.758 meters, and 2.871 meters, respectively.Figure 10
CDF versus localization error.As shown in Figure11, the localization error of the four algorithms decreases as the number of nodes increases. Among them, PF-AKF has the highest localization accuracy, as the mean value of RMSE of PF-AKF is 7.390%, 9.350%, and 13.770% lower than K-IMM, KF, and AKF, respectively.Figure 11
RMSE versus number of beacon nodes.As shown in Figure12, performance of the four algorithms degrades with the increase of umax. The localization accuracy of PF-AKF is higher than K-IMM and KF algorithms. In the case of umax=6, the localization accuracy of PF-AKF is 10.520%, 13.960%, and 14.190% higher than the K-IMM, AKF, and KF algorithms, respectively. Under small umax condition, RMSE of the three methods is similar except AKF whose RMSE is much larger than the other three. Therefore, PF-AKF plays an advantage when the mean of NLOS error is large.Figure 12
RMSE versus parameterumax.


### 4.3. The NLOS Errors Obey an Exponential Distribution

The default parameter values in the simulation are shown in Table3.Table 3
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
Eλ
E4
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure13 shows the cumulative distribution function of the localization error in the case where the NLOS noise obeys the exponential distribution. 90% of localization error of the PF-AKF, AKF, K-IMM, and KF is less than 1.647 m, 1.744 m, 1.879 m, and 2.348 m, respectively.Figure 13
CDF versus localization error.The variation of the RMSE of the four algorithms when changingλ is shown in Figure 14. As λ increases, the RMSE of the four algorithms increases. The localization error of PF-AKF and K-IMM is rather close. Moreover, in most cases, the localization error of both PF-AKF and K-IMM is smaller than the KF algorithm. When λ is less than 7.04, PF-AKF has the highest localization accuracy. However, when λ>7.04, PF-AKF performs poorer than K-IMM.Figure 14
RMSE versus parameterλ.


### 4.4. Discussion on Performance

In most cases, the proposed method shows superior performance compared to KF, AKF, and K-IMM. However, there exists some cases such as large standard deviation of measurement noise and small parameter value in exponential distribution where the proposed method performs badly, even inferior to K-IMM. Therefore, in general, the proposed method possesses higher localization accuracy than KF, AKF, and K-IMM, particularly under large NLOS noise cases. This is because PF-AKF’s NLOS identification could better handle large NLOS noise than KF, AKF, and K-IMM. There exists downturn of the proposed method due to low accuracy of polynomial fitting caused by large NLOS noise in the first few measurement times. Therefore, a large standard deviation of measurement noise might be a breakthrough for further study.


### 4.5. Experimental Results

To further testify the performance of PF-AKF, an experiment is implemented in a real environment. In order to obtain experimental data, a chirp spread spectrum (CSS) positioning system is employed. As shown in Figure15, 7 beacon nodes are deployed in a 9.9m×16m room. The beacon node and the mobile node are about 1.1 m above the ground. The distance measurement frequency of CSS nodes is set to 20 Hz. The mobile node moves along a U-shape trajectory at a uniform speed around the obstacle. In order to compare the proposed PF-AKF algorithm with the KF, AKF, and K-IMM algorithms, the same set of measured data is applied to all of them with the help of MATLAB2017a.Figure 15
The practical environment.In Figure16, the positioning errors of the PF-AKF, KF, AKF, and K-IMM algorithms at each sample point are compared. The average localization errors of PF-AKF, KF, AKF, and K-IMM algorithms are 0.5240 m, 0.6564 m, 0.6170 m, and 0.6343 m, respectively. The localization accuracy of PF-AKF algorithm is 20.17%, 15.07%, and 17.39% higher than that of KF, AKF, and K-IMM, respectively.Figure 16
Localization error of each sample point.Figure17 shows the cumulative distribution function of the localization error of the PF-AKF, KF, AKF, and K-IMM algorithms. 90% localization error of the PF-AKF, KF, AKF, and K-IMM algorithms is less than 0.6807 m, 0.8494 m, 0.7701 m, and 0.7641 m, respectively. When the localization error of the PF-AKF, KF, AKF, and K-IMM algorithms reaches 0.9162 m, 1.1390 m, 1.1310 m, and 1.062 m, respectively, the cumulative distribution function is close to 1.Figure 17
CDF versus localization error.


### 4.6. Computational Time

As shown in Table4, runtimes of the four algorithms are compared. The four algorithms are programmed under aforementioned configuration. The time duration of a single mobile target’s one (time) step localization is calculated and compared to that of KF, AKF, and K-IMM. The proposed method runs faster than the K-IMM algorithm. Since fitting curves consume some time, the proposed method is relatively more complex than the KF algorithm.Table 4
Running times of computation time.
Method used
Running times (s)
KF
0.000638
AKF
0.005433
K-IMM
0.018271
PF-AKF
0.006633


## 4.1. The NLOS Errors Obey a Gaussian Distribution

The default parameter values in the simulation are shown in Table1.Table 1
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
NμNLOS,σNLOS2
N2,42
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure5 shows the effect of changing the number of beacon nodes on the localization accuracy of PF-AKF, K-IMM, KF, and AKF. As the number of beacon nodes increases, the localization accuracy of the four algorithms increases. The localization accuracy of PF-AKF is always better than that of the other three, and AKF always performs better than KF. In the case of M=5, the localization accuracy of PF-AKF is 12.190%, 17.780%, and 22.300% higher than K-IMM, AKF, and KF, respectively.Figure 5
RMSE versus number of beacon nodes.Figure6 shows the effect of changing the probability of LOS propagation α on the localization accuracy. As α increases, RMSE decreases rapidly. The localization accuracy of PF-AKF is higher than K-IMM and KF. When α=0.6, the localization accuracy of PF-AKF is 11.500%, 18.020%, and 21.18% higher than that of K-IMM, AKF, and KF, respectively.Figure 6
RMSE versus probability of LOS propagationα.Figure7 shows the effect of changing the deviation of measurement noise σm on the localization accuracy. With the increase of σm, the localization accuracy of the four algorithms gradually degrades. In the case of small σm, PF-AKF performs better than the other three algorithms; however, under a large σm condition specifically σm greater than 3.5 m, the localization accuracy of K-IMM is higher than PF-AKF. When σm=1, the PF-AKF localization accuracy is 9.570%, 15.560%, and 19.960% higher than K-IMM, AKF, and KF, respectively. Meanwhile, when σm=4, the PF-AKF localization accuracy is 3.510% lower than K-IMM and 10.770% and 15.42% higher than AKF and KF, respectively. When σm reaches a certain value, performance of PF-AKF becomes poorer than K-IMM, which results from relatively low accuracy of NLOS detection.Figure 7
RMSE versus standard deviation of measurement noise.The effect of varying NLOS errors’ meanμNLOS and standard deviation σNLOS on the four algorithms is presented in Figures 8 and 9, respectively. As shown in Figure 8, as the mean of NLOS errors increases, the RMSE of the four algorithms increases. The average RMSE of PF-AKF is 2.7335 meters, and those of KF, AKF, and K-IMM algorithms are 3.0411 meters, 3.2599 meters, and 3.4464 meters, respectively. As shown in Figure 9, with the increase of the standard deviation of NLOS error, the localization accuracy of the four algorithms degrades. When σNLOS=2, the localization accuracy of PF-AKF improves about 13.480%, 17.180%, and 19.190% in comparison with K-IMM, AKF, and KF algorithms, respectively.Figure 8
RMSE versus mean of NLOS error.Figure 9
RMSE versus standard deviation of NLOS error.


## 4.2. The NLOS Errors Obey a Uniform Distribution

The default parameter values in the simulation are shown in Table2.Table 2
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
Uumin,umax
U2,5
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure10 shows the cumulative distribution function of the localization error. 97.5% of the localization error of PF-AKF, K-IMM, AKF, and KF goes below 1.960 meters, 2.067 meters, 2.758 meters, and 2.871 meters, respectively.Figure 10
CDF versus localization error.As shown in Figure11, the localization error of the four algorithms decreases as the number of nodes increases. Among them, PF-AKF has the highest localization accuracy, as the mean value of RMSE of PF-AKF is 7.390%, 9.350%, and 13.770% lower than K-IMM, KF, and AKF, respectively.Figure 11
RMSE versus number of beacon nodes.As shown in Figure12, performance of the four algorithms degrades with the increase of umax. The localization accuracy of PF-AKF is higher than K-IMM and KF algorithms. In the case of umax=6, the localization accuracy of PF-AKF is 10.520%, 13.960%, and 14.190% higher than the K-IMM, AKF, and KF algorithms, respectively. Under small umax condition, RMSE of the three methods is similar except AKF whose RMSE is much larger than the other three. Therefore, PF-AKF plays an advantage when the mean of NLOS error is large.Figure 12
RMSE versus parameterumax.


## 4.3. The NLOS Errors Obey an Exponential Distribution

The default parameter values in the simulation are shown in Table3.Table 3
Default parameter values.
Parameters
Symbol
Default values
Number of beacon nodes
M
6
Probability of LOS propagation
α
0.6
Standard deviation of measurement noise (m)
σm
1
NLOS error
Eλ
E4
Number of sample points
T
100
Number of Monte Carlo runs
tn
1000Figure13 shows the cumulative distribution function of the localization error in the case where the NLOS noise obeys the exponential distribution. 90% of localization error of the PF-AKF, AKF, K-IMM, and KF is less than 1.647 m, 1.744 m, 1.879 m, and 2.348 m, respectively.Figure 13
CDF versus localization error.The variation of the RMSE of the four algorithms when changingλ is shown in Figure 14. As λ increases, the RMSE of the four algorithms increases. The localization error of PF-AKF and K-IMM is rather close. Moreover, in most cases, the localization error of both PF-AKF and K-IMM is smaller than the KF algorithm. When λ is less than 7.04, PF-AKF has the highest localization accuracy. However, when λ>7.04, PF-AKF performs poorer than K-IMM.Figure 14
RMSE versus parameterλ.


## 4.4. Discussion on Performance

In most cases, the proposed method shows superior performance compared to KF, AKF, and K-IMM. However, there exists some cases such as large standard deviation of measurement noise and small parameter value in exponential distribution where the proposed method performs badly, even inferior to K-IMM. Therefore, in general, the proposed method possesses higher localization accuracy than KF, AKF, and K-IMM, particularly under large NLOS noise cases. This is because PF-AKF’s NLOS identification could better handle large NLOS noise than KF, AKF, and K-IMM. There exists downturn of the proposed method due to low accuracy of polynomial fitting caused by large NLOS noise in the first few measurement times. Therefore, a large standard deviation of measurement noise might be a breakthrough for further study.


## 4.5. Experimental Results

To further testify the performance of PF-AKF, an experiment is implemented in a real environment. In order to obtain experimental data, a chirp spread spectrum (CSS) positioning system is employed. As shown in Figure15, 7 beacon nodes are deployed in a 9.9m×16m room. The beacon node and the mobile node are about 1.1 m above the ground. The distance measurement frequency of CSS nodes is set to 20 Hz. The mobile node moves along a U-shape trajectory at a uniform speed around the obstacle. In order to compare the proposed PF-AKF algorithm with the KF, AKF, and K-IMM algorithms, the same set of measured data is applied to all of them with the help of MATLAB2017a.Figure 15
The practical environment.In Figure16, the positioning errors of the PF-AKF, KF, AKF, and K-IMM algorithms at each sample point are compared. The average localization errors of PF-AKF, KF, AKF, and K-IMM algorithms are 0.5240 m, 0.6564 m, 0.6170 m, and 0.6343 m, respectively. The localization accuracy of PF-AKF algorithm is 20.17%, 15.07%, and 17.39% higher than that of KF, AKF, and K-IMM, respectively.Figure 16
Localization error of each sample point.Figure17 shows the cumulative distribution function of the localization error of the PF-AKF, KF, AKF, and K-IMM algorithms. 90% localization error of the PF-AKF, KF, AKF, and K-IMM algorithms is less than 0.6807 m, 0.8494 m, 0.7701 m, and 0.7641 m, respectively. When the localization error of the PF-AKF, KF, AKF, and K-IMM algorithms reaches 0.9162 m, 1.1390 m, 1.1310 m, and 1.062 m, respectively, the cumulative distribution function is close to 1.Figure 17
CDF versus localization error.


## 4.6. Computational Time

As shown in Table4, runtimes of the four algorithms are compared. The four algorithms are programmed under aforementioned configuration. The time duration of a single mobile target’s one (time) step localization is calculated and compared to that of KF, AKF, and K-IMM. The proposed method runs faster than the K-IMM algorithm. Since fitting curves consume some time, the proposed method is relatively more complex than the KF algorithm.Table 4
Running times of computation time.
Method used
Running times (s)
KF
0.000638
AKF
0.005433
K-IMM
0.018271
PF-AKF
0.006633


## 5. Conclusions

In this paper, a polynomial fitting-based adjusted Kalman filter (PF-AKF) method is proposed, which is robust to NLOS errors. The method does not require any prior statistical information of the NLOS noise. First, historical measurement data obtained by the beacon nodes are used to fit the polynomial curve by the least squares method and predict the measured value at the next moment by the curve. Then, the residual is obtained by subtracting the observed value from the predicted value. The residual is then employed to judge whether the propagation condition is LOS or NLOS. Kalman filtering is performed after approximating the primary coefficient and the quadratic coefficient obtained by the polynomial fitting to the velocity and acceleration of the mobile node into the state vector. In the presence of NLOS error, the KF results are weighted with the fitted curves to attenuate the effect of NLOS errors. The simulation results show that the method is robust to the NLOS errors with higher localization accuracy, and the method has better performance when the measurement noise is small. For future work, the computational complexity would be reduced through optimization of the algorithm, and we will also apply the method to TOA/RSS data fusion-based localization with unscented Kalman filter to enhance the performance.



---

*Source: 6787252-2020-01-07.xml*