Graph Optimization Model Fusing BLE Ranging with Wi-Fi Fingerprint for Indoor Positioning

To improve the user’s positioning accuracy of a Wi-Fi fingerprint-based positioning algorithm, this study proposes a graph optimization model based on the framework of g2o that fuses a Wi-Fi fingerprint and Bluetooth Low Energy (BLE) ranging technologies. In our model, the improvement in positioning can be formulated as a nonlinear least-squares optimization problem that a graph can represent. The graph regards users as nodes and our self-designed error functions between users as edges. In the graph, the nodes obtain the initial coordinates through Wi-Fi fingerprint positioning, and all error functions aggregate to a total error function to be solved. To improve the solution effect of the total error function and weaken the influence of measurement error, an information matrix, an edge selection principle, and a Huber kernel function are introduced. The Levenberg–Marquardt (LM) algorithm is used to solve the total error function and the affine transformation estimation is used for the drifting solution. Through experiments, the influence of the threshold in the Huber kernel function is explored, the relationship between the number of nodes in the graph and the optimization effect is analyzed, and the impact of the distribution of nodes is researched. The experimental results show improvements in the positioning accuracy of four common Wi-Fi fingerprint-matching algorithms: KNN, WKNN, GK, and Stg.


Introduction
With increasing demand for location-based services (LBSs), indoor positioning has become a research topic of much interest [1]. In recent years, some RF-based methods for precise indoor positioning have emerged, mainly in two categories: trilateration [2] and fingerprint [3]. The trilateration method requires the physical distance between the target point and at least three anchor points. The measurements used to calculate physical distance in trilateration include received signal strength indication (RSSI) [4], time difference of arrival (TDOA) [5], and time of arrival (TOA) [6]. Compared to TDOA and TOA, RSSI-based ranging avoids the need for expensive clock synchronization systems and additional hardware, and has the advantages of low cost, easy deployment, and simple calculation. The fingerprint method is regarded as superior to trilateration in more complex indoor environments, where the interference, reflection, and refraction of signal forces are superimposed onto the specific points [7].
Wi-Fi technology is widely used in indoor target intrusion detection and localization [8,9]. The Wi-Fi wireless access points (APs) have been widely deployed in indoor environments because of the improvement of wireless sensor technology and the growing demand for network communication, and they provide a physical device foundation for the application of Wi-Fi fingerprint positioning. The Wi-Fi fingerprint method is regarded as low-cost and simple because it requires no additional hardware and performs satisfactorily in complex indoor environments [10]. However, room exists for improvement in its positioning accuracy in complex indoor environments such as underground parking lots.
Neural network and deep-learning technologies are also applied in Wi-Fi finger positioning [21,22], but regardless of the system or algorithm, the Wi-Fi fingerp method always follows a general process consisting of offline and online phases, as sh in Figure 1. In the offline phase, a certain number of APs are deployed, several reference p (RPs) are set up evenly in the indoor environment, and the RSSIs from APs at the RP collected to establish a fingerprint database. A fingerprint consists of the physical co nates and Wi-Fi signal fingerprints of the RPs. Suppose the number of APs involve positioning is , the number of RPs is , and , represents the RSSI from ceived at . Since the RSSI is easily affected by environmental changes, multiple collections are carried out at each RP to build a fingerprint database; so, the kth finger at can be expressed as In the online phase, we estimate the positions of users by collecting the RSSIs o ferent APs in real-time and applying the fingerprint-matching algorithm.

Signal Propagation of BLE
Common propagation path-loss models include the free space propagation and arithmic distance path-loss models [23]. The latter is generally used in BLE ranging nology [24,25] and is expressed as where is the RSSI when the distance between the Bluetooth transmitting point receiving point is , is the reference distance, is the RSSI when the dist between the transmitting and receiving point is , is Gaussian random noise, an is the path loss index. For convenience of calculation, usually, = 1 m. To fit the si propagation model based on the collected RSSI, the simplified model can be expresse RSSI = − 10 lg , In the offline phase, a certain number of APs are deployed, several reference points (RPs) are set up evenly in the indoor environment, and the RSSIs from APs at the RPs are collected to establish a fingerprint database. A fingerprint consists of the physical coordinates and Wi-Fi signal fingerprints of the RPs. Suppose the number of APs involved in positioning is n, the number of RPs is m, and RSSI (i,j) represents the RSSI from AP i received at RP j . Since the RSSI is easily affected by environmental changes, multiple RSSI collections are carried out at each RP to build a fingerprint database; so, the kth fingerprint at RP j can be expressed as In the online phase, we estimate the positions of users by collecting the RSSIs of different APs in real-time and applying the fingerprint-matching algorithm.

Signal Propagation of BLE
Common propagation path-loss models include the free space propagation and logarithmic distance path-loss models [23]. The latter is generally used in BLE ranging technology [24,25] and is expressed as where P t (d) is the RSSI when the distance between the Bluetooth transmitting point and receiving point is d, d 0 is the reference distance, P l (d 0 ) is the RSSI when the distance between the transmitting and receiving point is d 0 , ε is Gaussian random noise, and λ is the path loss index. For convenience of calculation, usually, d 0 = 1 m. To fit the signal propagation model based on the collected RSSI, the simplified model can be expressed as where RSSI = P t (d) and A = P l (1). Both A and λ are empirical values that are related to the indoor environment. By collecting the RSSI of BLE beacons at different distances and fitting the parameters by least-squares, Equation (4) can be used for BLE ranging.

