Bias Reduction for Direct Frequency Estimator of Complex Sinusoid

Frequency estimation of single-tone complex sinusoid is a fundamental problem in signal processing. An efficient way to estimate the frequency is to directly calculate the frequency from three discrete Fourier transform coefficients in the neighborhood of peak magnitude spectrum. Because of the approximation used in the calculation, these direct estimators exhibit non-uniform biases. In this paper, we first derive a general analytical expression for the direct estimators. Applying Taylor series to the general expression, we find two new direct estimators that can eliminate the bias up to the fourth order approximation. Experiments show that the new estimators can reduce the bias by at least five orders at large frequency offset as compared to the current best direct estimator.


Introduction
Frequency estimation of a single-tone complex sinusoidal signal is a fundamental problem in digital signal processing and has applications in areas such as radar, sonar, measurement, wireless communication, and speech analysis.The signal to be estimated can be described in the following equation: where 0 A , 0 f , and 0  are the amplitude, frequency, and phase of the signal, respectively, ] [n w is an additive noise, S f is the sampling frequency, and N is the number of samples.
For a zero-mean white Gaussian noise, it has been proved that the maximum likelihood frequency estimation is the frequency that maximizes the magnitude of the periodogram of ] [n s [1]: The Cramér Rao lower bound (CRLB) of the mean squared error (MSE) is: Where  is the SNR of ] [n s .Since Eq. ( 2) is a search in a continuous domain f , it is difficult to implement as an efficient computational algorithm.Therefore, the task to estimate the frequency is usually divided into two steps: a coarse search and a fine search.The coarse search finds the peak magnitude of the discrete Fourier transform (DFT) of The above equation gives the relationship between the estimated frequency and the sampling frequency S f .If only DFT is used for the frequency estimation,  ˆ is zero and the frequency resolution is N f S / .The estimation of the frequency offset  ˆ makes the estimated frequency much more accurate.
The methods for fine search can be divided into two categories: iterative approaches and direct approaches.For the iterative approaches,  ˆ is repeatedly refined.They can be further divided into two categories.The first is to repeatedly narrow the search range using dichotomous search [2][3][4][5].The second is to iteratively refine a single solution toward the final solution [1,[6][7][8][9][10][11][12].
For the direct approaches,  ˆ is directly calculated from the three spectrum lines ( ) [13][14][15][16][17][18] or more [15,19].Most direct estimators using three spectrum lines [13,14,16,17] can be expressed as the following general equation: , and C are all constant real numbers.The values of i a , j b , and C for some popular estimators are shown in Table 1.In the table, the Quinn estimator is modified from its original form to fit the general equation but the result is equivalent.Among the direct estimators, the performance of the estimator proposed by Candan [17] has been proved to be the best experimentally.
Due to the inevitable approximation in calculating the estimated frequency, both direct and iterative approaches exhibit estimation bias between  ˆ and the true offset  .The bias is usually non-uniform across the range of  .
While CRLB is the limiting factor of the estimator performance when SNR is low, the bias of the estimator is the limiting factor when SNR is high.Since the bias diminishes as the number of iteration increases for iterative estimators, it is critical to reduce bias for direct estimator to improve its performance.In this regard, Candan used a constant scaling factor to correct the result of the Jacobsen estimator [16] and achieved significant reduction on the bias.In this paper, we first derive a general analytical expression for the direct estimators described by Eq. ( 6).After applying Taylor series into the analytical expression, we find two new estimators that eliminate the higher-order terms up to the fourth order approximation and, thereby, achieve even better bias reduction than Candan's approach.We also show that the Candan estimator removes some of the higher-order terms but not all and that is the reason the new estimators are better than the Candan estimator.In the experiments, we show that the new estimators can reduce the bias by at least three orders at large  as compared to the Candan estimator.
This paper is organized as follows: In the second section, we present the theoretical analysis of the direct estimators.Through Taylor series, we derive the expressions for the new estimators and analyze the Candan estimator.In the third section, we show the experimental results.In the fourth section, we conclude.

Theoretical Analysis
The first step to derive an analytical expression for Eq. ( 6) is to derive the analytical expressions for the three   1).Plug Eq. ( 1) and Eq. ( 5) into Eq.( 4), we get: The equation is the sum of a geometric series and the three spectrum lines can be expressed as: Plugging the above three equations into Eq.( 6) and after simplification, we get: where the terms 1 T , 2 T , and 3 T are: Then, Eq. ( 11) can be written as: For Eq. ( 15) to converge, 2 r must be greater than 2 q .
However, if , only the first term in Eq. ( 15) remains and it is: As can be seen from the above analysis, the complex exponential term ) / exp( N j and its conjugate in Eq. ( 11) require an additional layer of approximation in Therefore, it will be beneficial to remove them from Eq. (11).To achieve this, we can multiply a constant exponential term After this simple correction, the new general equation becomes: We will call this estimator as "phase corrected" estimator and the estimator in Eq. ( 6) as "uncorrected" estimator.After the correction, the phases of the three spectrum lines are synchronized.Thus, the corresponding analytical expression for the phase corrected estimator is: It is worth noting that the above equation is the exact analytical expression without any approximation.
For the following analysis, we will apply Taylor series of sinusoidal function and terms higher than the fourth order in the series are neglected.Plugging Taylor series into Eq.( 12)-( 14), we get: For the uncorrected general equation (11), we also For both the estimators in Eq. ( 11) and Eq. ( 24) to be a good estimator, the numerator should approximate  and the denominator should approximate a constant.It is especially important that the denominator contains a constant term because the denominator will be zero at 0   without the constant term and the estimator becomes unstable.An additional benefit to maintain a constant term in the denominator is that it ensures the convergence constraint:  24) because they remove both the second and fourth order terms and maintain the first order  term.The second is that 2 T must be present in the denominator to maintain a constant term, i.e., 2 b must not be zero.
Therefore, Eq. ( 16) for the uncorrected case can be defined as: We then define Eq. ( 17) as: The above equation is directly proportional to Eq. (30).Therefore,   2 1 / r r up to the fourth order approximation.In addition, the fourth order approximations of 1 q and 2 q up to the fourth order approximation.In other words, the remaining terms after 2 1 / r r in Eq. ( 15) are all zeroes up to the fourth order approximation.Therefore, the following estimator: ) is an estimator without distortion up to the fourth order approximation for the uncorrected case.
Similarly, for the phase corrected case, the numerator in Eq. ( 24) can be defined as ] [ 3 1

T T 
and it is: The denominator in Eq. ( 24) can be defined as: ) Combining Eq. ( 33) and Eq. ( 34) yields an unbiased estimator up to the fourth order approximation.
Therefore, the new estimator for the phase corrected case is:

Experimental Results
The first set of the experiments verifies the performance on estimation bias and is conducted with offset and  ˆ is the estimated frequency offset.The absolute value of biases with respect to  for N=128 and N=1024 are shown in Fig. 1.For these experiments,  varies from -0.49 to 0.49 with step size 0.01.Since the bias is symmetric, the figures only show the part for positive  .Fig. 2 shows the absolute value of biases with respect to the number of samples N at two different  .Because the bias for the Quinn estimator [13] is the same as the bias of the Jacobsen estimator, we only show the results of the Jacobsen estimator to represent the results for both Quinn and Jacobsen estimators.In the figures, we can clearly see that both of the proposed estimators provide significant reduction on bias.As N gets larger, the difference between the proposed estimator and the Candan estimator gets larger.In Fig. 1, the minimum gaps between the proposed estimator and the Candan estimator are 10 2.1979 and 10 2.2649 and the maximum gaps are 10 5.1408 and 10 5.4383 for uncorrected and phase corrected cases respectively.In Fig. 2(a), the minimum gaps are 10 2.2550 and 10 2.3990 and the maximum gaps are 10 5.2683 and 10 5.4063 for uncorrected and phase corrected cases respectively.In Fig. 2(b), the minimum gaps are 10 2.6325 and 10 2.8936 and the maximum gaps are 10 5.6433 and 10 5.9044 for uncorrected and phase corrected cases respectively.This means that the proposed estimators can achieve an improvement of at least two orders at small  and at least five orders at large  .f is calculated as in Eq. ( 5) and 0 f is the true frequency.Each experiment is repeated 10,000 times and MSE of the estimated frequency is calculated.Fig. 3 shows MSE with respect to SNR.From the figure, we can see that the new estimators follow CRLB closely in the entire range of SNR while the Candan estimator is limited by its estimation bias and becomes flat at a threshold which depends on  and N.Only at N=128 and  =0.4,we can see that the new estimators become flat at very high SNR and the phase corrected version is better than uncorrected version because of the lower bias.

Conclusions
In this paper, we use Taylor series to derive two new frequency estimators which eliminate the estimation bias up to the fourth order approximation.The biases of the new estimators are significantly better than that of the currently best Candan estimator because the Candan estimator eliminates some of the higher-order terms in the bias but not all of them.Because of the reduced bias, we show in the experiments that the new estimators follow the CRLB closely up to very high SNR while the Candan estimator becomes flat at a much lower SNR threshold.Overall, we can conclude that the new estimators are much better than the Candan estimator at high SNR and are the same as the Candan estimator at moderate and low SNR.Therefore, they are superior to the Candan estimator.

k.
be the index of the peak magnitude in Eq. (4).The fine search finds an offset  ˆ in the vicinity of p k After the two searches, the estimated frequency is calculated as:

2 2 q
r  .From these observations, we can deduce two results.The first is that the best

Fig. 4 shows
Fig. 4 shows MSE with respect to  at two different SNR for N=128.At SNR = 100, we can see that the two new estimators have exactly the same MSE and are much closer to CRLB than the Candan estimator.At 49 .0    , the MSE of the Candan estimator is 835 times of the MSE of the new estimators.At SNR = 50, the two new estimators are the same as the Candan estimator.

Table 1
Coefficients Used in the General Expression for Some Popular Direct Frequency Estimators