CN105488341A - Denoising method based on hybrid EMD (Empirical Mode Decomposition) - Google Patents

Abstract

The invention discloses a denoising method based on hybrid EMD (Empirical Mode Decomposition), and belongs to the technical field of digital signal processing. According to the denoising method disclosed by the invention, EEMD (Ensemble Empirical Mode Decomposition) and MEMD (Multivariate Empirical Mode Decomposition) algorithms are combined, and aiming at a masking signal frequency used in the MEMD algorithm, parameter optimization is carried out, so that mode mixing can be more effectively eliminated, and a denoising effect is better; and according to the denoising method, targeted optimization is further carried on an additive noise amplitude in the EEMD, filtering thresholds of noise related IMF (Spurious modes) and filtering thresholds of noise IMF, so that an integral denoising effect is further improved. Compared with the prior art, the method disclosed by the invention can more effectively eliminate mode mixing, and has a better denoising effect.

Description

A kind of denoising method based on mixing empirical mode decomposition

Technical field

The present invention relates to a kind of denoising method, particularly relate to a kind of denoising method based on mixing empirical mode decomposition (EmpiricalModeDecomposition is called for short EMD), belong to digital signal processing technique field.

Background technology

EEMD (the EnsembleEmpiricalModeDecomposition of EMD and derivation, set empirical mode decomposition) and maskingsignalEMD (mask signal method empirical mode decomposition, after be called for short MEMD in literary composition) be the signal processing method grown up in the last few years.It decomposes based on signal time yardstick, is applicable to non-linear and unstable signal process, owing to not needing to determine the subjective experience optimum configurations such as the basis function Sum decomposition number of plies, in some cases, has better discomposing effect than wavelet transformation.However, EMD also has some problems and the defect of oneself.The problem of EMD mainly contains end effect, modal overlap.Although EEMD can solve the problem of EMD to a certain extent, it is uncertain for adding noise amplitude.EEMD in actual applications, especially when signal to noise ratio (S/N ratio) is lower, modal overlap problem still exists, and the noise that now EEMD adds can be used as the noise of signal as useful signal and decompose, thus brings modal overlap and more false IMF (Spuriousmodes).In addition, the amplitude of the ensemblenoise of EEMD is selected relevant with the signal to noise ratio (S/N ratio) of signal, and accurately select noise amplitude to have certain difficulty, great many of experiments shows, when signal to noise ratio (S/N ratio) is low, EMD decomposes will be advantageously.MEMD solves the problem of extreme's distribution by adding the high frequency mask signal determined, therefore for modal overlap, MEMD has obvious advantage and counting yield.The amplitude of the mask signal of MEMD and frequency parameter are selected to be a committed step, although give reference value in literary composition, but in actual applications, reference frequency is often lower than the requirement that maskingsignal is required, thus the ideal effect eliminating modal overlap can not be reached, therefore, for different signals, the result that the parameter of recommendation calculates non-optimal.IMF filtering mainly contains the interval filtering (EMD-IT) of wavelet filter, the direct filtering of EMD (EMD-DT) and EMD and related derivative filtering method, and wherein wavelet threshold system errors is comparatively large, is not easy to find accurately in practice; EMD-DT can not ensure the continuity of signal, for the signal filtering poor signal of multiple amplitude range; Although EMD-IT can ensure the continuity of signal, the filtering performance of several value signals of low order IMF can not be ensured.Generally speaking, EMD-IT filtering mainly contains two features, the first, after finding relative (noise be correlated with) IMF, directly abandon " noise " IMF of low order, not have consideration wherein with effective constituent; The second, filtering threshold is self-adaptation not, and single threshold strategies exists comparatively big error for several value signal filtering.

Summary of the invention

Technical matters to be solved by this invention is to overcome prior art deficiency, a kind of denoising method based on mixing empirical mode decomposition is provided, EEMD and MEMD algorithm combines by the method, and carry out parameter optimization for the mask signal frequency used in MEMD algorithm, more effectively can eliminate modal overlap, denoising effect is better.

Technical solution of the present invention is specific as follows:

Based on a denoising method for mixing empirical mode decomposition, comprise the following steps:

Step 1, utilize set empirical mode decomposition EEMD method to decompose original signal, obtain a series of natural mode of vibration component IMF and discrepance, and from all IMF, determine noise IMF, noise is correlated with IMF, useful IMF;

Step 2, threshold filter denoising is carried out to noise IMF IMF relevant with noise, and utilize mask signal method to carry out the elimination of modal overlap to the signal that useful IMF is formed; Wherein, the mask signal frequency that mask signal method uses is determined in accordance with the following methods: first carry out first derivation to the signal that useful IMF is formed, obtain first derivative signal; Then utilize Hilbert transform method to obtain the instantaneous frequency of described first derivative signal, then ask for the average of instantaneous frequency; Finally be multiplied by the average of described instantaneous frequency the weights that a span is (1,2), gained product is mask signal frequency;

Step 3, be correlated with to the noise IMF after threshold filter denoising, noise IMF, the signal that the useful IMF eliminating modal overlap is formed, and discrepance carries out signal reconstruction.

In order to improve the performance of EEMD, present invention further proposes following improvement project, the interpolation noise amplitude in EEMD be optimized:

The amplitude of the interpolation noise used in described EEMD method pre-determines by the following method: first carry out filtering to original signal, obtain filtered signal, and adds the noise of different amplitude respectively to original signal, obtains one group of noisy signal; Then respectively following process is carried out to each noisy signal: empirical mode decomposition is carried out to this noisy signal, and define with IMF from the IMF that empirical mode decomposition obtains; Calculate the summation of each useful IMF and the relevant root-mean-square error between filtered signal, and the summation of related coefficient often between adjacent two useful IMF; The amplitude of the interpolation noise finally used in EEMD method using the amplitude of the interpolation noise corresponding to the minimum noisy signal of two summation sums.

In order to improve the threshold filter effect of the IMF that to be correlated with to noise, present invention further proposes following improvement project:

The noise IMF that is correlated with is carried out the filtering threshold that threshold filter denoising uses and determined in accordance with the following methods: the IMF that is first correlated with by each noise carries out decorrelative transformation with noise IMF respectively; Be divided into two parts according to the positive and negative IMF that to be correlated with by noise after each decorrelative transformation of amplitude, its amplitude average and energy filter threshold value are obtained respectively to every part, and will both weighted sums, obtain the filtering threshold of this part.

In order to more effectively extract the useful information in noise IMF, the present invention proposes following improvement project further:

Carry out to noise IMF the filtering threshold that threshold filter denoising uses to determine in accordance with the following methods: the amplitude average of calculating noise IMF, and time-sequencing is carried out to each amplitude being greater than amplitude average in noise IMF, obtain the time series of this amplitude; Time series for obtained each amplitude carries out periodic authentication respectively, and therefrom selects periodically higher part-time sequence; Finally using the minimum value of amplitude corresponding to selected part-time sequence as described filtering threshold.

Compared to existing technology, the present invention and further improvement opportunity scheme thereof have following beneficial effect:

EEMD and MEMD algorithm of the present invention combines, and carries out parameter optimization for the mask signal frequency used in MEMD algorithm, and more effectively can eliminate modal overlap, denoising effect is better;

The present invention is optimized for the interpolation noise amplitude in EEMD further, using the RelativeRMSE value of useful IMF and the minimum choice criteria of adding noise amplitude as EEMD of related coefficient sum, thus more effectively can solve the problem that end effect brings;

The present invention is optimized for the be correlated with filtering threshold of IMF of noise further, adopts positive and negative signal amplitude separately strategy, by decorrelation, gets average and determine last filtering threshold according to energy three aspects, better to the be correlated with filter effect of IMF of noise;

The present invention utilizes the useful information in the method extraction noise IMF of periodic authentication further, avoids the loss directly abandoning brought useful information, further increases denoising effect.

Accompanying drawing explanation

Fig. 1 determines that EEMD adds the schematic flow sheet of noise amplitude method;

