GPS Sparse Multipath Signal Estimation Based on Compressive Sensing

A GPS sparse multipath signal estimation method based on compressive sensing is proposed. A new 0 norm approximation function is designed, and the parameter of the approximate function is gradually reduced to realize the approximation of 0 norm. The sparse signal is reconstructed by a modified Newton method. The reconstruction performance of the proposed algorithm is better than several commonly reconstruction algorithms at different sparse numbers and noise intensities. The GPS sparse multipath signal model is established, and the sparse multipath signal is estimated by the proposed reconstruction algorithm in this paper. Compared with several commonly used estimation methods, the estimation error of the proposed method is lower.


Introduction
In wireless communication system, due to the influence of scattering or reflection, radio wave may reach the receiver through multiple different paths, so that the received signal at the receiver is the superposition of multiple signals. The multipath interference is inevitable, which is one of the main reasons that affect the performance of GNSS receiver. It will lead to the distortion of the autocorrelation waveform of the pseudocode and affect the detection of the coherent peak, thus reducing the positioning accuracy of the receiver. It is necessary to detect the multipath signals [1,2].
Candès et al. [3] and Donoho [4] proposed the theory of compressive sensing (CS) in 2006, which can sample the sparse signal at low sampling rate. Currently, the research and application fields of CS include radar signal arrival angle estimation [5,6], satellite navigation signal processing [7][8][9], and sparse radar design [10][11][12][13]. The CS signal reconstruction needs to recover the original high-dimensional signal from a small amount of measurement data. There are several kinds of reconstruction algorithms. The relaxation algorithm transforms 0 norm into 1 norm. It converts the problem into a convex optimization problem which is easy to be solved [14,15]. The greedy algorithm has the advantages of simple calculation and easy implementation. However, the reconstruction accuracy is not high and the sparse number is needed [16,17]. In recent years, some 0 norm approximate solution algorithms have been developed, which transform the 0 norm into mathematical analyzable approximation functions. The sparse signal is reconstructed by a gradient method, Newton method, or other optimization methods [18,19]. This method can reduce the difficulty of signal reconstruction and obtain better results.
A new 0 norm approximation function is proposed. The sparse reconstruction problem of CS is translated into a Lagrange multiplier problem. It gradually reduces the parameters and reconstructs the original signal by the modified Newton method. A GPS sparse multipath signal model is established, and the estimation of GPS sparse multipath signal is realized by using the approximate function which is proposed in this paper. The rest of this paper is organized as follows. In Section 2, the GPS sparse multipath model is introduced. Section 3 analyzes the proposed 0 norm approximate function and reconstruction method. Sections 4 and 5 give the reconstruction and estimation results of the new proposed method. Finally, the paper is concluded in Section 6.

GPS Sparse Multipath Model
Assuming that the GPS receiver has accurately tracked the carrier frequency of the receiver signal and stripped the carrier, the received multipath signal in one navigation message can be expressed as In the formula, A k , τ k , and ϕ k are the amplitude, delay time, and phase of the multipath signal, respectively. sðtÞ is the spread spectrum code, and nðtÞ is the noise [20]. The corresponding signals of all sampling points can be obtained by cyclic shift, so the upper formula can be represented as a matrix form: where S = ½s τ 1 , s τ 2 , ⋯, s τ N and a k = A k e jϕ k . When there are multipath signals, the vector values equal to the amplitude value of the multipath signal, and the other values are 0. The vector a can be represented as The length of a GPS spread spectrum chip is 1/1023 ms. Suppose the sampling frequency of the GPS signal is 51.15 MHz. So there are 50 sampling points in a chip. Assuming that the number of multipath in a chip is 3; that is to say, three elements are not 0 and the others are 0. r is a sparse vector. The GPS multipath signal represented by formula (3) is a sparse model.
There are many GPS multipath estimation methods. In reference [21], a TK operator method is easy to implement but the estimation effect is not ideal in strong noise environment. The method based on maximum likelihood (ML) estimation is proposed in reference [22]. The ML method is improved, and the Generalized Likelihood Ratio (GLRT) is used to improve the performance of the estimation [23]. The MEDLL is used to analyze the signal, and the effect of MEDLL is good, but the method is complex [24]. There are other wavelet-based analysis methods and so on [25,26]. These methods estimate GPS multipath signals without considering the sparse characteristics of GPS multipath. If sparse constraints are added to the estimation, the performance of multipath estimation can be improved theoretically.  [27]. Instead of 0 norm, 1 norm is used to reconstruct the MRI images. The 1 norm is changed to the following formula:

New 0 Norm Approximation Function
In the formula, p is a small positive constant. The approximate formula is used to transform the nondifferential 1 norm into a differentiable function. Then, the MRI image is reconstructed by using the optimization algorithm. A CS reconstruction algorithm based on 0 norm approximation function is proposed [18]. The approximation function is a Gaussian function: The function realizes the approximate expression of 0 norm by decreasing the parameter σ gradually. Then, the sparse signal is reconstructed by the gradient method. Both methods convert a nonconvex 0 norm problem to a convex function which has some deviation in principle. Furthermore, the parameter p of formula (4) remains unchanged, and the reconstruction may have some deviation.
Based on the approximation of 0 norm given by Donoho and Gaussian function, a new approximate function is proposed. The coefficient 1/2 is changed to p/2 by modifying the approximation of 1 norm given by Donoho. An approximate 0 norm can be obtained when the coefficient p decreases gradually to 0. The designed 0 norm approximation function is The model has only one parameter, and the approximate function is a smooth differentiable function. So the gradient and Hessian matrix of the function can be obtained, which facilitates the solution of the model.

Sparse Signal Reconstruction Based on Approximate 0
Norm. Based on the proposed approximation function, the CS reconstruction problem can be described as the following mathematical constraint expression: Φ ∈ R m×n is a measurement matrix. In reference [18], the gradient method is used to reconstruct the sparse signal. The descending direction is "zigzag," which affects the efficiency of reconstruction. In this paper, the Newton method is used to reconstruct the sparse signal for improving the efficiency [28,29]. Newton's method requires positive definite Hessian 2 Wireless Communications and Mobile Computing matrix, so the Hessian matrix needs to be modified in iteration. The CS reconstruction problem is transformed into a Lagrange solution model: y is Lagrange multiplier. The gradient of formula (8) is : The Hessian matrix of formula (8) is obtained by derivation of expression (9): : ð10Þ When p = 1, u i is a positive value, so the Hessian matrix is a positive definite matrix. That is to say, the function is a convex function. The value p = 1 can be used as the starting point of iteration. According to Newton's method, the iterative formula can be expressed as In the formula ( (11) and (12)), α k is the scale factor and d k is the direction of gradient descent. g k−1 is the gradient of the k − 1 iteration which can be obtained by formula (9). In the iteration, u i may be a negative value in formula (10). Formula (12) requires a positive definite Hessian matrix. It needs to modify the value ofu i to guarantee a positive definite matrix in the iteration. The modification formula is δ is a very small positive value which ensures the Hessian matrix is a positive definite matrix. The empirical value is 10 -5 . The initial value of the parameter p in the approximate function is 1. The value is decreased in iteration according to the following formula: It can obtain better reconstruction effect when c belong to ð0:5, 1Þ.
The general idea of sparse reconstruction algorithm in this paper is given as below. The reconstruction has two search processes. For each given value p k in the external loop, an optimal value is obtained by the modified Newton method. An upgrading value is gotten by decreasing p gradually. Then, p k in the outer loop is reduced to the value p k+1 , and the new optimal value is searched again according to the previous optimal solution. The sparse solution is found iteratively. The implementation steps of this reconstruction Algorithm 1 are as follows: 3.3. Convergence Analysis of Approximate Function. In this section, the convergence of the algorithm is analyzed. In order to facilitate the analysis, a lemma is first introduced.