Graph Optimization Based on g 2 o
Many robotics and computer vision studies involve the minimization of nonlinear error functions, a typical example being simultaneous localization and mapping (SLAM) [26][27][28]. These studies aim to find parameters or state variables that best explain a set of measured values affected by Gaussian noise, where the LM algorithm is used to find the optimal solution [29]. g 2 o is a general back-end optimization framework suitable for nonlinear least-squares problems [30], which can be represented graphically in the framework of g 2 o, as shown in Figure 2, where each node represents a state variable to be optimized, and an edge between nodes represents a paired observation of the two nodes.
where RSSI = and = 1 . Both and are empirical values that to the indoor environment. By collecting the RSSI of BLE beacons at differen and fitting the parameters by least-squares, Equation (4) can be used for BLE

Graph Optimization Based on 2 o
Many robotics and computer vision studies involve the minimization o error functions, a typical example being simultaneous localization and mapp [26][27][28]. These studies aim to find parameters or state variables that best exp measured values affected by Gaussian noise, where the LM algorithm is used optimal solution [29]. g 2 o is a general back-end optimization framework suitable for nonl squares problems [30], which can be represented graphically in the framewor shown in Figure 2, where each node represents a state variable to be optimiz edge between nodes represents a paired observation of the two nodes. Graph optimization based on g 2 o represents a nonlinear least-squares p graph and visually displays the structural relationship between the variables can be directed or undirected, and the corresponding optimization problem scribed as = , Ω , , + , Ω , , + , Ω , , + , Ω , , + , Ω , , + ⋯ where is the final total error function to be solved, , is the error f tween nodes and , and Ω , is the information matrix of edge , . ranging data to optimize Wi-Fi fingerprint location results can be abstracted optimization problem. Users can be defined as nodes in the graph, and the d edges between nodes is introduced in Section 3.2.

Model Overview
Wi-Fi fingerprint positioning can provide position results at the current m are not dependent on the previous one, so there is no cumulative error as the u However, coordinates obtained in this way are not accurate enough due to n timize the positioning results of a Wi-Fi fingerprint, physical distances bet measured by BLE ranging technology are added to establish redundant constr Graph optimization based on g 2 o represents a nonlinear least-squares problem as a graph and visually displays the structural relationship between the variables. The graph can be directed or undirected, and the corresponding optimization problem can be described as 3) e (2,3) + e T (2,4) Ω (2,4) e (2,4) + e T (3,5) Ω (3,5) e (3,5) + e T (4,5) Ω (4,5) e (4,5) + . . . , where F(x) is the final total error function to be solved, e (i,j) is the error function between nodes x i and x j , and Ω (i,j) is the information matrix of edge e (i,j) . Using BLE ranging data to optimize Wi-Fi fingerprint location results can be abstracted as a graph optimization problem. Users can be defined as nodes in the graph, and the definition of edges between nodes is introduced in Section 3.2.

Model Overview
Wi-Fi fingerprint positioning can provide position results at the current moment that are not dependent on the previous one, so there is no cumulative error as the user moves. However, coordinates obtained in this way are not accurate enough due to noise. To optimize the positioning results of a Wi-Fi fingerprint, physical distances between users measured by BLE ranging technology are added to establish redundant constraints. Then, graph optimization based on g 2 o can be adopted to adjust the coordinates of users in the graph and improve positioning accuracy. The construction of the model is shown in Figure 3.  We define the users as nodes of the graph and obtain their initial positioning coordinates by the construction of an offline Wi-Fi fingerprint database and online fingerprint matching. Users receive BLE RSSIs from other users' devices and input them to the BLE ranging model to obtain physical distances to other users. We define the error function and regard it as the edges between pairs of nodes. The total error function is constructed by all edges and grows in size with the number of users. The optimal solution of the total error function is obtained by the solver, i.e., the LM algorithm and the drifting solution introduced in Section 3.3.
There will inevitably be several edges with large relative errors in the graph that affect optimization. The robust kernel function, the method of reconstructing the information matrix introduced in Section 3.2, and the filtering conditions of the edges determined by the experiments in Section 4.2 are introduced to weaken the effects of these edges.
When we experimented with the proposed model, we used a laptop (Lenovo, Beijing, China) running Windows 10 as the central server for user positioning calculation. The users need to provide the RSSIs received from APs to the central server, and the initial coordinates of the users can be obtained based on the Wi-Fi fingerprint positioning. The RSSIs received from other users' devices with BLE function also need to be input into server to run the graph optimization model, and users can obtain the optimized position from the server.