Fig. 2 is the schematic flow sheet of MEMD mask signal frequency determination methods;

Fig. 3 is that noise is correlated with the filtering threshold defining method schematic flow sheet of IMF.

Embodiment

Below in conjunction with accompanying drawing, technical scheme of the present invention is described in detail:

Basic ideas of the present invention are combined by EEMD and MEMD algorithm, and carry out parameter optimization for the mask signal frequency used in MEMD algorithm, and more effectively can eliminate modal overlap, denoising effect is better; The present invention optimizes targetedly the be correlated with filtering threshold of IMF, the filtering threshold of noise IMF of the interpolation noise amplitude in EEMD, noise, thus further increases global de-noising effect further.

For the ease of public understanding, with a preferred embodiment, technical solution of the present invention is described in detail below.The present invention is based on the denoising method of mixing empirical mode decomposition, specifically comprise the following steps:

Step 1, utilize set empirical mode decomposition EEMD method to decompose original signal, obtain a series of natural mode of vibration component IMF and discrepance, and from all IMF, determine noise IMF, noise is correlated with IMF, useful IMF.

Traditional E EMD method can be adopted to the decomposition of original signal in this step or improve EEMD method.EEMD is prior art, its detailed content can see document [Wu.Zhaohua, HuangNE, Ensembleempiricalmodedecomposition:anoiseassisteddataana lysismethod, Adv.Adapt.DataAnal.1 (1) (2009) 1 – 41.].

As a kind of innovatory algorithm of EMD, although EEMD can solve the problem of EMD to a certain extent, it adds noise amplitude is uncertain.EEMD in actual applications, especially when signal to noise ratio (S/N ratio) is lower, modal overlap problem still exists, and the noise that now EEMD adds can be used as the noise of signal as useful signal and decompose, thus brings modal overlap and more false IMF (Spuriousmodes).In addition, the amplitude of the ensemblenoise of EEMD is selected relevant with the signal to noise ratio (S/N ratio) of signal, accurately selects noise amplitude to have certain difficulty.For this reason, the present invention improves existing EEMD method, interpolation noise amplitude is optimized, using the RelativeRMSE value of useful IMF and the minimum choice criteria of adding noise amplitude as EEMD of related coefficient sum, thus more effectively can solves the problem that end effect brings.Particularly, as shown in Figure 1, the amplitude of the interpolation noise used in described EEMD method pre-determines by the following method:

First rough filtering (for the sake of simplicity, preferably adopting median filter method) is carried out to original signal, obtain filtered signal, and the noise of different amplitude is added respectively to original signal, obtain one group of noisy signal; Preferably, add the amplitude interval of noise for 0.001-0.4 original signal standard deviation doubly;

Then respectively following process is carried out to each noisy signal: EMD decomposition is carried out to this noisy signal, and decompose the IMF obtained from EMD and define with IMF; Calculate the summation of each useful IMF and the relevant root-mean-square error RelativeRMSE between filtered signal, and the summation of related coefficient often between adjacent two useful IMF; The amplitude of the interpolation noise finally used in EEMD method using the amplitude of the interpolation noise corresponding to the minimum noisy signal of two summation sums (representing this parameter with Rel-CorrRMSE).The mathematical expression of Rel-CorrRMSE is specific as follows:

$\mathrm{Re}l-CorrRMSE=\sqrt{\frac{\underset{i=1}{\overset{n}{\Σ}}\underset{k=1}{\overset{N}{\Σ}}{\left({\hat{x}}_0\right(k)-c_{i,rel}(k\left)\right)}^2}{\underset{k=1}{\overset{N}{\Σ}}{\hat{x}}_0^2\left(k\right)}}+\underset{i=1}{\overset{n-1}{\Σ}}corr\left(c_{i,rel}\right(k),c_{i+1,rel}(k\left)\right)---\left(1\right)$

Wherein, the signal of original signal after mean filter, c _{i, rel}k () is i-th useful IMF, corr () expression ask related coefficient to the signal of two in bracket.

Determine add noise amplitude after, EEMD can be utilized to decompose original signal, and determine noise IMF, noise is correlated with IMF and useful IMF.

