WO2021042589A1

WO2021042589A1 - Extremum energy decomposition-based electrocardiogram signal analysis method

Info

Publication number: WO2021042589A1
Application number: PCT/CN2019/121414
Authority: WO
Inventors: 周作建; 宁新宝; 曾彭; 姜晓东; 王�华; 刘红星; 王斌斌
Original assignee: 江苏华康信息技术有限公司
Priority date: 2019-09-06
Filing date: 2019-11-28
Publication date: 2021-03-11
Also published as: CN110558974A; CN110558974B

Abstract

An extremum energy decomposition-based electrocardiogram signal analysis method, comprising: acquiring an ECG signal x(t) in an unknown state at a given time and at a given sampling frequency; performing de-noising preprocessing on the ECG signal x(t); taking the de-noised ECG signal x(t) as an original signal, decomposing the original signal x(t) into n extremum mode function components and a margin, the n extremum mode function components, obtained by decomposing the original signal x(t), representing components of different frequency bands of the original signal; and determining, according to the n extremum mode function components, whether the ECG signal is an abnormal electrocardiogram signal. The ECG signal is analyzed by means of extremum energy decomposition, the original signal is decomposed into a plurality of components, i.e. extremum component functions, and the energy of each component is calculated to obtain the energy distribution thereof.

Description

An ECG signal analysis method based on extreme energy decomposition method

Technical field

The invention relates to an electrocardiogram signal analysis, in particular to an electrocardiogram signal analysis method based on an extreme value energy decomposition method.

Background technique

Physiological signals are generated by the interaction of multiple systems in a living body, and the time and intensity of the action of different systems are different, leading to the complexity of time and space in physiological signals. The ECG (electrocardiogram) signal reflects the electrical activity process of the heart beat, and has important guiding significance for the basic functions of the heart and the diagnosis of diseases. The ECG signal includes several different wave groups, such as P wave, QRS wave, and T wave. Each wave group contains different frequency components, and the energy distribution ratios of different frequency components in the ECG signal are different. Studies have shown that 99% of the energy of the ECG signal is concentrated in the range of 0 to 40 Hz. For different wave groups, the frequency distribution range and energy ratio have certain rules. Heart disease can cause changes in the energy distribution of ECG signals. Studying the energy distribution of ECG has important clinical guiding significance for revealing the functional changes caused by heart diseases.

In the energy (spectrum) research of ECG, there are many traditional classical methods, such as spectrum analysis, time domain analysis, etc. However, ECG signals are non-stationary and nonlinear signals. Traditional frequency domain analysis methods are more suitable for analyzing stationary signals and can only Gives global frequency information. In addition, when segmenting the frequency domain, the traditional method often divides according to several fixed frequency points, without considering the fluctuation characteristics of different signals and the differences between the signals.

Therefore, it is urgent to solve the above-mentioned problems.

Summary of the invention

Object of the invention: The object of the present invention is to provide an ECG signal analysis method based on extreme energy decomposition method that can use less data and can intuitively reflect the true law of ECG energy distribution and signal fluctuation.

Technical solution: In order to achieve the above objectives, the present invention discloses an ECG signal analysis method based on extreme energy decomposition method, which includes the following steps:

(1) Obtain the ECG signal x(t) in an unknown state at a given time and a given sampling frequency;

(2) Perform denoising preprocessing of the ECG signal x(t); the specific preprocessing method is: pass the ECG signal through a 40Hz zero-phase FIR low-pass filter to remove high-frequency noise, and then pass a median filter to remove the baseline drift;

(3) Take the denoised ECG signal x(t) as the original signal, find all the local extreme points of the original signal, and then use the spline curve for all the maximum points and all the minimum points of the original signal They are connected to form the upper envelope e _max and the lower envelope e _{min respectively} , and the envelope mean signal m(t)=(e _max +e _min )/2 of the upper and lower envelopes is obtained;

(4) Subtract the envelope mean signal m(t) from the original signal x(t) to obtain h(t)=x(t)-m(t); then judge whether h(t) satisfies the extreme mode If the judgment condition of the function is not met, return h(t) as the original signal to step (3) until h _k (t) meets the judgment condition of the extreme modal function, then write c ₁ (t) = h _k (t), as the first extreme value modal function component;

