A Novel Parameter Estimation Method Based on LSU-EKF for Polynomial Phase Signal

The parameter estimation problem for polynomial phase signals (PPSs) arises in a number of fields, including radar, sonar, biology, etc. In this paper, a fast algorithm of parameter estimation for monocomponent PPS is considered. We propose the so-called LSU-EKF estimator, which combines the least squares unwrapping (LSU) estimator and the extended Kalman filter (EKF). First, the coarse estimates of the parameters of PPS are obtained by the LSU estimator using a small number of samples. Subsequently, these coarse estimates are used to initial the EKF. Monte-Carlo simulations show that the computation complexity of the LSU-EKF estimator is much less than that of the LSU estimator, with little performance loss. Similar to the LSU estimator, the proposed algorithm is able to work over the entire identifiable region. Moreover, in the EKF stage, the accurate estimated results can be output point-by-point, which is useful in real applications.


Introduction
Parameter estimation of Polynomial phase signals (PPSs) are common in many fields, including radar, sonar, geophysics, radio communication and biology [1][2][3].In radio communication, continuous-phase signals are used to generate signals with attractive spectral efficiency [4,5].The parameter estimation for PPSs is well indicated for the underwater monitoring of vessels and marine fauna [6,7].It is known that frequency modulated (FM) signals are frequently used in radar and sonar [8], where the phase shift caused by the relative target motion can be modeled as a continuous function of time.In radar applications, the performance of target detection and motion compensation can be significantly improved by using the a priori knowledge of motion parameters.Therefore, precise estimations of range, velocity, and acceleration obtained from the echoes of high-speed targets, such as rockets and missiles, are of great importance for radar detection and imaging [9,10].In recent years, parameter estimation of PPS has attracted considerable attention and many methods have been proposed [11][12][13][14][15][16][17][18][19][20][21][22][23][24][25].Among these methods, a typical estimator of the unknown coefficients is the least squares estimator [11,12], which is an effective approach, being both computationally efficient and statistically accurate when the noise is white and Gaussian.However, this method becomes extremely complex for PPS with second and higher orders.
To estimate the parameters of higher-order PPS, researchers have developed more efficient methods based on "multilinear transforms" or "phase unwrapping".Examples of multilinear transform techniques include the cubic phase function (CPF) [13,14], the high order ambiguity function (HAF) [15,16] and the product HAF (PHAF) [17], which have good performance with less computational complexity, and can be adapted to process multicomponent signals.Nevertheless, for the HAF-based approaches, the phase differentiation (PD) is used to decrease the order of PPSs, which increases the number of noise terms.As a result, each PD increases the SNR threshold and the MSE.Another drawback of the HAF-based approach and the CPF estimator, is that the region of coefficients searched by these estimators is much smaller than the size of the identifiable region.One of the alternative methods is phase unwrapping (PU), which utilizes the fact that the phase of a complex sinusoid is a linear function which is wrapped into (0, 2π).A classical estimator based on PU is the least squares unwrapping estimator (LSU) [18,19].The significant advantage of this approach is that it works for PPS contained anywhere inside the identifiable region.The estimator also appears to perform well when the low signal-to-noise ratio (SNR) threshold is reached.The major drawback of the LSU estimator is that, in general, computing a nearest lattice point leads to a large computational cost.Especially, the algorithm becomes prohibitively complex when the number of samples is larger than 1000 [19].
The parameter estimation by using the extended Kalman filter (EKF) has been investigated in the case of PPS affected by Gaussian noise [20].A disadvantage is that the convergence of EKF depends strongly on the initial values.In this paper, a novel method which combines the LSU estimator and the EKF is proposed.A small number of initial samples are used to estimate the parameters of PPS.Then these coarse estimates are used to initiate the EKF.Meanwhile, a theoretical derivation is conducted on the proper number of samples used in the initialization step.Numerical simulations show that the proposed method can work with low SNR threshold and exhibit high performance with small computational complexity.Moreover, this method is accurate over a much wider range of parameters than other state-of-the-art existing estimators.
The rest of this paper is organized as follows: The signal model is presented in Section 2 and the LSU estimator in introduced in Section 3.Then, the state-space representation of PPS, and the EKF based parameter estimation algorithm are introduced in Section 4. In Section 5, the proper number of samples used in the initialization phase is investigated, and the arithmetic complexity of the proposed method is analyzed.To validate the effectiveness and advantages of the proposed algorithm, Monte-Carlo simulations are conducted in Section 6.Finally, conclusions are given in Section 7.