Construction and Solution of Total Error Function
Constructing a reasonable total error function is the key to solving the graph optimization problem. As the number of user nodes in the graph increases, the number of variables and in the total error function expand twice as fast as nodes. We find that the number of nodes grows linearly, while the number of edges grows exponentially in the graph as the number of users joining the positioning increases. Because the total error function can easily fall into a locally optimal solution during iteration, assigning each node a good initial value can effectively solve the above problems, that is, obtaining the user's initial coordinates through the Wi-Fi fingerprint positioning. The graph, nodes, and edges are defined as follows.
Suppose there are user nodes in the graph, user node is defined as = , , 1 ≤ ≤ , and its initial value as obtained by Wi-Fi fingerprint positioning is denoted as = , , 1 ≤ ≤ . The distance obtained by the BLE ranging model between nodes and is defined as , , 1 ≤ , ≤ , ≠ , and the distance between nodes and in the graph is defined as We define the users as nodes of the graph and obtain their initial positioning coordinates by the construction of an offline Wi-Fi fingerprint database and online fingerprint matching. Users receive BLE RSSIs from other users' devices and input them to the BLE ranging model to obtain physical distances to other users. We define the error function and regard it as the edges between pairs of nodes. The total error function is constructed by all edges and grows in size with the number of users. The optimal solution of the total error function is obtained by the solver, i.e., the LM algorithm and the drifting solution introduced in Section 3.3.
There will inevitably be several edges with large relative errors in the graph that affect optimization. The robust kernel function, the method of reconstructing the information matrix introduced in Section 3.2, and the filtering conditions of the edges determined by the experiments in Section 4.2 are introduced to weaken the effects of these edges.
When we experimented with the proposed model, we used a laptop (Lenovo, Beijing, China) running Windows 10 as the central server for user positioning calculation. The users need to provide the RSSIs received from APs to the central server, and the initial coordinates of the users can be obtained based on the Wi-Fi fingerprint positioning. The RSSIs received from other users' devices with BLE function also need to be input into server to run the graph optimization model, and users can obtain the optimized position from the server.

Construction and Solution of Total Error Function
Constructing a reasonable total error function is the key to solving the graph optimization problem. As the number of user nodes in the graph increases, the number of variables (x i and y i ) in the total error function expand twice as fast as nodes. We find that the number of nodes grows linearly, while the number of edges grows exponentially in the graph as the number of users joining the positioning increases. Because the total error function can easily fall into a locally optimal solution during iteration, assigning each node a good initial value can effectively solve the above problems, that is, obtaining the user's initial coordinates through the Wi-Fi fingerprint positioning. The graph, nodes, and edges are defined as follows.
Suppose there are n user nodes in the graph, user node i is defined as , and its initial value as obtained by Wi-Fi fingerprint positioning is denoted as The distance obtained by the BLE ranging model between nodes i and j is defined as d b (i,j) , (1 ≤ i, j ≤ n, i = j), and the distance between nodes i and j in the graph is defined as Since the distance between users measured by the BLE ranging model adds constraints on nodes, in our proposed model, d (i,j) will be close to d b (i,j) by adjusting the locations of nodes. So, the edge between p i and p j is defined as The total error function of the entire graph is composed of all of the error functions (called edges) between pairs of nodes. It can be defined as where Ω (i,j) is the information matrix of e p i , p j , and p = {p 1 , p 2 , . . . , p n } and p * = p * 1 , p * 2 , . . . , p * n are respectively a set of user nodes before and after optimization.
is the final position of p i after optimization, because p i keeps changing during the iterative solution of the total error function.
After defining the nodes and edges, we can define the size of the information matrix as 1 × 1, which can be regarded as a weighting coefficient. The distance error between nodes based on the Wi-Fi fingerprint positioning cannot be predicted accurately; so, we do not consider adding the statistical indicators of Wi-Fi fingerprint positioning to calculate the information matrix.
Previous experiments with the BLE ranging model have shown that the RSSI of BLE can be approximated to a Gaussian distribution for simplicity [31]. As the distance increases, the BLE RSSI fluctuation increases, which directly leads to an increase in the fluctuation of the ranging error. The relationship between true distance and BLE ranging distance can be defined as where d r (i,j) is the real distance between p i and p j , and l (i,j) is the error between d b (i,j) and d r (i,j) . We use d b (i,j) as a substitute for d r (i,j) because d r (i,j) cannot be predicted. The experiments of Section 4.2 show that the distance error of BLE ranging is not only related to the physical distance, but is also closely related to the standard deviation (SD) of the ranging results. Thus, the SD of BLE ranging is added to calculate the information matrix: For simplicity of notation, in the rest of this paper, we will express the error function in a new form: The error function can be approximated by its first-order Taylor expansion around the current initial guessp: where J (i,j) is the Jacobian of e (i,j) (p), computed from the current estimatep. Since there are errors in both BLE ranging and Wi-Fi fingerprint positioning, the graph will contain wrong edges in many cases. Nonlinear least-squares optimization is greatly affected by errors; thus, it is difficult to guarantee the convergence and performance of the total error function in the iterative process, even if the information matrix has been introduced to weaken the effect of errors. To reduce the occurrence of the above situation, the Huber kernel function [32] is introduced to adjust the edges.
When the error is small, the error function is approximated as a quadratic function with faster growth. When the error is large, it is approximated as a linear function with slower growth, which still ensures that the function is continuous and differentiable. The Huber kernel function is and the total error function with the Huber kernel function can be rewritten as To preserve the characteristics of the original total error function and the calculated parameters, we can simplify the processing of the Huber kernel function by adjusting the information matrix, where Substituting the updated information matrix in the original total error function, we obtain a new total error function with the original form: Substituting Equation (13) in Equation (19), F H (p + ∆p) can be approximated as For simplicity of notation, we set can be minimized in ∆p by solving the linear system based on the LM algorithm, where H is the information matrix of this system, I is an identity matrix, and α is a damping factor. This is useful to control the step size in case of nonlinear surfaces and to control the convergence. The solution is obtained by adding the increments ∆p * to the initial guess: The popular LM algorithm iterates the linearization in Equation (20), the solution in Equation (21), and the update step in Equation (22). In every iteration, the previous solution is used as the linearization point and the initial guess until a given termination criterion is met. The derivation of the LM algorithm is not the work of this paper and is not repeated here. The implementation process of the LM algorithm is shown in Algorithm 1. solution is used as the linearization point and the initial guess until a given termination criterion is met. The derivation of the LM algorithm is not the work of this paper and is not repeated here. The implementation process of the LM algorithm is shown in Algorithm 1.
Algorithm 1: Pseudo code of LM algorithm

