Line spectrum extraction method based on hidden Markov model

To solve the difﬁcult problem of line spectrum detection with a low signal-to-noise ratio, a line spectrum extraction algorithm based on a hidden Markov model (HMM) is proposed. The forward–backward al- gorithm was used to track single and multiple spectral lines in a lofargram, and a new algorithm for the state transition probability matrix was proposed to handle the large amount of HMM calculation and the uncertainty difﬁculty in processing real signals. Finally, a median ﬁtting algorithm was used to correct the few outlier points in the line spectrum estimation process. Simulation and sea trial data showed that the algo- rithm had good line spectrum tracking ability for a low signal-to-noise ratio.

✉ Email: makainue@163.com To solve the difficult problem of line spectrum detection with a low signal-to-noise ratio, a line spectrum extraction algorithm based on a hidden Markov model (HMM) is proposed. The forward-backward algorithm was used to track single and multiple spectral lines in a lofargram, and a new algorithm for the state transition probability matrix was proposed to handle the large amount of HMM calculation and the uncertainty difficulty in processing real signals. Finally, a median fitting algorithm was used to correct the few outlier points in the line spectrum estimation process. Simulation and sea trial data showed that the algorithm had good line spectrum tracking ability for a low signal-to-noise ratio.
Introduction: The radiated noise signals [1] of underwater targets often contain rich components of line spectra. These line spectra are usually caused by mechanical and propeller noise [2], so they are rich in target information.
The commonly used line spectrum extraction methods [3][4][5] often make use of the fact that the spectral level of a line spectrum is higher than that of a continuous spectrum by 10 dB or more [6]. However, for the case of weak targets, this method performs poorly. Because the line spectrum may be submerged in the noise background, its spectrum level may be lower than that of the continuous spectrum. Li et al.'s [7] proposal of a fusion method to extract frequency features proved to be an effective algorithm. Paris et al. [8] proposed a series of methods for applying a hidden Markov model (HMM) to model line spectrum extraction, but the amount of calculation was large and the probability of state transition was difficult to determine when processing actual signals. To solve the problem of the difficulty of extracting the line spectrum under a low signal-to-noise ratio, the large amount of calculation of the HMM algorithm and the difficulty of determining the probability of state transition, a line spectrum extraction method based on an HMM is proposed.

Basic theory of HMM:
HMM expression: The HMM [9,10] algorithm is a relatively mature algorithm that has been studied for many years. The state sequence of an HMM can be described as X k = (x 1 , x 2 , x 3 , …, x k ), and the observation sequence is Y k = (y 1 , y 2 , y 3 , …, y k ), where x k is the frequency number of the lofargram and y k can usually be expressed by the power of the lofargram; i.e. y k = P k .
The HMM model is usually represented by a triple λ = (π , A, B), where A is the state transition matrix, A = (a e (g)) M × M , M is the number of frequency points, and a e (g) is where g is the frequency number at the current moment, e is the frequency number at the previous moment, c is a normalization factor, R is the covariance matrix that determines the frequency distribution of the line spectrum, H is the state matrix, det indicates a determinant, and ζ is the signal-to-noise ratio at time k, σ is the standard deviation of noise, and I o represents the Bessel function. However, under normal circumstances, the signal-to-noise ratio is unknown, so b kj can be expressed as where P k,j represents the amplitude of the jth frequency at time k.
Forward-backward algorithm: The forward-backward (FB) algorithm is a local optimization algorithm that is the maximum posterior probability estimation of the state. The FB algorithm can therefore use the observation data Y k to estimate the possible statex k of the line spectrum, expressed asx The FB algorithm can implement the recursive calculation in the forward and backward directions using α and β, respectively. The forward recursion formula can be expressed as where i ∈ I, k = 2, 3, …, K.
The backward recursion formula can be expressed as Then, according to (5) and (6), the FB algorithm can be expressed as The line spectrum estimated according to the FB algorithm is Similarly, when there are multiple spectral lines, the forward recursion formula of the HMM algorithm evolves into where k = 2, 3, …, K. The backward recursion formula is At the same time, to eliminate the coupling between the forward and backward probabilities, their weights must be obtained: In Equation (12), ω α,l k and ω β,l k are the weighting factors of the forward and backward probabilities, respectively. To further eliminate the mutual influence of multiple spectral lines, the state probability around the Lth spectral line must be integrated to obtain where v α and v β are the numbers of the frequency points that must be integrated around the Lth spectral line in the forward and backward recursive algorithms. Substituting Equation (13) into Equation (11), the equation becomes so the FB algorithm can finally be expressed as The estimated value of the Lth spectral line is where a and b are integer constants and a ≥ b.
Improved algorithm proposed in the present paper: For the HMM algorithm proposed above, the following shortcomings exist.
1. Either in the recursion processing or in the state transition probability calculation, the calculations must be carried out point-by-point. When the sampling rate is high and the number of frequency points is large, the calculations can be large. 2. When processing real signals, the state matrix H in the state transition probability a e (g) [9] and the covariance matrix R that determines the line spectrum frequency are unknown and difficult to determine.
To address these issues, the following improvements were made.
1. When processing real signals, it was usually assumed that the line spectrum existed in a local maximum of the power spectrum. We therefore added a pre-processing step to the FB algorithm to extract the local maximum of the power at any given time. This could greatly reduce the amount of calculations of the processing algorithm. 2. In the calculation of the above state transition probability, many parameters could not be determined. Hence, a truncated Gaussian probability model was adopted. Then, Equation (1) became where is the line spectrum offset and δ is the maximum allowable line spectrum offset.
When the signal-to-noise ratio is low, the estimated line spectrum values might contain outliers. In such cases, the median value fitting method was used to eliminate the outliers. This process is shown in Figure 1 and proceeds as follows.
1. First, the HMM algorithm is used to obtain the estimated line spectrum value x k . 2. The median value M 0 of x k is calculated, as is the deviation d k between each x k and M 0 .