Step 2, threshold filter denoising is carried out to noise IMF IMF relevant with noise, and utilize mask signal method to carry out the elimination of modal overlap to the signal that useful IMF is formed; Wherein, the mask signal frequency that mask signal method uses is determined in accordance with the following methods: first carry out first derivation to the signal that useful IMF is formed, obtain first derivative signal; Then utilize Hilbert transform method to obtain the instantaneous frequency of described first derivative signal, then ask for the average of instantaneous frequency; Finally be multiplied by the average of described instantaneous frequency the weights that a span is (1,2), gained product is mask signal frequency.

Utilize the MEMD algorithm of improvement to carry out the elimination of modal overlap to the signal that useful IMF is formed in the present embodiment, and utilize the filtering threshold optimized to carry out threshold filter to noise IMF IMF relevant with noise respectively.Be described respectively below.

(1) optimization of mask signal frequency

In actual applications, reference frequency, often lower than the requirement that maskingsignal is required, thus can not reach and eliminate the ideal effect of modal overlap existing MEMD algorithm, therefore, for different signals, and the result that the parameter of recommendation calculates non-optimal.For this reason, the present invention improves existing MEMD algorithm, promotes by optimizing mask signal frequency the effect eliminating modal overlap.Particularly, as shown in Figure 2, the mask signal frequency that mask signal method uses is determined in accordance with the following methods:

First first derivation is carried out to the signal that useful IMF is formed, obtain first derivative signal; Suppose that the signal that useful IMF is formed is x ₁=x ₂+ x ₃, x in composition ₂, x ₃corresponding frequency is f ₁, f ₂, and hypothesis f ₂>f ₁, then have:

x ₁′＝x′ ₂+x ₃′(2)

First order derivative is asked in apostrophe representative above.

Then utilize Hilbert transform method to obtain the instantaneous frequency of described first derivative signal, and ask for the average of instantaneous frequency; Its mathematical expression is as follows:

$h\[\left(x_1^{\′}\right)\]=\frac{1}{\mathrm{\π}}P{\∫}_{-\mathrm{\∞}}^{+\mathrm{\∞}}\frac{x_2^{\′}\left(\mathrm{\τ}\right)+x_3^{\′}\left(\mathrm{\τ}\right)}{t-\mathrm{\τ}}d\mathrm{\τ}---\left(3\right)$

$\stackrel{\‾}{f}=\frac{\underset{i=1}{\overset{k}{\Σ}}{f}_1^{\′}\left(i\right)}{k}\≈f_2---\left(4\right)$