(5) Subtract the first extremum modal function component c ₁ (t) from the _{original signal x(t) to obtain the margin r 1} (t) = x(t)-c ₁ (t), and then judge Whether h _k (t) meets the stopping criterion, if not, take r ₁ (t) as the new original sequence x(t), and return to steps (3) and (4) until h _k (t) meets the stopping criterion , Get the 2, 3,..., n extremum modal function components and margin r _n (t), then decompose the original signal x(t) into n extremum modal function components and a margin, namely

(6) Perform frequency spectrum analysis on the extreme value modal function components c _i (t), i=1, 2,..., n to obtain the center frequency of each extreme value modal function component;

(7) The n extremum modal function components obtained by decomposing the original signal x(t) represent the components of different frequency bands of the original signal, and then calculate the energy of each component

E _i ＝∫|c _i (t)| ² dt,i＝1,2,...,n

Normalize each energy value to get a normalized energy distribution vector

p _i ＝E _i /E,i＝1,2,...,n

among them,

The first component p ₁ represents the energy of the highest frequency band, which represents the proportion of the energy distribution of the signal in the highest frequency range, and the last component p _n represents the proportion of the energy distribution of the signal in the lowest frequency range; according to the normalized energy distribution vector Draw a normalized energy distribution map, where the abscissa represents the component level, the ordinate represents the normalized energy distribution vector value, the curve represents the average value, and the error bar represents the standard deviation;

Probability on each extremum mode function component normalization energy distribution vector p _i and a significance detector (8), the normalized unknown state ECG signal energy distribution vector p _i and the standard ECG signals to obtain Value P _i , judge whether the probability value P _{1 on the} first extreme value modal function component and the probability value P ₂ on the second extreme value modal function component are less than the probability standard value, if not, return to step (1) again Get the signal;

(9) Under the condition that P ₁ and P _{2 are} less than the probability standard value, if the energy distribution vector p _{1 of the} first extreme value modal function component is smaller than the energy distribution vector p of the second extreme value modal function component ₂ , it is determined that the ECG signal is an abnormal electrocardiogram signal.

Wherein, the minimum amount of data required for the original signal x(t) is N=2 ⁿ⁺¹ , where n is the number of decomposed extreme value modal function components.

Preferably, the determination condition of the extreme value modal function in the step (4) is: (a). In the entire data sequence, the number of extreme points is equal to or one difference from the number of zero crossing points; (b). At any moment, the upper and lower envelopes are symmetrical to the time axis.

_{Furthermore, the formula for h k} (t) satisfying the stopping criterion in the step (5) is:

ε represents the screening threshold, which is between 0.2 and 0.3.

Further, the probability standard value in the steps (8) and (9) is 0.05.

Furthermore, the repetition frequency of the waveform of the sixth extremum modal function component obtained in the step (5) is the human heart rate, that is, the cardiac cycle.

Beneficial effects: Compared with the prior art, the present invention has the following significant advantages: the present invention uses an extreme energy decomposition method (EED) to analyze the ECG signal, and decomposes the original signal into multiple components, that is, extreme components Function, calculate the energy of each component, obtain its energy distribution, and determine whether the signal is normal or not according to the comparison of the energy values of the components at different levels; the present invention can decompose the signal into subordinates according to the fluctuation law of the biomedical signal itself Signals with different time levels from high frequency to low frequency; and for ECG signals, the repetition frequency can be determined from a special level, that is, the heart rate value of the ECG signal (including the instantaneous heart rate value and the average heart rate value). The method of obtaining the heart rate value is better than The traditional method is more accurate; the data length obtained by extreme value decomposition at all levels is the same, so the data length will not be reduced, so that it can be used for short-term data analysis, that is, it requires a small amount of data to analyze and get accurate results ; EED is not easy to be affected by noise for different levels of component energy analysis.

Description of the drawings

Figure 1 is a schematic diagram of the original signal in the present invention;

Fig. 2 is a schematic diagram of obtaining an envelope of an original signal in the present invention;

Fig. 3 is a schematic diagram of the original signal minus the envelope mean value signal in the present invention;

Fig. 4 is a schematic diagram of the first extremum modal function component obtained in the present invention;

Figure 5 is a schematic flow chart of the extreme value energy decomposition method in the present invention;

6 is a schematic diagram of EED decomposition of ECG signals in Embodiment 1 of the present invention;