Drifting Solution Based on Affine Transformation Estimation
The optimized graph has a better structure, but the nodes may drift together because of the missing absolute coordinate constraints. Since the affine transformation estimation can relatively better preserve the optimized structure, we use it to solve the drifting problem. We define the nodes calculated by the LM algorithm as = , , … and the nodes estimated by affine transformation as = , , … . We redefine the coordinates of each node as = , , 1 , 1 ≤ ≤ , = , , 1 , 1 ≤ ≤ and = , , 1 , 1 ≤ ≤ . The transformation matrix is calculated by where represents the transformation matrix. The optimal estimates of parameters , , , , and in are calculated by the least-squares method. Our final optimization results of the nodes are

BLE Signal Filter and Ranging Model
During actual signal collection, the received RSSI from the BLE beacon fluctuates in the indoor environment. BLE ranging can affect the accuracy of the edges in the graph thus, BLE ranging errors should be controlled within a certain range to build the network structure of nodes in the graph. In our experiments, the transmit power of the BLE beacon was +4 dBm, and its broadcast period was set to 80 ms. Figure 4a shows the RSSI data collected at 3, 6, and 9 m between the BLE beacon and the signal receiving point, and it can be observed that the original RSSI fluctuates greatly.

Drifting Solution Based on Affine Transformation Estimation
The optimized graph has a better structure, but the nodes may drift together because of the missing absolute coordinate constraints. Since the affine transformation estimation can relatively better preserve the optimized structure, we use it to solve the drifting problem. We define the nodes calculated by the LM algorithm as p l = p l 1 , p l 2 , . . . p l n and the nodes estimated by affine transformation as p t = p t 1 , p t 2 , . . . p t n . We redefine the coordinates of each node as p l The transformation matrix is calculated by where M represents the transformation matrix. The optimal estimates of parameters a, b, c, d, e, and f in M are calculated by the least-squares method. Our final optimization results of the nodes are p t = p l M.

BLE Signal Filter and Ranging Model
During actual signal collection, the received RSSI from the BLE beacon fluctuates in the indoor environment. BLE ranging can affect the accuracy of the edges in the graph; thus, BLE ranging errors should be controlled within a certain range to build the network structure of nodes in the graph. In our experiments, the transmit power of the BLE beacon was +4 dBm, and its broadcast period was set to 80 ms. Figure 4a shows the RSSI data Since the broadcast period of the BLE beacon was set to 80 ms, the average number of RSSIs that could be collected every second was 12.5. To balance the positioning time and a reasonable filter window size, we adopted mean filtering and set the filter window to 10, that is, 10 raw BLE RSSIs were processed each time before input to the BLE ranging model. Obviously, this strategy can effectively reduce the signal fluctuation, as shown in Figure 4b.   Table 1 presents the comparison of the SD of RSSIs collected at 3, 6, and 9 m before and after filtering. The SD can be reduced by nearly 50% after filtering. This paper uses the simplified BLE ranging model (Equation (4)) and estimates the values of and by the filtered RSSI at different distances (0.25 m interval) and leastsquare curve fitting, as shown in Figure 5.
Then, we can calculate the distance between user nodes based on this mode, RSSI = −45.688 − 20.835 lg .
Through the linear fitting to the SD of the RSSI at different distances, we can obtain the signal propagation model applicable to our experiments,  Since the broadcast period of the BLE beacon was set to 80 ms, the average number of RSSIs that could be collected every second was 12.5. To balance the positioning time and a reasonable filter window size, we adopted mean filtering and set the filter window to 10, that is, 10 raw BLE RSSIs were processed each time before input to the BLE ranging model. Obviously, this strategy can effectively reduce the signal fluctuation, as shown in Figure 4b. Table 1 presents the comparison of the SD of RSSIs collected at 3, 6, and 9 m before and after filtering. The SD can be reduced by nearly 50% after filtering. This paper uses the simplified BLE ranging model (Equation (4)) and estimates the values of A and n by the filtered RSSI at different distances (0.25 m interval) and least-square curve fitting, as shown in Figure 5.
Then, we can calculate the distance between user nodes based on this mode, Through the linear fitting to the SD of the RSSI at different distances, we can obtain the signal propagation model applicable to our experiments,