Wherein, H [(x ₁')] be exactly ask for Hilbert conversion, represent that Hilbert converts the average of the instantaneous frequency of trying to achieve, f ₁' (i) represents instantaneous frequency, f ₂represent the frequency values of the composition that frequency is the highest in aliasing signal.

Finally to the average of instantaneous frequency be multiplied by the weights that a span is (1,2), gained product is mask signal frequency f _mask, that is:

$1.1\stackrel{\&OverBar;}{f}<f_{mask}<2\stackrel{\&OverBar;}{f}---\left(5\right)$

(2) noise is correlated with the filtering threshold optimization of IMF

Noise to be correlated with the threshold value problem identificatioin of IMF, to consider the asymmetry of signal in noise IMF, need the positive and negative separate computations to original signal point amplitude.First noise is correlated with IMF (IMF _kk=2,3...) and noise IMF (IMF ₁) carry out decorrelation operation, remove part noise, as shown in formula (6).Then the average of signal amplitude is asked for, as shown in formula (7), as comprehensively measuring of noise and useful signal; Again according to energy threshold formula determine corresponding energy threshold, finally process is weighted to average threshold value and energy threshold, obtains final noise threshold, as formula (7), (8).

${\mathrm{IMF}}_K=1{\mathrm{MF}}_K-{\mathrm{corr}}_{{\mathrm{IMF}}_K-{\mathrm{IMF}}_1}\&CenterDot;{\mathrm{IMF}}_1---\left(6\right)$

$T_{mean}=\frac{median\left(\right|c_i\left|\right)}{N}---\left(7\right)$

${\hat{E}}_k=\frac{E_1^2}{\mathrm{\&beta;}}{\mathrm{\&rho;}}^{-k},k=2,3,4,\mathrm{...};T_{\widehat{\mathrm{\&sigma;}}}=C\&CenterDot;\sqrt{{\hat{E}}_{k}2ln\left(N\right)}---\left(8\right)$

$T_{uni}=a_1\&CenterDot;T_{mean}+a_2\&CenterDot;T_{\widehat{\mathrm{\&sigma;}}}---\left(9\right)$

Wherein, IMF _kfor noise is correlated with IMF, IMF ₁for noise IMF, | c _i| be the absolute value of IMF each point amplitude, T _meanfor mean filter threshold value, N is signal number. for IMF ₁energy square, for the energy estimators of a kth IMF, β, ρ, C are constant, β=0.719, ρ=2.01, C ∈ [0.1,1.4], T _unifor the filtering threshold of final utilization.Two weights a during weighted sum ₁, a ₂usually be all 0.5.

(3) the filtering threshold optimization of noise IMF

Noise IMF normally directly abandons by prior art, likely loses useful information wherein.Utilize the method determination filtering threshold of periodic authentication in the present embodiment, fully to extract the useful information in noise IMF for this reason.Particularly, carry out to noise IMF the filtering threshold that threshold filter denoising uses to determine in accordance with the following methods: the amplitude average of calculating noise IMF, and time-sequencing is carried out to each amplitude being greater than amplitude average in noise IMF, obtain the time series of this amplitude; Time series for obtained each amplitude carries out periodic authentication respectively, and therefrom selects periodically higher part-time sequence; Finally using the minimum value of amplitude corresponding to selected part-time sequence as described filtering threshold.

Periodic authentication method in the present embodiment is specific as follows: to the time series of this amplitude, first the in chronological sequence order arrangement of the mistiming of adjacent amplitude, obtains very first time difference sequence, and calculates the average k of these mistimings; Then difference done to adjacent time difference in very first time difference sequence again and by difference in chronological sequence order arrangement, obtain the second mistiming sequence; To each data in the second mistiming sequence, its business divided by average k gained and predetermined threshold value α are compared, be less than α, assignment is 1, otherwise assignment is 0, thus obtains a binary sequence; Calculate the Lempel-Ziv complexity of this binary sequence, complexity is lower, then the seasonal effect in time series of the amplitude corresponding to this binary sequence is periodically higher.

Its algorithm is expressed as follows: supposing that the time series of similar magnitude does poor the second mistiming sequence obtained through twice is { c (n ₁), c (n ₂) ..., c (n _m), by following formula, binary conversion treatment is carried out to the second mistiming sequence, obtain a new binary sequence { Bc (n ₁), Bc (n ₂) ..., Bc (n _m) }:

$BC\left(n_i\right)=\left\{\begin{array}{cc}1,& c\left(n_i\right)/k\&GreaterEqual;\mathrm{\&alpha;}\\ 0,& c\left(n_i\right)/k<\mathrm{\&alpha;}\end{array}\right.,$

Then binary sequence { Bc (n is calculated ₁), Bc (n ₂) ..., Bc (n _m) Lempel-Ziv complexity C _lZN, and judge whether the time series of this amplitude has periodically according to Lempel-Ziv complexity:

$\left\{\begin{array}{cc}non-period& C_{LZN}\&GreaterEqual;\mathrm{\&beta;}\\ period& C_{LZN}<\mathrm{\&beta;}\end{array}\right.$

Wherein, c (n _i) be i-th element in the second mistiming sequence; Bc (n _i) be the element of i-th in binary sequence; C _lZNfor Lempel-Ziv complexity; β is the threshold value of Lempel-Ziv complexity, is to judge whether signal has periodic final step, generally gets the value (the smaller the better) of less than 0.5.

Step 3, be correlated with to the noise IMF after threshold filter denoising, noise IMF, the signal that the useful IMF eliminating modal overlap is formed, and discrepance carries out signal reconstruction.

Claims (7)

1., based on a denoising method for mixing empirical mode decomposition, it is characterized in that, comprise the following steps:

Step 1, utilize set empirical mode decomposition EEMD method to decompose original signal, obtain a series of natural mode of vibration component IMF and discrepance, and from all IMF, determine noise IMF, noise is correlated with IMF, useful IMF;

Step 2, threshold filter denoising is carried out to noise IMF IMF relevant with noise, and utilize mask signal method to carry out the elimination of modal overlap to the signal that useful IMF is formed; Wherein, the mask signal frequency that mask signal method uses is determined in accordance with the following methods: first carry out first derivation to the signal that useful IMF is formed, obtain first derivative signal; Then utilize Hilbert transform method to obtain the instantaneous frequency of described first derivative signal, then ask for the average of instantaneous frequency; Finally be multiplied by the average of described instantaneous frequency the weights that a span is (1,2), gained product is mask signal frequency;

Step 3, be correlated with to the noise IMF after threshold filter denoising, noise IMF, the signal that the useful IMF eliminating modal overlap is formed, and discrepance carries out signal reconstruction.

2. denoising method as claimed in claim 1, it is characterized in that, the amplitude of the interpolation noise used in described EEMD method pre-determines by the following method: first carry out filtering to original signal, obtain filtered signal, and original signal is added respectively to the noise of different amplitude, obtain one group of noisy signal; Then respectively following process is carried out to each noisy signal: empirical mode decomposition is carried out to this noisy signal, and define with IMF from the IMF that empirical mode decomposition obtains; Calculate the summation of each useful IMF and the relevant root-mean-square error between filtered signal, and the summation of related coefficient often between adjacent two useful IMF; The amplitude of the interpolation noise finally used in EEMD method using the amplitude of the interpolation noise corresponding to the minimum noisy signal of two summation sums.

3. denoising method as claimed in claim 2, is characterized in that, is 0.001-0.4 original signal standard deviation doubly to the amplitude range of the noise of the different amplitudes that original signal is added.

4. denoising method as claimed in claim 1, is characterized in that, the noise IMF that is correlated with is carried out to the filtering threshold that threshold filter denoising uses and determines in accordance with the following methods: the IMF that is first correlated with by each noise carries out decorrelative transformation with noise IMF respectively; Be divided into two parts according to the positive and negative IMF that to be correlated with by noise after each decorrelative transformation of amplitude, its amplitude average and energy filter threshold value are obtained respectively to every part, and will both weighted sums, obtain the filtering threshold of this part.

5. denoising method as claimed in claim 4, it is characterized in that, when being weighted summation to amplitude average and energy filter threshold value, the weights of the two are 0.5.

6. denoising method as claimed in claim 1, it is characterized in that, carry out to noise IMF the filtering threshold that threshold filter denoising uses to determine in accordance with the following methods: the amplitude average of calculating noise IMF, and time-sequencing is carried out to each amplitude being greater than amplitude average in noise IMF, obtain the time series of this amplitude; Time series for obtained each amplitude carries out periodic authentication respectively, and therefrom selects periodically higher part-time sequence; Finally using the minimum value of amplitude corresponding to selected part-time sequence as described filtering threshold.

7. denoising method as claimed in claim 6, it is characterized in that, for the time series of any one amplitude, described periodic authentication is specific as follows: to the time series of this amplitude, first the in chronological sequence order arrangement of the mistiming of adjacent amplitude, obtain very first time difference sequence, and calculate the average k of these mistimings; Then difference done to adjacent time difference in very first time difference sequence again and by difference in chronological sequence order arrangement, obtain the second mistiming sequence; To each data in the second mistiming sequence, its business divided by average k gained and predetermined threshold value α are compared, be less than α, assignment is 1, otherwise assignment is 0, thus obtains a binary sequence; Calculate the Lempel-Ziv complexity of this binary sequence, complexity is lower, then the seasonal effect in time series of the amplitude corresponding to this binary sequence is periodically higher.