FIG. 7 is a spectrum diagram of the ECG signal in Embodiment 1 of the present invention;

8 is a schematic diagram of EED decomposition of ECG signals in Embodiment 2 of the present invention;

FIG. 9 is a time-frequency diagram of ECG of a healthy person in implementation 2 of the present invention;

Fig. 10 is a time-frequency diagram of the ECG of CHF patients in the second embodiment of the present invention.

detailed description

The technical scheme of the present invention will be further described below in conjunction with the accompanying drawings.

As shown in Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5, an ECG signal analysis method based on the extreme value energy decomposition method of the present invention includes the following steps:

(1) Get the ECG signal x(t) of unknown state at a given time and a given sampling frequency; the minimum amount of data required for the original signal x(t) is N=2 ⁿ⁺¹ , where n is the decomposed pole The number of value modal function components;

(4) Subtract the envelope mean signal m(t) from the original signal x(t) to obtain h(t)=x(t)-m(t); then judge whether h(t) satisfies the extreme mode If the judgment condition of the function is not met, return h(t) as the original signal to step (3) until h _k (t) meets the judgment condition of the extreme modal function, then write c ₁ (t) = h _k (t), as the first extreme value modal function component; among them, the judging condition of the extreme value modal function is: (a). In the whole data sequence, the number of extreme points is equal to the number of zero crossing points or one difference ; (B). At any moment, the upper and lower envelopes are symmetrical to the time axis;

The formula for h _k (t) to satisfy the stopping criterion is:

ε represents the screening threshold, which is between 0.2 and 0.3; the extremum modal decomposition that meets the stopping criterion meets the following two conditions: (a) The final extremum modal function component c _n (t) or margin r _n (t) is less than the preset threshold; (b) the residual signal r _n (t) becomes a monotonic signal, and the extreme modal function signal cannot be extracted from it;

E _i ＝∫|c _i (t)| ² dt,i＝1,2,...,n

Normalize each energy value to get a normalized energy distribution vector

p _i ＝E _i /E,i＝1,2,...,n

among them,

Probability on each extremum mode function component normalization energy distribution vector p _i and a significance detector (8), the normalized unknown state ECG signal energy distribution vector p _i and the standard ECG signals to obtain Value P _i , judge whether the probability value P _{1 on the} first extreme value modal function component and the probability value P ₂ on the second extreme value modal function component are less than the probability standard value, if not, return to step (1) again Obtain the signal; the standard value of probability is selected as 0.05;

The repetition frequency of the sixth extreme value modal function component waveform obtained by the present invention is the human heart rate, that is, the cardiac cycle.

The extreme energy decomposition method (Extremum Energy Decomposition, EED) adopted by the present invention is a method based on the concept of extreme value modal function, which is a type of signal with a single frequency that satisfies the following two conditions at the same time , The two conditions are:

(a) In the entire data sequence, the number of extreme points (including maximum and minimum) and the number of zero-crossing points must be equal or differ by at most one;

(b) At any time, the upper envelope formed by the local maximum point and the lower envelope formed by the local minimum point have a mean value of zero, that is to say, the local upper and lower envelopes are locally symmetric with respect to the time axis;

The above two conditions, condition (a) is similar to the requirement of Gaussian normal stationary process for traditional narrowband, condition (b) guarantees that the instantaneous frequency calculated by the extreme value modal function has physical meaning.

The standard selection for the termination of the extreme value modal function decomposition of the present invention should be moderate. Too strict conditions will cause the last few extreme value modal function components to lose meaning; too loose conditions will cause useful components to be lost; in practical applications, It is also possible to set the number of extreme modal function component levels that need to be decomposed according to requirements, and the calculation is terminated when the decomposition level is satisfied.

Example 1

An ECG signal analysis method based on the extreme value energy decomposition method of the present invention includes the following steps:

(1) Obtain the ECG signal x(t) at a given time and a given sampling frequency from the normal sinus database nsrdb of physionet of healthy people; the data length is 8s, and the minimum amount of data required for the original signal x(t) is N=2 ⁿ⁺¹ =2 ⁹ , where n is the number of components of the extremum modal function decomposed, n=8;