Signal Model
The received PPS can be written as where A is the signal amplitude, and w(n) denotes the complex additive white Gaussian noise (AWGN), with mean zero and variance σ 2 .θ(n) is a polynomial of order m, which is expressed by where u i , i = 0, . . ., m, are real-valued and unknown coefficients, and their identifiable regions are written as [26] −

The Least Squares Unwrapping Estimator
The LSU estimator attempts to estimate the coefficients by performing linear regression on the phase of received signal.The angle of y(n), denoted y(n), is written as where φ n is a random variable representing the phase noise induced by w(n), and it can be calculated as [28] φ n = According to the fundamental of least square method, the parameters u 0 , u 1 , ..., u m can be evaluated by where N represents the number of samples.Due to the phase ambiguity in y(n), we may define where x denotes the nearest integer to x.By considering ϕ n = θ n − ∑ m k=0 u k n k as a nuisance parameters, SS can be rewritten as The LSU estimator is the minimizer of SS(u) over the identifiable region, which can be represented as a nearest lattice point problem [29].When the number of samples is less than 60, the sphere decoder [30] can exactly compute a nearest lattice point.On the other hand, if the number of samples is larger than 60, another approximate nearest point algorithm called the K-best algorithm [31] is applied, which requires O(N 3 log N) operations.As a result, the LSU estimator needs much more operations as the number of samples, i.e., N increases.
According to the asymptotic normality and strong consistency of LSU estimator [18,19], the mean square error (MSE) of this estimator can be obtained by where P signal denotes the power of the polynomial phase signal, f (•) is the probability density function of φ n and A./B means element-by-element division between matrix A and matrix B. The superscript -1 denotes the inverse operator.Matrix D is a m + 1 by m + 1 Hilbert matrix with elements and the matrix G can be defined as According to (9), more samples for the LSU estimator should be chosen to make the EKF work effectively at lower SNR.However, this will increase the computational cost conversely.

EKF based Parameter Estimation for PPSs
The state-space model associated with PPSs is described by two equations: the state transition equation and the observation equation: where X n is the discrete-time state vector and Y n is the observation vector, υ n−1 represents the state noise vector and ω n is the vector of noise in the measured signal.

The Observation Equation
Since the signal ( 1) is complex, it is necessary to express the measured signal as a 2 × 1 vector in terms of its real and imaginary part In this sense, the observation equation takes the form where the 2 × 1 nonlinear observation function h(X n ) is written as The observation noise vector is defined as In this paper, we assume that the process-noise υ n and the observation noise ω n are not correlated, i.e., Consequently, the observation model of PPSs is nonlinear.

The State Transition Equation
The state-model of PPSs can be described by a set of continuous-time differential equations, that is, where the (m + 1) × (m + 1) constant matrix A and the (m + 1) × m constant matrix L, which characterize the behavior of the model, are written as and υ(t) is the process noise with a process noise matrix The discrete-time state vector X n in the case of PPSs is given by the phase and its derivatives where θ (m) (n) denotes the m order derivative of θ(n).In order to use the discrete-time state (21), the model ( 18) must be discretized somehow [27], and can be described in the form where F n is the transition matrix of the dynamic model, which can be represented as where T s is the step-size of the discretization.On the other hand, υ n−1 with variance Q n−1 is the discretized form of υ(t) and Q n−1 satisfies Therefore, the state-model of PPSs is linear.