Experimental Environment
Our experiments used the underground parking lot of North China Electric Power University as a testbed, where the Wi-Fi fingerprint database was constructed. The testbed was 58 m long and 42 m wide, and Wi-Fi APs were installed on the wall 2 m high relative to the parking lot surface. Figure 6a is a plan view of the testbed, where black dots mark the 10 AP locations. Test points (TPs) represent users, and RPs and TPs were set in the middle of the passage. The distance between two adjacent RPs was 2 m, and adjacent TPs were 1 m apart. We deployed 86 RPs and 84 TPs, where an "x" represents an RP and a four-color "o" represents a TP. Figure 6b shows the coordinates of RPs and TPs. Our device was set 1 m high relative to the parking lot surface, collecting Wi-Fi data at each RP. We collected 10 sets of Wi-Fi data of 10 APs continuously at each RP to construct the Wi-Fi fingerprint database, and we used the Wi-Fi fingerprint-matching algorithm to obtain the initial location of each TP. To verify the validity of the proposed model

Experimental Environment
Our experiments used the underground parking lot of North China Electric Power University as a testbed, where the Wi-Fi fingerprint database was constructed. The testbed was 58 m long and 42 m wide, and Wi-Fi APs were installed on the wall 2 m high relative to the parking lot surface. Figure 6a is a plan view of the testbed, where black dots mark the 10 AP locations. Test points (TPs) represent users, and RPs and TPs were set in the middle of the passage. The distance between two adjacent RPs was 2 m, and adjacent TPs were 1 m apart. We deployed 86 RPs and 84 TPs, where an "x" represents an RP and a four-color "o" represents a TP. Figure 6b shows the coordinates of RPs and TPs.

Experimental Environment
Our experiments used the underground parking lot of North China Electric Power University as a testbed, where the Wi-Fi fingerprint database was constructed. The testbed was 58 m long and 42 m wide, and Wi-Fi APs were installed on the wall 2 m high relative to the parking lot surface. Figure 6a is a plan view of the testbed, where black dots mark the 10 AP locations. Test points (TPs) represent users, and RPs and TPs were set in the middle of the passage. The distance between two adjacent RPs was 2 m, and adjacent TPs were 1 m apart. We deployed 86 RPs and 84 TPs, where an "x" represents an RP and a four-color "o" represents a TP. Figure 6b shows the coordinates of RPs and TPs. Our device was set 1 m high relative to the parking lot surface, collecting Wi-Fi data at each RP. We collected 10 sets of Wi-Fi data of 10 APs continuously at each RP to construct the Wi-Fi fingerprint database, and we used the Wi-Fi fingerprint-matching algorithm to obtain the initial location of each TP. To verify the validity of the proposed model Our device was set 1 m high relative to the parking lot surface, collecting Wi-Fi data at each RP. We collected 10 sets of Wi-Fi data of 10 APs continuously at each RP to construct the Wi-Fi fingerprint database, and we used the Wi-Fi fingerprint-matching algorithm to obtain the initial location of each TP. To verify the validity of the proposed model and control the cost, we conducted experiments with simulated BLE ranging data and randomly selected points from 84 TPs to perform the graph optimization.

Edge Selection
Edges with larger errors in the graph are often formed by two user nodes that are far apart. Giving these edges a reasonable information matrix can weaken the influence of errors. However, when we calculate the information matrix Ω (i,j) (Equation (11)), we need to determine the value of σ (i,j) . We collected the BLE RSSIs every 0.5 m within a range of 20 m from the BLE beacon and input them to the BLE ranging model (Equation (25)) to calculate the error and SD of ranging in each collection point, as shown in Figure 7.
Sensors 2022, 22, x FOR PEER REVIEW 11 and control the cost, we conducted experiments with simulated BLE ranging data randomly selected points from 84 TPs to perform the graph optimization.

Edge Selection
Edges with larger errors in the graph are often formed by two user nodes that ar apart. Giving these edges a reasonable information matrix can weaken the influenc errors. However, when we calculate the information matrix Ω , (Equation (11)) need to determine the value of , . We collected the BLE RSSIs every 0.5 m with range of 20 m from the BLE beacon and input them to the BLE ranging model (Equa (25)) to calculate the error and SD of ranging in each collection point, as shown in Figur   Figure 7. Error and standard deviation (SD) of BLE ranging at different distances.
As the distance increases, the error and SD of ranging generally show a linear ward trend. The greater the distance, the more frequent the fluctuations of the error SD. Experiments show that most BLE ranging errors within 15 m can be stably con trated within 5 m, and more are concentrated within 2 m. Table 2 shows the ranging e in different distances within 15 m at 1 m intervals. In the experiment, we collected R at different distances and calculated the physical distance between the transmitting p and the receiving point, but the results show the discontinuity of the error. This is ma due to the fact that the influence of the multipath effect on the ranging error is more nificant than that of signal attenuation in the range of 20 m. Considering the relation between ranging error and distance, we decided to only add edges within 15 m to graph.  As the distance increases, the error and SD of ranging generally show a linear upward trend. The greater the distance, the more frequent the fluctuations of the error and SD. Experiments show that most BLE ranging errors within 15 m can be stably concentrated within 5 m, and more are concentrated within 2 m. Table 2 shows the ranging error in different distances within 15 m at 1 m intervals. In the experiment, we collected RSSIs at different distances and calculated the physical distance between the transmitting point and the receiving point, but the results show the discontinuity of the error. This is mainly due to the fact that the influence of the multipath effect on the ranging error is more significant than that of signal attenuation in the range of 20 m. Considering the relationship between ranging error and distance, we decided to only add edges within 15 m to the graph.