Lemma 1.
Suppose that A ∈ R m×n is a Unique Representation Property (URP) matrix. All m × m submatrixes have inverse matrix [30]. A vector s = ½s 1 , s 2 , ⋯, s n T belongs to zero space of A. If more than n − m elements converge to 0 in the vector, then the vector converges to 0.
Proof. Suppose that A has been normalized by a column vector. For any γ > 0, I γ is a subscript set that all the corresponding s i are larger than γ.Â is a submatrix which is composed of the columns corresponding to the subscript set.
Because A is a URP matrix, each column satisfies the linear independence. The matrixÂ has a left inverse matrix. M = max fkA∧ −1 kg; jI γ j is defined as the number of Input: random measurement matrix Φ, measurement value y; Output: reconstruction signal xR; (1) Initialization: Set the initial parameter: p 0 , λ, c, α k , Total iteration times J, inner iteration L k , reconstruction error threshold ε, the signal initial value x 0 = 0, external iteration count k, internal iteration count t; (2) The external iteration: if k > J, turn to step 5. Otherwise, x k is calculated by Newton method; (3) The internal iteration: t = t + 1, if t > L k , the internal iteration is over and turn to the next step. Otherwise, calculate x t k by x k and p t k which decrease gradually. If E k t = norm ðH k t Þ < ε, turn to the next step; (4) p k+1 = p k , k = k + 1, and turn to step 2, searching the optimal value based on x t k ; (5) All the iteration is over and get the sparse solution xR = x t k ; Output2: There is no satellite in the line of sight.
Algorithm 1: GPS sparse multipath signal estimation based on compressive sensing. 3 Wireless Communications and Mobile Computing elements in the subscript set. From the known conditions, There are 2 norms for the upper formulas: rewrite the upper formula to a matrix form: And Proof of lemma 1 is over.

Corollary 2.
Suppose that A ∈ R m×n is a URP matrix and the equation y = As has a sparse solution. If there is another solutions = ðs 1 ,s 2 , ⋯,s n Þ T and ksk 0 = k < m/2, the vector converges to 0. The proof of corollary can be drawn from the lemma. Because 0 ≤ ks − s 0 k 0 ≤ 2k < m, more than n − m elements are 0, ands − s 0 belongs to the zero space of A. The matrix is a URP matrix, so the difference vectors − s 0 converges to 0. That is to say, the sparse solution is unique.
For the approximate function proposed in this paper, suppose that ζ is a very small positive number. If γ = ±ð1 − ζÞ 1/p , for a given jγj ≤ js i j ≤ 1, the following equations hold: Because the range of p is ð0, 1, so j1 − ζj 1/p < j1 − ζj < 1. For a very small coefficient s i in the vector, the value of the proposed function is approximately 1.
For any given γ = ±ð1 + ζÞ 1/p , the following equations hold: For all coefficients which are greater than 1, the function values are approximately 1. Therefore, the proposed approximate function in this paper realizes the function of 0 norm.
The above analysis theoretically ensures that the approximate function can find the sparse solution of the signal. The constraint condition k < m/2 is a sufficient rather than a necessary condition, and the constraint condition is strong. If this condition is not satisfied in practical application, the algorithm can reconstruct the sparse signal also, but the probability of reconstruction is very low.

Simulation Experiment of Reconstruction Algorithm
Three experiments are designed to verify the feasibility and effectiveness of the algorithm. In experiment 1, the measurement matrix is fixed. The reconstruction performance is analyzed at different sparse numbers. Antinoise performance is analyzed in experiment 2. The running time of the algorithm is also an important standard to test the algorithm. Experiment 3 analyzes the running time of algorithm.

Different Sparse Numbers.
According to the proposed approximate function, a CS reconstruction model is defined as In order to compare the reconstruction effect, the proposed algorithm is compared with OMP, Donoho method (1 norm) [27], and the reference [18] algorithm (SL0). The noise type is zero mean Gaussian white noise, and the noise variance is 0.01 in experiment 1. The measurement matrix is Φ ∈ R 200×500 . The length of measurement value and sparse signal is fixed at 200 and 500, respectively. The simulation is run 1000 times at each sparse number, and the total reconstruction times are counted if the signal-to-noise ratio (SNR) of the reconstructed signal is greater than 25 dB. The experimental results are shown in Table 1.
As can be seen from Table 1, each algorithm can reconstruct the signal completely when the sparse number is less than 50. When the sparse number is greater than 130, the reconstruction SNR of all algorithms are less than 25 dB. When the sparse number increases, the successful reconstruction times of the proposed algorithm in this paper are greater than other algorithms. For example, when the number of nonzero data reaches 120, there are still 7 times reconstructed signals with SNR above 25 dB, but SL0 algorithm only has 3 times and 1 norm and the OMP algorithm can no longer reconstruct the signal above 25 dB. Comparing the data in Table 1, the reconstruction effect algorithm proposed in this paper is better.  Table 2.

Robust
From Table 2, it can be seen that the reconstruction SNR is above 25 dB at low sparse number when the noise variance is 0.01 and 0.02. The proposed algorithm is effective. The reconstruction SNR can reach 25.272 dB when the sparse number is 100 and noise variance is 0.01. The proposed algorithm shows strong robustness. However, when the noise variance is 0.04, the reconstruction effect of this algorithm is not ideal.