Parameter Estimation of PPSs using EKF Algorithm
Since the observation model is nonlinear, we use the EKF algorithm to estimate the parameters of PPSs: Predict Equations: Update Equations: The relationship between the parameters u i , i = 0, • • • , m, of PPSs and the initial state vector X 0 is described by where the matrix C is defined as Therefore, the parameters of PPSs is estimated in the form Using a Taylor series based transformation, the EKF extends the scope of Kalman filter to nonlinear optimal filtering problems by forming a Gaussian approximation to the joint distribution of state X and measurements Y.In this sense, the EKF is sensitive to the initial conditions.If the filter is improperly initialized, the filter may fail.In this paper, in order to make the filter work effectively, we utilize the estimated results of the LSU estimator to initialize the EKF.Considering the computational cost of the LSU estimator, only a few samples are used to obtain the initial estimates.

Performance Analysis of The Proposed Method
As the convergence of EKF depends strongly on the initial values of coefficients, in the proposed LSU-EKF estimator, the LSU estimator is used to obtain the initial values for the EKF and make the filter work effectively.As described above, the computational cost of LSU estimator increases significantly as the number of samples increases.In this section, we first discuss about the appropriate number of samples used in the initialization phase.Meanwhile, the computational complexity of the proposed method is analyzed.

The Number of Samples Used in The Initialization Phase
In order to obtain the appropriate number of samples used in the initialization phase, it is necessary to analyze the stability and convergence of the EKF.Due to the fact that the EKF equations are approximate nonlinear, it is difficult to conduct the convergence analysis for these systems.Many research have been conducted to deal with this problem [32][33][34][35].In general, the sufficient conditions to ensure convergence of the EKF in this paper can be described as follows: where and m+1) are unknown time-varying diagonal matrices, i.e., and which satisfy where λ max (∼) represents the maximum eigenvalue, and On the other hand, the state-space modelization of polynomial phase signals is linear, we can obtain Therefore, the appropriate number of samples used in the initialization phase can be obtained when satisfying (35).

Complexity Analysis
The main process of the EKF composes of two stages: prediction and updating.Assume that the discrete-time state vector X n of PPSs is an i × 1 vector and the observation vector Y n is a j × 1 vector.The arithmetic complexity of the proposed method can be derivated as follows.
The number of arithmetic operations of the predict equations ( 26) and the updating equations ( 27) are 2i 3 + 3i 2 + ij and i 3 + j 3 + 4ji 2 + ij 2 + i 2 + j 2 + 3ij + j + i, respectively.Thus, the number of arithmetic operations of the EKF is In the initialization process, the LSU estimator using the K-best algorithm [31] requires O(N 1 3 log N 1 ) (N 1 means the number of samples used for the LSU estimator) arithmetic operations.
Denoting the number of arithmetic operations in the initialization process as T LSU , then the computational complexity of the proposed algorithm is This shows that the computational complexity of the proposed method is O(N).

Simulation
In this section, to find the minimum number of samples used in the initialization step, the proposed method is first evaluated with different number of initial samples.Then we compare the the performance of the proposed method with other existing estimator, in the cases of big and small coefficients, respectively.In each experiment, 2000 times of Monte Carlo simulations have been carried out.