Influence of Threshold δ in Huber Kernel Function
The Huber kernel function can reduce the effect of large error edges in the graph. If the value of δ is set too small, the weight of most edges will be reduced, which will weaken the role of precise edges. If it is set too large, the weight of large error edges will not be reduced, which will cause the Huber kernel function to lose its effect. So, the δ value should be chosen appropriately.
The experimental design of this part was as follows. We obtained the initial coordinates of 84 TPs by KNN matching algorithm and set the Huber kernel function threshold δ from 0.5 to 3.5. Seven control groups were set up with 7,9,11,13,15,17, and 19 nodes per group. Each group randomly selected a corresponding number of points from 84 TPs for optimization and repeated the experiment 1000 times to obtain sufficient data. We also added a comparison experiment without the Huber kernel function in each group. The results in Figure 8 show that the mean positioning error and 75th percentile error were reduced by adjusting the value of threshold δ from 0.5 to 2. Adjusting the value of threshold δ from 2 to 3.5, the mean positioning error and 75th percentile error increased. Obviously, setting δ to 2 is reasonable.

Influence of Threshold in Huber Kernel Function
The Huber kernel function can reduce the effect of large error edges in the graph. If the value of is set too small, the weight of most edges will be reduced, which will weaken the role of precise edges. If it is set too large, the weight of large error edges will not be reduced, which will cause the Huber kernel function to lose its effect. So, the value should be chosen appropriately.
The experimental design of this part was as follows. We obtained the initial coordinates of 84 TPs by KNN matching algorithm and set the Huber kernel function threshold from 0.5 to 3.5. Seven control groups were set up with 7,9,11,13,15,17, and 19 nodes per group. Each group randomly selected a corresponding number of points from 84 TPs for optimization and repeated the experiment 1000 times to obtain sufficient data. We also added a comparison experiment without the Huber kernel function in each group. The results in Figure 8 show that the mean positioning error and 75th percentile error were reduced by adjusting the value of threshold from 0.5 to 2. Adjusting the value of threshold from 2 to 3.5, the mean positioning error and 75th percentile error increased. Obviously, setting to 2 is reasonable.

Impact of Number of Nodes in Graph
In the experiment, we obtained the initial coordinates of 84 TPs based on the KNN algorithm and set the Huber kernel function threshold to 2. Nine control groups were set up with 7,9,11,13,15,17,19,21, and 23 nodes per group. Each group randomly selected a corresponding number of points from 84 TPs for optimization and repeated the experiment 1000 times to obtain sufficient positioning coordinates. Figure 9 shows that the optimization effect improves as the number of nodes added to the graph for optimization increases from 7 to 23, and the magnitude of the improvement gradually decreases and the error curve tends to be steady.
When the number of nodes reaches 19, the reduction of the mean error and 75th percentile error is less than 2 cm for every two additional nodes. When the number of nodes is 7, the average distance between nodes is often greater than the threshold of 15 m, and too few edges can be added to the graph. Therefore, insufficient constraints and large errors will cause the failure of the total error function to converge and an improvement of positioning accuracy that is not obvious. When the number of nodes reaches 9, the total error function can converge in most cases, further improving the positioning accuracy. The experimental results show that under the KNN algorithm, we can obtain relatively stable positioning results when the number of nodes reaches 19 and is set to 2. Figure  10 shows the error cumulative distribution function (CDF) curves of different numbers of nodes.

Impact of Number of Nodes in Graph
In the experiment, we obtained the initial coordinates of 84 TPs based on the KNN algorithm and set the Huber kernel function threshold δ to 2. Nine control groups were set up with 7,9,11,13,15,17,19,21, and 23 nodes per group. Each group randomly selected a corresponding number of points from 84 TPs for optimization and repeated the experiment 1000 times to obtain sufficient positioning coordinates. Figure 9 shows that the optimization effect improves as the number of nodes added to the graph for optimization increases from 7 to 23, and the magnitude of the improvement gradually decreases and the error curve tends to be steady.
When the number of nodes reaches 19, the reduction of the mean error and 75th percentile error is less than 2 cm for every two additional nodes. When the number of nodes is 7, the average distance between nodes is often greater than the threshold of 15 m, and too few edges can be added to the graph. Therefore, insufficient constraints and large errors will cause the failure of the total error function to converge and an improvement of positioning accuracy that is not obvious. When the number of nodes reaches 9, the total error function can converge in most cases, further improving the positioning accuracy. The experimental results show that under the KNN algorithm, we can obtain relatively stable positioning results when the number of nodes reaches 19 and δ is set to 2. Figure 10 shows the error cumulative distribution function (CDF) curves of different numbers of nodes.  It can be observed that the spacing between the curves decrease nodes increases, and the curves with 19, 21, and 23 nodes are tightly pa gap between the curves is especially obvious at 2-5 m, indicating that zation model can significantly improve the positioning accuracy at 2of nodes increases.  It can be observed that the spacing between the curves decrease nodes increases, and the curves with 19, 21, and 23 nodes are tightly pa gap between the curves is especially obvious at 2-5 m, indicating that zation model can significantly improve the positioning accuracy at 2of nodes increases.