(2) Perform denoising preprocessing of the ECG signal x(t); the specific preprocessing method is: Since the ECG energy is mainly concentrated in 0-40Hz, the ECG signal is filtered through a 40Hz zero-phase FIR low-pass filter to eliminate high frequency Noise, and then pass the median filter to remove baseline drift;

(5). Subtract the first extremum modal function component c ₁ (t) from the _{original signal x(t) to get the margin r 1} (t) = x(t)-c ₁ (t), and replace r ₁ (t) as the new original sequence x(t), return to steps (3) and (4) to obtain the second, third, ..., eighth extreme value modal function components and the margin r ₈ (t), Then the original signal x(t) is decomposed into 8 extreme value modal function components and a margin, namely

As shown in Figure 6, according to the characteristics of the waveform, the ECG signal can be divided into three wave groups: P, QRS, and T waves; according to the decomposition of 8 extreme value modal function components, c ₁ (t) and c ₂ (t ) With the highest frequency represents the decomposition component of the QRS complex with the highest frequency in the ECG signal, c ₃ (t) represents the decomposition component of the high frequency QRS complex and the P wave in the _{ECG signal, and c 4} (t) represents the high frequency in the ECG signal The decomposed components of the superposition of QRS complex, P wave and T wave, c ₅ (t) represents the decomposed component of the low-frequency part of QRS complex, P wave and T wave in the ECG signal, and c ₆ (t) represents the heartbeat rhythm The cardiac cycle, c ₇ (t) and c ₈ (t) represent the cardiac physiological adjustment rhythm on a larger time scale representing the long-term rhythm of the heart; the amplitude of the observed signal, the highest frequency component has a higher amplitude Value, the energy is the highest; the component with the lowest frequency has the lower amplitude and lower energy;

(6) For the extreme value modal function components c _i (t), i = 1, 2, ..., 8, perform spectrum analysis to obtain the center frequency of each extreme value modal function component, as shown in Table 1; where c ₁ Perform spectrum analysis to obtain the spectrum diagram shown in Figure 7. It can be concluded that _{the center frequency of c 1} is about 20 Hz, and the main frequency is concentrated in the range of 15 to 25 Hz; as the existing research shows, the spectrum range of P wave is 0 to 18 Hz ( ±3Hz), the energy is mainly concentrated in 5～12Hz; the spectrum range of QRS wave is 0～37Hz (±5Hz), the energy is mainly concentrated in 6～18Hz; the spectrum range of T wave is 0～8Hz (±2Hz), the energy is mainly Focus on 0～8Hz. Comparing Table 1, it can be seen that the frequency band of the QRS complex mainly contains _{two components, c 1} and c ₂ , the P wave mainly contains _{two components c 3} and c ₄ , and the T wave mainly contains the components _{c 4} ～c ₈ section. It should be noted that the inclusion here does not mean that each component is only determined by a specific ECG wave group (P, QRS, T wave), or that each ECG wave group is only contained in a specific component. The upper wave group ——Energy relationship is a major correspondence, not all. _{For example, the c 5} to c ₈ components representing the low frequency part are the result of the superposition of the low frequency parts of each ECG wave group, rather than a specific ECG wave group. The extreme value modal function components of the present invention and the extreme value modal function components of each level all show that the extreme value modal function components of ECG can represent a certain ECG wave group fluctuation situation, and reflect the fluctuation law of ECG at different levels; Compared with the traditional frequency domain analysis method, the EED method can directly observe the fluctuation of ECG at various levels, which is very intuitive. The repetition frequency of the waveform of the sixth extremum modal function component obtained is the human heart rate, that is, the cardiac cycle.

Table 1 The center frequency of each extreme value modal function component

Example 2

The EED analysis method is used to analyze the energy distribution of ECG in healthy people and CHF patients at different levels.

An ECG signal analysis method based on extreme energy decomposition method for healthy people includes the following steps:

(1). Obtain the ECG signal x(t) with a data length of 10s and a sampling frequency of 128Hz from the normal sinus database nsrdb of physionet; the nsrdb database contains 18 healthy people (age 34.3±8.4), the original signal x(t ) The minimum amount of data required is N=2 ⁿ⁺¹ =2 ⁹ , where n is the number of decomposed extreme value modal function components, n=8;

(6) Perform spectrum analysis on the extreme value modal function components c _i (t), i = 1, 2, ..., 8, and obtain the center frequency of each extreme value modal function component, and obtain the frequency domain analysis result diagram;