Numerical Analysis on the Number of Samples in the Initialization
Step First, we conduct some simulations at different SNRs.When the SNR is 12 dB, and the number of samples used in the simulations is 499, the sample MSEs of the four coefficients (u 0 , u 1 , u 2 , u 3 ) for PPSs of order 3 are shown in Fig. 1  LSU-EKF estimator performs well when the number of initial samples is over 12. Considering the complexity of LSU estimator, the optimal number of initial samples for the LSU-EKF estimator is 12 when the SNR is 12 dB.In the same way, the optimal number of initial samples for PPSs of order m = 3 and m = 4 as shown in Table 1 and Table 2 at different SNRs, respectively.Then, the performances of the LSU estimator [18], the HAF estimator [16], the PHAF estimator [17], the Zhou and Wang (ZW) estimator [36] are compared with the proposed LSU-EKF estimator under the same experimental parameters.In practice, one would prefer to keep the number of samples used for the LSU estimator constant rather than changing the value with different conditions whether the optimal value is known or not.Here the number of samples used in the simulations is set to 199 and the number of initial samples used in our method is set to 64.The MSEs of the LSU estimator, the HAF estimator, the PHAF estimator, the ZW estimator, and the proposed estimator are shown in Fig. 2 and Fig. 3, where the values of coefficients are u = 1 4 , 1 4 , 1 8 , 1 24 , 1 96 , 1 480 T for the 5-order PPSs and for the 3-order PPSs.The volume of the coefficients suitable for the HAF estimator, the PHAF estimator and the CPF estimator is much smaller than the identifiable region [14,15,17].Hence, in this case, these estimators will fail.Alternatively, the ZW estimator can work on the entire identifiable region by extending the Euclidean algorithm, which comes with the cost of performance, including MSE and the threshold SNR.The proposed LSU-EKF estimator can reach a threshold of 6 dB for the 3-order PPSs and of 17 dB for the 5-order PPSs.The corresponding threshold SNRs for the LSU estimator are 2 dB for the 3-order PPSs and 12 dB for the 5-order PPSs, respectively.The reason is that the number of samples used for the LSU estimator is larger than that for the LSU-EKF estimator in the initialization step.Table 3 shows the time costs of all the estimators.It is seen that the computational complexity of the LSU estimator is approximately 24 times as much as that of the LSU-EKF estimator for the 3-order PPSs and 15 times for the 5-order PPSs.This simulation is performed on a high-performance computer with Intel(R) Core(TM) i5-4590 CPU 3.30 GHz and 8.0 GB of RAM.In fact, the SNR thresholds for the LSU-EKF estimator can be lower by using more samples for initialization, with increased time cost.Therefore, it is a tradeoff between the performance and the computational complexity.the 3-order PPSs.In Fig. 4, it can be seen that the HAF estimator reaches a threshold of 6 dB, which is the same as the proposed method for the 3-order PPSs.For the 5-order PPSs, the threshold of the HAF estimator is about 18 dB, slightly higher than that of the proposed method.Also, the PHAF estimator has a good performance on estimation of the coefficients u 3 for the 3-order PPSs and u 5 for the 5-order PPSs.As expected, the performance of the proposed method remains unchanged for the set of small coefficients.Finally, the proposed method is investigated using real radar data.Fig. 5(a) shows the estimated velocity in comparison with the real velocity of the target, while the errors between the estimated velocity and the real velocity are shown in Fig. 5(b).It is seen that at the beginning of the velocity estimation, the estimated errors are relatively high, as there is an error in the initial values obtained by the LSU estimator.Nevertheless, after the initial transient period (about 100 pulses), an accurate velocity estimation is established by the proposed LSU-EKF estimator.The average error of the estimated velocity by the proposed LSU-EKF estimator is nearly 0 m/s and the RMSE of the estimated velocity is about 0.0076 m/s in the steady state.

Figure 1 .
Figure 1.Estimation for PPSs of order 3 by the LSU-EKF estimator with different initial samples.(a)MSE for the coefficient u 0 .(b)MSE for the coefficient u 1 .(c)MSE for the coefficient u 2 .(d)MSE for the coefficient u 3 .

Table 1 .
The optimal number of samples used in the initialization phase at different SNRs for PPSs of order 3.

Table 2 .
The optimal number of samples used in the initialization phase at different SNRs for PPSs of order 4.

Table 3 .
Time cost of different methods for PPSs of order 3 and PPSs of order 5. (in seconds)In order to make the values of the coefficients suitable for the HAF and PHAF estimator, we set Tfor Preprints (www.