Impact of Distribution of Nodes in Graph
The 84 TPs were labeled in U-shaped order with markers ranging It can be observed that the spacing between the curves decreases as the number of nodes increases, and the curves with 19, 21, and 23 nodes are tightly packed together. The gap between the curves is especially obvious at 2-5 m, indicating that our graph optimization model can significantly improve the positioning accuracy at 2-5 m as the number of nodes increases.

Impact of Distribution of Nodes in Graph
The 84 TPs were labeled in U-shaped order with markers ranging from 1 to 84, which were sequentially divided into three groups. Each group contained 28 TPs, representing one of the three passages in the experimental area. We defined a good distribution as a random selection of nodes from only one of the three groups for graph optimization, and a bad distribution as a random and even selection of nodes from each of the three groups for graph optimization. In the experiment, we set the Huber kernel function threshold δ to 2 and separately conducted group experiments for two distributions. The number of nodes in each group was different and the experiment was repeated 1000 times per group. Figure 11a shows the positioning errors under different numbers of nodes in good distribution. It can be found that the mean errors and 75th percentile errors of the good distribution are lower than the random distribution (as shown in Figure 9) under the same number of nodes. The reduction of the mean error is less than 2 cm for every two additional nodes in 15 nodes, and the reduction of 75th percentile error is less than 2 cm for every two additional nodes in 17 nodes. Figure 11b shows the positioning errors under different numbers of nodes in bad distribution. In addition to the fact that the errors are higher than the random distribution, the reduction of the mean error and 75th percentile error is less than 2 cm for every two additional nodes in 23 nodes. In the bad distribution, the positioning accuracy decreases compared with the KNN algorithm when the number of nodes is 7, and the positioning error is only slightly lower than the KNN algorithm when the number of nodes reaches 9, which indicates that the number of nodes in our proposed model needs to be greater than or equal to 9 for better work in our experimental area. additional nodes in 15 nodes, and the reduction of 75th percentile error is less than 2 cm for every two additional nodes in 17 nodes. Figure 11b shows the positioning errors under different numbers of nodes in bad distribution. In addition to the fact that the errors are higher than the random distribution, the reduction of the mean error and 75th percentile error is less than 2 cm for every two additional nodes in 23 nodes. In the bad distribution, the positioning accuracy decreases compared with the KNN algorithm when the number of nodes is 7, and the positioning error is only slightly lower than the KNN algorithm when the number of nodes reaches 9, which indicates that the number of nodes in our proposed model needs to be greater than or equal to 9 for better work in our experimental area.
(a) Error in good distribution (b) Error in bad distribution The experimental results show that the error tends to level off with the node increase under both good and bad distributions, and the errors of a good distribution tend to level off faster than those of a bad distribution. This shows that our graph optimization model is influenced by the node distribution. We find that the better the node distribution, the higher the positioning accuracy for the same number of nodes. The positioning accuracy of a good distribution is more likely to achieve relatively stable optimization results.

Comparison of Wi-Fi Fingerprint-Matching Algorithms after Optimization
The proposed graph optimization model can be used to optimize Wi-Fi fingerprintmatching algorithms other than the KNN algorithm. In the experiments, we constructed the user graph with 19 random nodes and set the Huber kernel function threshold to 2 to optimize four Wi-Fi fingerprint-matching algorithms. The error CDF curves of the KNN, W-KNN, GK, and Stg algorithms before and after optimization are shown in Figure  12. Obviously, the proposed graph optimization model can improve the positioning accuracy. We find that the GK algorithm has the highest positioning accuracy, and the ranking of positioning effects before and after optimization has not changed.
The cumulative error probability of different algorithms under a fixed accuracy limit is shown in Table 3. It can be seen that the optimized algorithm (KNN, WKNN, and GK) does not dominate within 1 m, but the positioning effect is gradually improved from 1 m to 3 m, being most obvious at 3 m. This shows that our proposed algorithm can effectively correct node coordinates with large errors. The optimized Stg algorithm is slightly differ- The experimental results show that the error tends to level off with the node increase under both good and bad distributions, and the errors of a good distribution tend to level off faster than those of a bad distribution. This shows that our graph optimization model is influenced by the node distribution. We find that the better the node distribution, the higher the positioning accuracy for the same number of nodes. The positioning accuracy of a good distribution is more likely to achieve relatively stable optimization results.