Average CPU Time of Algorithm.
The reconstruction time is an important index to estimate the algorithm. The experimental parameter is designed as follows: the measurement matrix Φ ∈ R m×n , m/n = 0:4, and sparse number m = 4k. According to the previous experiments, this numerical setting can effectively reconstruct sparse signals in a noiseless environment. The length of the signal is constantly adjusted, and the CPU times for each of the algorithms averaged over 100 trials are tabulated in Table 3.
From Table 3, it can be found that the running time of several algorithms is close when the total length of the signal is less than 2000 sampling points. However, when the signal length reaches 4000, the reconstruction CPU time of OMP is almost twice that of the proposed algorithm. For highdimensional matrix, the pseudoinverse calculation of the matrix in the OMP algorithm takes too much time. SL0 is solved by the optimal gradient method which appears "zigzag" when finding the optimal value. This reduces the efficiency of finding the optimal solution. In this paper, the modified Newton method is used to improve the reconstruction speed. According to Table 3, the algorithm is faster than the other algorithms, so the algorithm is more suitable for high-dimensional signal processing.

GPS Sparse Multipath Signal Estimation
The GPS sparse multipath model is given in formula (3). The sparse multipath signal is estimated by CS theory in this section. First, the GPS signal is measured by measurement matrix Φ ∈ R m×n : Then, sparse multipath signal is estimated by using the 0 norm approximation function which is proposed in this paper: The multipath model discussed in this paper is based on the following assumptions [31,32]: (1) The transmission delay of the nondirect signal is greater than that of the direct signal, which is in line with the actual model. The nondirect signal is often formed by refraction or reflection which takes more time (2) The amplitude of nondirect signal is lower than that of the direct signal. In the shelter or indoor environment, this assumption is generally not true which needs to estimate multipath with an additional method. But this situation is beyond the scope of this paper which is not analyzed here (3) In this paper, the multipath delay in on chip is only considered. That is to say, the multipath delay is limited to one chip (4) Multiple sparse values may be estimated due to noise interference. Assuming that the number of sparse multipath is known, the maximum sparse value is selected as the multipath signal. The small sparse values is considered as noise and excluded

Comparison of Different Sparse Multipath Estimation
Methods. The TK operator [21] and ML estimate method [22] are selected to estimate the multipath signal. Assuming that the GPS signal has two nondirect signals. The direct signal delay is 0, and the normalized amplitude is 1. The delay time of the nondirect signal changes randomly in a chip. The ranges of absolute amplitude vary randomly between 0.2 and 0.8. The GPS signal length is 1 ms and sampling rate is 51.15 MHz. The SNR of correlation result is 35 dB. The compression rate is m/n = 0:4. The estimation results of different methods are shown in Figure 1. The amplitude of two random nondirect signals is (0.6, 0.4), and the position is (0.2, 0.58) chip. From Figure 1, three methods can correctly estimate the multipath signal delay with some deviation in magnitude. The TK operator and ML estimate method can correctly find multipath signals. Due to noise interference, the estimation amplitude deviation of multipath signals is large, and the  Sparse  number  50  60  70  80 90 100 110 120 130   OMP  1000 998 981 944 729 228 25 0  0  1 norm  1000 1000 995 978 864 403 96 0  0  SL0  1000 1000 1000 996 947 429 185 3  0  This paper  1000 1000 1000 1000 981 506 227 7  0 5 Wireless Communications and Mobile Computing estimated values of nonsparse positions fluctuate greatly which result in more fake multipath signals. The CS estimation method proposed in this paper adds sparse prior information, the estimated multipath signal is consistent with the simulation signal position, and the estimated multipath signal amplitudes are close to the real multipath signal.
ML and TK operators are not effective because of noise interference. After the multipath estimation of the signal, some detection algorithm is used to detect it, and the false multipath signal is excluded in practical application. For example, the GLRT criterion is used to analyze the estimated multipath signal in reference [23]. The proposed CS method adds sparse constraints and selects the largest signal as a sparse multipath signal. It can be seen from Figure 1 that the performance is better than that of ML and TK operators. It is unfair to compare the CS estimation method with the common methods, so it only compares the method of this paper with the method based on the CS theory in the subsequent performance analysis.