(7). The 8 extreme modal function components decomposed from the original signal x(t) represent the components of the original signal in different frequency bands, and then the energy of each component is calculated

E _i ＝∫|c _i (t)| ² dt,i=1，2，…，8

Normalize each energy value to get a normalized energy distribution vector

p _i ＝E _i /E,i＝1,2,…,8

among them,

The first component p ₁ represents the energy of the highest frequency band, and represents the proportion of the energy distribution of the signal in the highest frequency range, and the last component p _n represents the proportion of the energy distribution of the signal in the lowest frequency band.

An ECG signal analysis method based on extreme energy decomposition method for CHF patients includes the following steps:

(1) Get the ECG signal x(t) with a data length of 10s and a sampling frequency of 128Hz from the CHF database chfdb. The chfdb database contains 15 CHF patients (age 58.8±9.1), and the original signal x(t) requires the least The amount of data N=2 ⁿ⁺¹ =2 ⁹ , where n is the number of decomposed extreme value modal function components, n=8;

((2), the ECG signal x(t) is subjected to denoising preprocessing; the specific method of preprocessing is: Since the ECG energy is mainly concentrated in 0-40Hz, the ECG signal is filtered through a 40Hz zero-phase FIR low-pass filter to eliminate high Frequency noise, and then pass the median filter to remove baseline drift;

(6) Perform spectrum analysis on the extreme value modal function components c _i (t), i = 1, 2, ..., 8, and obtain the center frequency of each extreme value modal function component, and obtain the frequency domain analysis result graph, such as As shown in Figure 9, the energy distribution of healthy people is in the range of 0-40Hz, and the high-frequency energy above 20Hz has a higher proportion. As shown in Figure 10, the CHF energy is mainly concentrated below 20Hz, and the high-frequency energy is significantly reduced ；

E _i ＝∫|c _i (t)| ² dt,i=1，2，…，8

Normalize each energy value to get a normalized energy distribution vector

p _i ＝E _i /E,i＝1,2,…,8

among them,

The first component p ₁ represents the energy of the highest frequency band, which represents the proportion of the energy distribution of the signal in the highest frequency range, and the last component p _n represents the proportion of the energy distribution of the signal in the lowest frequency range; according to the return of healthy people and CHF patients A normalized energy distribution map is drawn with the unified energy distribution vector, where the abscissa represents the component level, the ordinate represents the normalized energy distribution vector value, the curve represents the average value, and the error bar represents the standard deviation;

As shown in Figure 8, at level 1, the energy of healthy people is higher than that of patients with CHF. As the level increases, the energy of healthy people gradually decreases; while patients with CHF gradually increase at levels 1 to 3 and reach the highest value at level 3. When the level is greater than 3, the energy gradually decreases as the level increases; the energy of healthy people is mainly concentrated on levels 1 to 4, which is in the high frequency range, indicating that the short-term regulation ability of the heart of healthy people is stronger; the low-level energy of CHF patients The relative decrease indicates that the short-term adjustment ability is reduced, and the main power is concentrated on the middle level. At the level 6 that reflects the heart rhythm, CHF patients are higher than healthy people, indicating that CHF patients have a higher level of regulating heart rhythm. The proportion of energy.

(8), the normalized ECG signal CHF patients normalized energy in the energy distribution of the vectors p _i and healthy human ECG signal distribution vector p _i for a significance test to get on each extremum mode function component the probability value P _i, P _i when the probability value less than 0.05 indicates a significant difference between the two, as shown in the table, at the level of 1,2,3,5,6 significant difference was found between the two, CHF patients Decrease in energy at low levels of decomposition may indicate that diseases cause a decrease in the heart's ability to regulate on small time scales; while healthy people have higher energy at low levels of decomposition, indicating that healthy people's hearts have better short-term adjustment capabilities and have better short-term adjustment capabilities to the external environment. Have a better ability to adapt to changes in the physical environment;

Table 2 Energy vector T test of healthy people and CHF patients

(9) Under the condition that P ₁ and P _{2 are} _{less than 0.05, the energy distribution vector p 1} of the first extreme value modal function component of CHF patients is smaller than the energy distribution vector p of the second extreme value modal function component ₂ .

Claims

An ECG signal analysis method based on extreme energy decomposition method is characterized in that it comprises the following steps:

(1) Obtain the ECG signal x(t) in an unknown state at a given time and a given sampling frequency;

(2) Perform denoising preprocessing of the ECG signal x(t); the specific preprocessing method is: pass the ECG signal through a 40Hz zero-phase FIR low-pass filter to remove high-frequency noise, and then pass a median filter to remove the baseline drift;

(3) Take the denoised ECG signal x(t) as the original signal, find all the local extreme points of the original signal, and then use the spline curve for all the maximum points and all the minimum points of the original signal They are connected to form the upper envelope e max and the lower envelope e min respectively , and the envelope mean signal m(t)=(e max +e min )/2 of the upper and lower envelopes is obtained;

(4) Subtract the envelope mean signal m(t) from the original signal x(t) to obtain h(t)=x(t)-m(t); then judge whether h(t) satisfies the extreme mode If the judgment condition of the function is not met, return h(t) as the original signal to step (3) until h k (t) meets the judgment condition of the extreme modal function, then write c 1 (t) = h k (t), as the first extreme value modal function component;

(5) Subtract the first extremum modal function component c 1 (t) from the original signal x(t) to obtain the margin r 1 (t) = x(t)-c 1 (t), and then judge Whether h k (t) meets the stopping criterion, if not, take r 1 (t) as the new original sequence x(t), and return to steps (3) and (4) until h k (t) meets the stopping criterion , Get the 2, 3,..., n extremum modal function components and margin r n (t), then decompose the original signal x(t) into n extremum modal function components and a margin, namely

(6) Perform frequency spectrum analysis on the extreme value modal function components c i (t), i=1, 2,..., n to obtain the center frequency of each extreme value modal function component;

(7) The n extremum modal function components obtained by decomposing the original signal x(t) represent the components of different frequency bands of the original signal, and then calculate the energy of each component

E i ＝∫|c i (t)| 2 dt,i＝1,2,...,n

Normalize each energy value to get a normalized energy distribution vector

p i ＝E i /E,i＝1,2,...,n

among them,
The first component p 1 represents the energy of the highest frequency band, which represents the proportion of the energy distribution of the signal in the highest frequency range, and the last component p n represents the proportion of the energy distribution of the signal in the lowest frequency range; according to the normalized energy distribution vector Draw a normalized energy distribution map, where the abscissa represents the component level, the ordinate represents the normalized energy distribution vector value, the curve represents the average value, and the error bar represents the standard deviation;

Probability on each extremum mode function component normalization energy distribution vector p i and a significance detector (8), the normalized unknown state ECG signal energy distribution vector p i and the standard ECG signals to obtain Value P i , judge whether the probability value P 1 on the first extreme value modal function component and the probability value P 2 on the second extreme value modal function component are less than the probability standard value, if not, return to step (1) again Get the signal;

(9) Under the condition that P 1 and P 2 are less than the probability standard value, if the energy distribution vector p 1 of the first extreme value modal function component is smaller than the energy distribution vector p of the second extreme value modal function component 2 , it is determined that the ECG signal is an abnormal electrocardiogram signal.
An ECG signal analysis method based on the extreme value energy decomposition method according to claim 1, characterized in that: the minimum amount of data required for the original signal x(t) is N=2 n+1 , where n is the decomposed extreme The number of components of the value modal function.
An ECG signal analysis method based on the extreme value energy decomposition method according to claim 1, characterized in that the determination condition of the extreme value modal function in the step (4) is: (a). In the entire data sequence, The number of extreme points is equal to or one difference from the number of zero-crossing points; (b). At any moment, the upper and lower envelopes are symmetrical to the time axis.
An electrocardiogram signal analysis method based on extreme energy decomposition method according to claim 1, characterized in that: in the step (5), the formula that h k (t) satisfies the stopping criterion is:

ε represents the screening threshold, which is between 0.2 and 0.3.
An ECG signal analysis method based on extreme energy decomposition method according to claim 1, characterized in that: the probability standard value in step (8) and step (9) is 0.05.
An ECG signal analysis method based on the extreme value energy decomposition method according to claim 1, characterized in that: the repetition frequency of the sixth extreme value modal function component waveform obtained in the step (5) is the human heart rate, namely Cardiac cycle.