Comparison of Wi-Fi Fingerprint-Matching Algorithms after Optimization
The proposed graph optimization model can be used to optimize Wi-Fi fingerprintmatching algorithms other than the KNN algorithm. In the experiments, we constructed the user graph with 19 random nodes and set the Huber kernel function threshold δ to 2 to optimize four Wi-Fi fingerprint-matching algorithms. The error CDF curves of the KNN, W-KNN, GK, and Stg algorithms before and after optimization are shown in Figure 12. Obviously, the proposed graph optimization model can improve the positioning accuracy. We find that the GK algorithm has the highest positioning accuracy, and the ranking of positioning effects before and after optimization has not changed.   Table 4 shows the optimization results of four Wi-Fi fingerprint-m rithms. In terms of the mean positioning error, the KNN algorithm improve ing accuracy by 23.1%, WKNN by 21.5%, GK by 19.4%, and Stg by 27.2%. Th that the 75th percentile error and SD are also reduced after optimization, wh effectiveness of mixed BLE ranging and Wi-Fi fingerprint positioning based optimization model at improving the positioning accuracy. We find that initial localization accuracy, the more obvious the improvement of localiza which indicates that the proposed model works better for Wi-Fi matching al larger errors.  The cumulative error probability of different algorithms under a fixed accuracy limit is shown in Table 3. It can be seen that the optimized algorithm (KNN, WKNN, and GK) does not dominate within 1 m, but the positioning effect is gradually improved from 1 m to 3 m, being most obvious at 3 m. This shows that our proposed algorithm can effectively correct node coordinates with large errors. The optimized Stg algorithm is slightly different at 1 m and 1.5 m, which is caused by the uneven error distribution of the Stg algorithm in the range of 0-2 m, as shown in Figure 12. The obvious optimized results at 2 m, 2.5 m and 3 m of the Stg algorithm again validated our conclusion.  Table 4 shows the optimization results of four Wi-Fi fingerprint-matching algorithms. In terms of the mean positioning error, the KNN algorithm improves the positioning accuracy by 23.1%, WKNN by 21.5%, GK by 19.4%, and Stg by 27.2%. The results show that the 75th percentile error and SD are also reduced after optimization, which shows the effectiveness of mixed BLE ranging and Wi-Fi fingerprint positioning based on the graph optimization model at improving the positioning accuracy. We find that the lower the initial localization accuracy, the more obvious the improvement of localization accuracy, which indicates that the proposed model works better for Wi-Fi matching algorithms with larger errors. From the experimental results, it can be seen that the more accurate the initial positioning results output by Wi-Fi positioning, the more accurate the optimized positioning results, which indicates that the positioning results of our graph optimization model are sensitive to the initial positioning.

Conclusions and Further Work
We integrated Wi-Fi fingerprint positioning and BLE ranging technology to locate users in indoor environments based on a graph optimization model. After deriving the relationship between the physical distance and the RSSI of BLE, a BLE ranging simulation was carried out based on the actual RSSI standard deviation. We obtained initial positioning results through Wi-Fi fingerprint positioning, built a graph of users, added edges between user nodes based on BLE ranging data, and summarized them in a large-scale total error function, which we solved to obtain a better structure of nodes. A drifting solution based on the affine transformation was applied to obtain the final optimized position of the users. The experiments showed that the proposed graph optimization model has a relatively good positioning result when the kernel function threshold δ was set to 2. The reduction of the average positioning error and 75th percentile error was less than 2 cm for every two additional nodes when the number of nodes in the graph reached 19. Our graph optimization model can be influenced by the node distribution. The positioning accuracy of a good distribution is more likely to achieve relatively stable and accurate optimization results. The application of the proposed model to KNN, WKNN, GK, and Stg resulted in improvements, which shows that the proposed model is versatile. The positioning results of the four optimized Wi-Fi fingerprint-matching algorithms showed that the initial positioning results can affect the optimized positioning results.
Although the proposed model can improve the positioning accuracy of traditional Wi-Fi fingerprint positioning, there is room for improvement. The experiments showed that the positioning accuracy rankings of other Wi-Fi algorithms before and after optimization are basically the same, and they verified the effect of the initial positioning results on the optimization results. In addition, we hypothesize that the accuracy of ranging between the user nodes in the graph can affect the final positioning results. Therefore, our next step will be to go beyond the use of BLE beacons for ranging. For example, wireless sensors with higher accuracy can be used to construct more reliable edges of the user graph.
Author Contributions: Conceptualization, R.Z. and J.T.; methodology, R.Z.; software, P.C. and F.M.; validation, P.C. and F.M.; formal analysis, J.T.; investigation, R.Z.; resources, R.Z.; data curation, P.C.; writing-original draft preparation, R.Z.; writing-review and editing, J.T. and P.C.; visualization, P.C.; supervision, R.Z.; project administration, R.Z.; funding acquisition, R.Z. All of the authors contributed to the interpretation and discussion of the results, as well as playing a role in drafting the final manuscript. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement:
The experiment uses an internal data set. The data presented in this study are available on request from the corresponding author.