Multipath Signal Estimation Based on CS.
It is mainly affected by three factors: noise, sparse number (number of multipath), and measurement matrix dimension m to estimating multipath signal based on CS. This section analyzes the estimated effect from these aspects. First, the sparse estimation results based on CS are analyzed by an experiment. 1 ms GPS signal is selected with sampling rate 40.92 MHz. The sparse number is set to 4, and the delay of nondirect signal varies randomly in a chip. The absolute amplitude varies randomly between 0.2 and 0.8. The delay and amplitude are uniformly distributed. SNR is set to 0 dB. Compression ratio m/n is 0.5. In the simulation, assuming that the number of sparse multipath is known, the estimated results of the proposed method are shown in Figure 2. In a strong noise environment, the amplitude of the multipath signal and the position of the multipath signal may be deviated. Multipath signal estimation error is used to evaluate the estimation results, which is defined as In the formula, N is the total number of sampling point of a chip. K is a sparse number. a i andâ i represent real and estimated multipath signals, respectively. All multipath signal amplitudes are normalized. It can be seen from Figure 2 that if the position and phase of the estimated signal are consistent with the real multipath signal, the estimation error is small. If the position or phase of estimated signal is different from the real signal, the real value and the estimated value are counted into error. Formula (25) can objectively evaluate the estimated performance.
The multipath estimation effect is analyzed at different multipath numbers and noise intensities in experiment 1. The measurement matrix is fixed and m/n = 0:3. The SNR and sparse number vary from -20 dB to 20 dB and 2-8, respectively. The other parameters are set as above. The experiment repeats 100 times at each sparse number. The average estimation error is shown in Figure 3. From the figure, it can be found that the estimation error increases when the sparse number increase. It is consistent with the CS theory that the reconstruction effect decreases with the increase of sparsity at the same measurement matrix. The best estimation effect is obtained when there is only one nondirect   When the multipath number is fixed, the influence of measurement matrix is analyzed in experiment 2. The multipath number is set to 6. The column of the matrix remains unchanged, and m/n range from 0.2 to 0.8. The other related parameters are set in the same experiment 1. The experiment repeats 100 times at each m/n. The average estimation error is shown in Figure 4. The estimation error will be relatively reduced when the number of random measurement is increased, which can retain more information of sparse signal. The estimation effect is not good when the SNR is less than 0 dB.
The reconstruction algorithm has some effect on sparse signal reconstruction in CS, which also has some influence on sparse multipath estimation. The OMP, 1 norm, and SL0 are selected to estimate the multipath signal. The estimation effect of different algorithms is analyzed in experiment 3. The sparse number of signal is set to 5. The m/n of the measurement matrix is 0.5. The GPS signal length is 1 ms with 40.92 MHz sampling rate. The experiment repeats 100 times at each SNR, and the average estimation error is shown in Figure 5. The estimation error of all methods is large at low SNR. When the SNR is higher than 0 dB, the estimation effect is good. As the SNR increases, three multipath estimation errors based on optimization algorithms are slightly lower than OMP methods.
From the three experiments, it can be seen that the sparse multipath signal can be estimated by the CS method effectively. When the SNR is low, the estimation effect is poor. For example, the CS method cannot estimate the multipath  7 Wireless Communications and Mobile Computing signal at -20 dB strong noise interference. As the SNR increases, the estimation error decreases. When the SNR is higher than 5 dB, the estimated multipath signal is basically consistent with the added multipath signal. In this paper, the length of GPS signal is 1 ms. In practical application, the SNR can be improved by increasing the GPS signal length because the time of navigation message is 20 ms, which can obtain lower estimation error than 1 ms.

Conclusions
In this paper, we proposed a new GPS multipath signal estimation method based on CS. A new 0 norm approximate function is proposed. The sparse signal can be reconstructed by using the modified Newton method. A large number of simulation results show that the new method can reconstruct the sparse signal effectively at different sparse numbers and SNR. The reconstruction performance of the proposed algorithm is better than OMP and other optimization algorithms. The GPS signal sparse multipath model is established. The GPS multipath signal can be estimated by the proposed method in this paper. The estimation result is better than the TK and ML method and other CS reconstruction methods. Because there are 20 periods of spread spectrum code in a navigation message, the SNR can be improved by accumulating for a long time, and then, the estimation effect can be improved. Our research team will continue to delve into how to improve the estimation effect of this method.

Data Availability
The simulation data used to support the research of this paper are available from the corresponding author.

Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.