b) Comparison with results obtained using the forward-backward (FB) algorithm; (c) Comparison with results obtained with proposed algorithm
3. The maximum offset e is set and the deviation value d k which is greater than e is set to 0. 4. The value of k when the deviation d k = 0 is determined, and x k = 0 is made. 5. In the line spectrum estimation x k , the nearest two non-zero values x k1 and x k2 for each x k that is equal to 0 are found, and the mean value x k0 of x k1 and x k2 is calculated. 6. Finally, in the lofargram, the values of the n frequency points centered around x k0 at the current time were found and the position x kmax of the maximum among the values of the n frequency points was used as the estimated value of the line spectrum for that time.
Simulation and processing of sea trial data: Simulation: Single spectral line: The signal was a single frequency sine wave, the frequency was 400 Hz, the sampling frequency was 2048 Hz, the one-time processing time was 1 s, the total duration was 500 s, the noise was Gaussian white noise, and the signal-to-noise ratio was −20 dB ( Figure 2). Figures 2a and 2b show the original lofargram and the comparison with the line spectrum. Figure 2a shows that the signal-to-noise was relatively low and the line spectrum was quite weak. It can be seen from Figure 2b that the deviation between the true and estimated values obtained from the FB algorithm was very small, and the result was very good.
To verify the performance of the improved algorithm at an even lower signal-to-noise ratio and in the presence of outliers in the estimated line spectrum values, the signal-to-noise ratio was changed to −22 dB. As shown in Figure 3a, the line spectrum trajectory was then weaker. Figure 3b shows a comparison of the estimated and true values of the line spectrum after processing with the FB algorithm. As shown in the figure, the performance of the FB algorithm degraded at the lower signalto-noise ratio and some discrete outlier points appeared. Figure 3c shows the results after applying median fitting to the results obtained with the proposed FB algorithm. The figure also shows that the discrete outliers were successfully removed after median fitting.
Multiple spectral lines: The case of multiple spectral lines was then considered. The signal was again a single-frequency sinusoidal signal, including three spectral lines at frequencies of 100, 400, and 700 Hz. The sampling frequency was 2048 Hz, the one-time processing duration was 1 s, the total duration was 500 s, the noise was Gaussian white noise, and the signal-to-noise ratio was −22 dB ( Figure 4). Figure 4a shows the lofargram for a signal-to-noise ratio of −22 dB. The figure shows that the line spectrum was weak for the lower signalto-noise ratio. Figure 4b shows the estimated value of the line spectrum processed with the FB algorithm proposed in this paper. As can be seen from the figure, the line spectrum was weak and some discrete outlier points appeared. Figure 4c shows the estimated line spectrum obtained from the median fitting to the data processed with the FB algorithm proposed in this paper. The figure also shows that the discrete points were successfully removed and that the line spectrum estimation worked well.
Verification of sea trial data: The performance of the algorithm was verified with data received from the sea trial test. Because the actual processing required the filtering of the algorithm, only part of the frequency band was extracted. Figure 5a shows the spectrogram of the first segment of the processed data. Only the line spectrum in the 100-400 Hz frequency band was studied, and the time duration of the data was 100 s. The figure shows many spectral lines, including some indistinct, weak spectral lines, e.g. at 332 Hz. Figure 5b shows the spectrogram after processing with the algorithm proposed in this paper. The spectral lines were extracted mostly intact, including the weak spectral line at 332 Hz that was successfully extracted after the algorithm processing. Figures 6a and 6b show the original spectrogram of the second segment of data and the line spectrogram extracted with the proposed algorithm, respectively. As can be seen from the figures, the proposed algorithm performed well in extracting the line spectrum. The algorithm performed well even for a line spectrum with weak energy. Summary: To address the difficulty of line spectrum extraction under a low signal-to-noise ratio, an improved HMM-based line spectrum extraction algorithm is proposed. By improving the calculation method of the state transition probability and the FB recursion calculation, the proposed algorithm greatly reduced the amount of calculation and improved the practicality of the algorithm. To solve the problem of discrete outliers that tended to appear in line spectrum estimates at low signal-to-noise ratio, a median fitting method was proposed that was validated with simulation and sea trial data. The line spectrum extraction algorithm still performs well when the signal-to-noise ratio is as low as −22 dB.