CN110464337A - A kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method - Google Patents
A kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method Download PDFInfo
- Publication number
- CN110464337A CN110464337A CN201910841605.6A CN201910841605A CN110464337A CN 110464337 A CN110464337 A CN 110464337A CN 201910841605 A CN201910841605 A CN 201910841605A CN 110464337 A CN110464337 A CN 110464337A
- Authority
- CN
- China
- Prior art keywords
- signal
- extreme value
- component
- energy
- mode function
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
- A61B5/0205—Simultaneously evaluating both cardiovascular conditions and different types of body conditions, e.g. heart and respiratory condition
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
- A61B5/024—Detecting, measuring or recording pulse rate or heart rate
- A61B5/02405—Determining heart rate variability
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
- A61B5/024—Detecting, measuring or recording pulse rate or heart rate
- A61B5/0245—Detecting, measuring or recording pulse rate or heart rate by using sensing means generating electric signals, i.e. ECG signals
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
- A61B5/318—Heart-related electrical modalities, e.g. electrocardiography [ECG]
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7203—Signal processing specially adapted for physiological signals or for diagnostic purposes for noise prevention, reduction or removal
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7235—Details of waveform analysis
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7235—Details of waveform analysis
- A61B5/725—Details of waveform analysis using specific filters therefor, e.g. Kalman or adaptive filters
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Surgery (AREA)
- General Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Medical Informatics (AREA)
- Molecular Biology (AREA)
- Physics & Mathematics (AREA)
- Animal Behavior & Ethology (AREA)
- Pathology (AREA)
- Public Health (AREA)
- Veterinary Medicine (AREA)
- Cardiology (AREA)
- Physiology (AREA)
- Signal Processing (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Psychiatry (AREA)
- Pulmonology (AREA)
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
Abstract
The invention discloses a kind of heart rate variability signals analysis methods based on extreme value Energy Decomposition method, and the ECG signal including obtaining the unknown state under given time and given sample frequency obtains RRI signal x (t) after denoising；It regard RRI signal x (t) as original signal, original signal x (t) is decomposed into n extreme value mode function component and a surplus, the n extreme value mode function component that original signal x (t) is decomposed, the component of original signal different frequency range is represented, determines whether the RRI signal is abnormal cardiac rate Variability Signals according to n extreme value mode function component.The present invention analyzes RRI signal using extreme value Energy Decomposition method, and original signal is decomposed into multiple components, that is, extreme value component function, the energy of each component is calculated, obtains its Energy distribution.
Description
Technical field
The present invention relates to a kind of analysis of ECG signal more particularly to a kind of heart rate variabilities based on extreme value Energy Decomposition method
Property signal analysis method.
Background technique
Physiological signal is generated by the multiple system interactions of life entity, and the time and intensity of different systemic effects is not
Together, lead to the complexity of physiological signal with the time and spatially.Heart rate variability (HRV) refer to measurement continuous cardiac cycle it
Between time-variance number, precisely, it should be the variance for measuring the difference between the normal P-P interphase continuously occurred.But
Since P wave is not as good as R wave is obvious or P wave crest end has that time width is blunt, thus to the research of heart rate variability signals usually with P-P interphase
Equal R -- R interval signal (RRI) replaces.Studies have shown that the Non-invaive examination that HRV can be used as activity of vegetative nervous system refers to
Mark, the important in inhibiting especially in terms of the prognosis for judging certain cardiovascular diseases.
Study in the prior art heart it is long when adjust the rhythm and pace of moving things (< 1Hz)) common heart-rate variability signal (HRV signal) make
To analyze object.It is a large amount of studies have shown that correlation and nonlinear kinetics complexity when human body HRV signal has long, and
Age and disease will lead to Dynamic complexity reduction.To heart rate variability (Heart Rate Varibility, HRV) signal
Time of the research the most commonly used is RR interphase (Interbeat Intervals, RRI) signal, i.e., between continuous RRI signal R wave
Blank signal.
The most common method of energy change for studying HRV signal is power spectrumanalysis (PSD).PSD passes through Fourier transformation
The power of HRV signal is converted to the function of frequency, studies the watt level of different frequency domains, usual HRV frequency spectrum is divided into height
Frequently several frequency ranges such as (HF), low frequency (LF) and very low frequencies (VLF).LF/HF ratio has important clinical value.Heart disease can draw
The change of HRV power spectrum is played, for example heart failure and myocardial infarction cause normalization HF to increase, LF and VLF reduction.However, PSD is not
A method of based on data-driven, and it is relatively rough to the segmentation of frequency domain, cause details to lack, divides also inflexible.
It would therefore be highly desirable to solve the above problems.
Summary of the invention
Goal of the invention: little data, which can be used, the object of the present invention is to provide one kind can intuitively reflect electrocardiogram energy point
The heart rate variability signals analysis method based on extreme value Energy Decomposition method of the true rule of cloth and signal fluctuation.
Technical solution: in order to achieve the above object, the invention discloses a kind of heart rate variabilities based on extreme value Energy Decomposition method
Property signal analysis method, includes the following steps:
(1), the ECG signal of given time and the unknown state under given sample frequency are obtained, then ECG signal is carried out
Noise suppression preprocessing therefrom extracts RRI signal, obtains the RRI signal x (t) of unknown state；
(2), it regard RRI signal x (t) as original signal, finds out all Local Extremums of original signal, it then will be original
All maximum points and all minimum points of signal are linked up using spline curve and are respectively formed coenvelope line emaxAnd lower envelope
Line emin, obtain envelope mean value signal m (t)=(e of upper and lower envelopemax+emin)/2；
(3), original signal x (t) is subtracted into envelope mean value signal m (t), obtains h (t)=x (t)-m (t)；Then judge h
(t) whether meet the decision condition of extreme value mode function, if conditions are not met, h (t) is back to step (3) as original signal,
Until hk(t) decision condition for meeting extreme value mode function, then remember c1(t)=hk(t), as first extreme value mode function point
Amount；
(4), original signal x (t) is subtracted into first extreme value mode function component c1(t), surplus r is obtained1(t)=x (t)-
c1(t), then judge hk(t) whether meet stopping criterion, if conditions are not met, by r1(t) it as new original series x (t), returns
It is back to step (2) and (3)；If hk(t) when meeting stopping criterion but n < 8, return step (1) reacquires original signal；Such as
Fruit hk(t) when meeting stopping criterion and n >=8, obtain the 2nd, 3 ..., n extreme value mode function component and surplus rn(t), then will
Original signal x (t) is decomposed into n extreme value mode function component and a surplus, i.e.,
(5), to extreme value mode function component ci(t), i=1,2 ..., n carry out spectrum analysis and obtain each extreme value mode letter
The centre frequency of number component；
(6), n extreme value mode function component for decomposing original signal x (t), represents original signal different frequency range
Component, then calculate the energy of its each component
Ei=∫ | ci(t)|2Dt, i=1,2 ..., n
Each energy value is normalized, normalized Energy distribution vector is obtained
pi=Ei/ E, i=1,2 ..., n
Wherein,One-component p1The energy for indicating highest frequency range, represents signal in highest frequency range model
The ratio of energy distribution is enclosed, the last one component pnIndicate signal in the ratio of peak low band range energy distribution；According to
Normalized Energy distribution vector-drawn normalized energy distribution map, wherein abscissa indicates component level, and ordinate expression is returned
The one Energy distribution vector value changed, curve indicate that average value, error bar indicate standard deviation；
(7), second component p is chosen2To the 7th component p7, energy difference calculated different value EDV, EDV=(p2+p3+p4)-
(p5+p6+p7), when EDV≤first standard value, then determine that the RRI signal is normal heart rate variability signals, when the first standard
When value < EDV the second standard value of <, then determine that the RRI signal is doubtful abnormal cardiac rate Variability Signals, when the standard of EDV >=second
When value, then determine the RRI signal for abnormal cardiac rate Variability Signals.
Wherein, minimal data amount N=2 needed for the original signal x (t)n+1, wherein n is the extreme value mode function decomposited
The quantity of component.
Preferably, noise suppression preprocessing in the step (1) method particularly includes: ECG signal is passed through into 40Hz zero phase FIR
High-frequency noise is eliminated in low-pass filter filtering, then removes baseline drift by median filter.
Furthermore the decision condition of extreme value mode function in the step (3) are as follows: (a), in entire data sequence, extreme value
Point quantity it is equal with the quantity of zero crossing or difference one；(b), at any time, upper and lower envelope is for time shaft pair
Claim.
Further, h in the step (4)k(t) meet stopping criterion formula are as follows:
ε indicates screening thresholding, takes between 0.2~0.3.
Preferably, the first standard value in the step (7) is -0.15, and the second standard value is 0.08.
The utility model has the advantages that compared with prior art, the present invention has following remarkable advantage:
The present invention analyzes RRI using extreme value Energy Decomposition method (Extremum Energy Decomposition, EED)
Original signal is decomposed into multiple components, that is, extreme value component function, calculates the energy of each component, obtain it by signal
Energy distribution；Of the invention by signal decomposition can be difference from high frequency to low frequency according to the fluctuation pattern of raw RRI signal itself
Time multi-level signal, it is more careful to the segmentation of frequency range；The data length that polar decomposition obtains on all levels is identical, thus
It not will lead to data length reduction, to allow to analyze for short time data, that is, need less data amount that can analyze
Obtain accurate result；EED is not easy the analysis of different levels component energy affected by noise.
Detailed description of the invention
Fig. 1 is the schematic diagram of original signal in the present invention；
Fig. 2 is the schematic diagram that original signal seeks envelope in the present invention；
Fig. 3 is the schematic diagram for subtracting envelope mean value signal of original signal in the present invention；
Fig. 4 is to obtain the schematic diagram of first extreme value mode function component in the present invention；
Fig. 5 is the flow diagram of extreme value Energy Decomposition method in the present invention；
Fig. 6 is the EED decomposition diagram of RRI signal in the embodiment of the present invention 1；
Fig. 7 is normalized energy distribution schematic diagram in RRI signal in the embodiment of the present invention 1.
Specific embodiment
Technical solution of the present invention is described further with reference to the accompanying drawing.
As shown in Figure 1, Figure 2, shown in Fig. 3, Fig. 4 and Fig. 5, a kind of heart rate variability based on extreme value Energy Decomposition method of the invention
Signal analysis method includes the following steps:
(1), the ECG signal of given time and the unknown state under given sample frequency are obtained, then ECG signal is carried out
Noise suppression preprocessing therefrom extracts RRI signal, obtains the RRI signal x (t) of unknown state；The wherein specific method of noise suppression preprocessing
Are as follows: since ECG energy is concentrated mainly on 0~40Hz, ECG signal is eliminated by the filtering of 40Hz zero phase FIR low pass filter
Then high-frequency noise removes baseline drift by median filter；
(2), it regard RRI signal x (t) as original signal, minimal data amount N=2 needed for original signal x (t)n+1, wherein n
For the quantity of the extreme value mode function component decomposited；All Local Extremums for finding out original signal, then by original signal
All maximum points and all minimum points link up using spline curve and be respectively formed coenvelope line emaxWith lower envelope line
emin, obtain envelope mean value signal m (t)=(e of upper and lower envelopemax+emin)/2；
(3), original signal x (t) is subtracted into envelope mean value signal m (t), obtains h (t)=x (t)-m (t)；Then judge h
(t) whether meet the decision condition of extreme value mode function, if conditions are not met, h (t) is back to step (2) as original signal,
Until hk(t) decision condition for meeting extreme value mode function, then remember c1(t)=hk(t), as first extreme value mode function point
Amount；The wherein decision condition of extreme value mode function are as follows: (a), in entire data sequence, the quantity of extreme point and the number of zero crossing
Measure equal or difference one；(b), at any time, upper and lower envelope is for time axial symmetry；
(4), original signal x (t) is subtracted into first extreme value mode function component c1(t), surplus r is obtained1(t)=x (t)-
c1(t), then judge hk(t) whether meet stopping criterion, if conditions are not met, by r1(t) it as new original series x (t), returns
It is back to step (2) and (3)；If hk(t) when meeting stopping criterion but n < 8, return step (1) reacquires original signal；Such as
Fruit hk(t) when meeting stopping criterion and n >=8, obtain the 2nd, 3 ..., n extreme value mode function component and surplus rn(t), then will
Original signal x (t) is decomposed into n extreme value mode function component and a surplus, i.e.,
Wherein hk(t) meet stopping criterion formula are as follows:
ε indicates screening thresholding, takes between 0.2~0.3；Meet stopping criterion
Extreme value mode decomposition then meets two following conditions: the extreme value mode function component c (a) finally obtainedn(t) or surplus rn(t)
Less than preset threshold value；(b) residue signal rn(t) become monotonic signal, cannot therefrom extract extreme value mode function again
Signal；
(5), to extreme value mode function component ci(t), i=1,2 ..., n carry out spectrum analysis and obtain each extreme value mode letter
The centre frequency of number component；
(6), n extreme value mode function component for decomposing original signal x (t), represents original signal different frequency range
Component, then calculate the energy of its each component
Ei=∫ | ci(t)|2Dt, i=1,2 ..., n
Each energy value is normalized, normalized Energy distribution vector is obtained
pi=Ei/ E, i=1,2 ..., n
Wherein,One-component p1The energy for indicating highest frequency range, represents signal in highest frequency range model
The ratio of energy distribution is enclosed, the last one component pnIndicate signal in the ratio of peak low band range energy distribution；According to
Normalized Energy distribution vector-drawn normalized energy distribution map, wherein abscissa indicates component level, and ordinate expression is returned
The one Energy distribution vector value changed, curve indicate that average value, error bar indicate standard deviation；
(7), second component p is chosen2To the 7th component p7, energy difference calculated different value EDV, EDV=(p2+p3+p4)-
(p5+p6+p7), when EDV≤first standard value, then determine that the RRI signal is normal heart rate variability signals, when the first standard
When value < EDV the second standard value of <, then determine that the RRI signal is doubtful abnormal cardiac rate Variability Signals, when the standard of EDV >=second
When value, then determine the RRI signal for abnormal cardiac rate Variability Signals；First standard value therein is -0.15, and the second standard value is
0.08。
The extreme value Energy Decomposition method (Extremum Energy Decomposition, EED) that the present invention uses is
The method of concept based on extreme value mode function, extreme value mode function be meet simultaneously following two condition have single-frequency
A kind of signal, two conditions are as follows:
(a), in entire data sequence, the quantity of extreme point (including maximum and minimum) and the quantity of zero crossing must
Must it is equal or it is most difference one；
(b) at any time, the lower envelope line that the coenvelope line and local minizing point that Local modulus maxima is formed are formed
Mean value be zero, that is to say, that part up and down envelope for time shaft Local Symmetric；
Both the above condition, condition (a) are similar to requirement of the Gaussian normal stationary process for traditional narrow, condition (b)
It ensure that the instantaneous frequency being calculated by extreme value mode function has physical significance.
Extreme value mode function of the invention decomposes the standard selection terminated and wants moderate, and condition is too stringent, will lead to last several
A extreme value mode function component loses meaning；Condition is too loose, will lead to useful component and loses；In practical applications, may be used
To set the extreme value mode function component number of plies for needing to decompose according to demand, terminates and calculate when meeting Decomposition order.
Embodiment 1
EED analysis method is used to analyze the Energy distribution of Healthy People and CHF patient ECG under different levels.
A kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method of Healthy People, includes the following steps:
(1), the ECG signal of healthy population is obtained from the RR interval data library nsr2db of physionet；Wherein data include
(at the age 28~76, average 61) then to carry out noise suppression preprocessing to ECG signal, therefrom extraction RRI signal, obtains 54 Healthy Peoples
To RRI signal x (t)；Wherein noise suppression preprocessing method particularly includes: since ECG energy is concentrated mainly on 0~40Hz, ECG is believed
Number by 40Hz zero phase FIR low pass filter filtering eliminate high-frequency noise, then by median filter remove baseline drift；
(2), it regard RRI signal x (t) as original signal, minimal data amount N=2 needed for original signal x (t)n+1=29,
Middle n is the quantity of the extreme value mode function component decomposited, n=8；All Local Extremums of original signal are found out, then will
All maximum points and all minimum points of original signal are linked up using spline curve and are respectively formed coenvelope line emaxWith under
Envelope emin, obtain envelope mean value signal m (t)=(e of upper and lower envelopemax+emin)/2；
(3), original signal x (t) is subtracted into envelope mean value signal m (t), obtains h (t)=x (t)-m (t)；Then judge h
(t) whether meet the decision condition of extreme value mode function, if conditions are not met, h (t) is back to step (2) as original signal,
Until hk(t) decision condition for meeting extreme value mode function, then remember c1(t)=hk(t), as first extreme value mode function point
Amount；The wherein decision condition of extreme value mode function are as follows: (a), in entire data sequence, the quantity of extreme point and the number of zero crossing
Measure equal or difference one；(b), at any time, upper and lower envelope is for time axial symmetry；
(5), original signal x (t) is subtracted into first extreme value mode function component c1(t), surplus r is obtained1(t)=x (t)-
c1(t), by r1(t) as new original series x (t), step (2) and (3) are back to, obtain the 2nd, 3 ..., 8 extreme value mode
Function component and surplus r8(t), original signal x (t) is then decomposed into 8 extreme value mode function components and a surplus, i.e.,
Obtain the extreme value mode decomposition schematic diagram of Healthy People as shown in FIG. 6, it can be seen that component 1 has highest frequency
Rate, signal fluctuate in shortest time scale, and as component serial number increases, frequency is gradually decreased；
(6), to extreme value mode function component ci(t), i=1,2 ..., 8, it carries out spectrum analysis and obtains each extreme value mode letter
The centre frequency of number component, obtains frequency-domain analysis result figure；
(7), 8 extreme value mode function components for decomposing original signal x (t), represent original signal different frequency range
Component, then calculate the energy of its each component
Ei=∫ | ci(t)|2Dt, i=1,2 ..., 8
Each energy value is normalized, normalized Energy distribution vector is obtained
pi=Ei/ E, i=1,2 ..., 8
Wherein,One-component p1The energy for indicating highest frequency range, represents signal in highest frequency range model
The ratio of energy distribution is enclosed, the last one component pnIndicate signal in the ratio of peak low band range energy distribution.
(8), second component p is chosen2To the 7th component p7, energy difference calculated different value EDV, EDV=(p2+p3+p4)-
(p5+p6+p7)=- 0.2092 ± 0.2940.
A kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method of CHF patient, includes the following steps:
(1), ECG signal is obtained from the RR interval data library chf2db database of physionet, wherein chfdb database
Comprising 29 congestive heart failures (Congestive Heart Failure, CHF) patient, (55), then the age 34~79 is averaged
Noise suppression preprocessing is carried out to ECG signal, RRI signal is therefrom extracted, obtains RRI signal x (t)；Wherein noise suppression preprocessing is specific
Method are as follows: since ECG energy is concentrated mainly on 0~40Hz, ECG signal is filtered by 40Hz zero phase FIR low pass filter
High-frequency noise is eliminated, then removes baseline drift by median filter；
(2), it regard RRI signal x (t) as original signal, minimal data amount N=2 needed for original signal x (t)n+1=29,
Middle n is the quantity of the extreme value mode function component decomposited, n=8；All Local Extremums of original signal are found out, then will
All maximum points and all minimum points of original signal are linked up using spline curve and are respectively formed coenvelope line emaxWith under
Envelope emin, obtain envelope mean value signal m (t)=(e of upper and lower envelopemax+emin)/2；
(3), original signal x (t) is subtracted into envelope mean value signal m (t), obtains h (t)=x (t)-m (t)；Then judge h
(t) whether meet the decision condition of extreme value mode function, if conditions are not met, h (t) is back to step (2) as original signal,
Until hk(t) decision condition for meeting extreme value mode function, then remember c1(t)=hk(t), as first extreme value mode function point
Amount；The wherein decision condition of extreme value mode function are as follows: (a), in entire data sequence, the quantity of extreme point and the number of zero crossing
Measure equal or difference one；(b), at any time, upper and lower envelope is for time axial symmetry；
(4), original signal x (t) is subtracted into first extreme value mode function component c1(t), surplus r is obtained1(t)=x (t)-
c1(t), by r1(t) as new original series x (t), step (2) and (3) are back to, obtain the 2nd, 3 ..., 8 extreme value mode
Function component and surplus r8(t), original signal x (t) is then decomposed into 8 extreme value mode function components and a surplus, i.e.,
(5), to extreme value mode function component ci(t), i=1,2 ..., 8, it carries out spectrum analysis and obtains each extreme value mode letter
The centre frequency of number component；
(6), 8 extreme value mode function components for decomposing original signal x (t), represent original signal different frequency range
Component, then calculate the energy of its each component
Ei=∫ | ci(t)|2Dt, i=1,2 ..., 8
Each energy value is normalized, normalized Energy distribution vector is obtained
pi=Ei/ E, i=1,2 ..., 8
Wherein,One-component p1The energy for indicating highest frequency range, represents signal in highest frequency range model
The ratio of energy distribution is enclosed, the last one component pnIndicate signal in the ratio of peak low band range energy distribution；According to
The normalized Energy distribution vector-drawn normalized energy distribution map of Healthy People and CHF patient, wherein abscissa indicates component
Level, ordinate indicate normalized Energy distribution vector value, and curve indicates that average value, error bar indicate standard deviation；
(7), second component p is chosen2To the 7th component p7, energy difference calculated different value EDV, EDV=(p2+p3+p4)-
(p5+p6+p7)=0.2642 ± 0.4070.
The present invention in order to further analyze Healthy People and CHF patient HRV signal different levels component mean center frequency
Rate obtains the result of table 1.
The mean center frequency of 1 Healthy People of table and CHF patient HRV signal difference component level
It can be concluded that as the centre frequency that component level increases HRV gradually decreases.
Typical way of traditional power spectral density (Power Spectral Density, PSD) method to Dividing in frequency domain
Are as follows: HF (0.15~0.4Hz), LF (0.04~0.15Hz), VLF (0.0033~0.04Hz).EED method of the present invention is to component layers
Secondary 1, frequency is higher than HF；Frequency range of the level 2 in HF；Range of the level 3,4 in LF；Frequency range of the level 5~7 in VLF；
Level 8 is lower than VLF.In addition it can be seen that the frequency of CHF patient is slightly above Healthy People when identical level, heart disease pair has been reacted
HRV fluctuates the influence of the rhythm and pace of moving things, and CHF patient HRV is fluctuated faster under same level.
The extreme value mode function component C that the present invention obtains the HRV signal decomposition of two groupsi(t), it calculates and is normalized
Energy distribution vector, and make EED curve graph, as shown in Figure 7.Fig. 7 is Healthy People and CHF patient RR interphase signal EED analysis knot
Fruit schematic diagram, wherein data length is 10000 points, and curve indicates mean value, and error bar indicates standard deviation.Symbol * above curve
Indicate that two groups of crowd's energy T examine p < 0.01.On hierarchy selection, level 1 and the level higher than 7, including surplus are eliminated.Layer
Secondary 1 is easily affected by noise, and leads to the biggish fluctuation of energy, causes result standard deviation excessive, and its frequency may be up to
Several kHz or more, because without specific physiological significance；Higher than the long Shi Jielv that 7 level has reacted signal, it is highly susceptible to
The influence of external environment, and its frequency is very low, and physiological significance is unknown.In Fig. 7 on component low level (level 2,3), CHF suffers from
The normalized energy value of person is apparently higher than Healthy People, and on component high-level (level > 5), opposite variation, Healthy People energy occurs
Amount is higher than CHF patient energy；Healthy People energy is relatively more steady in level 2~5, and when level is greater than 5, energy slowly rises, and CHF
Patient declines rapidly in level 2~4, tends towards stability later；On 4 levels (2,3,6,7), the two energy has conspicuousness difference
(p<0.01).In order to compare, invention increases the EED of alternate data (Healthy Surrogate, CHF Surrogate)
As a result, as shown in fig. 7, alternate data is generated by the method for being randomized initial data, the energy of alternate data divides for analysis
Cloth increases monotonic decreasing with scale, and compare CHF patient, and energy is higher in small time scale, the energy in long time scale
It measures lower.
In order to assess Energy distribution difference of the EED curve on low level and high-level decomposition, the energy of each crowd is calculated
Difference value EDV, high EDV value indicate the higher component low level Energy distribution of RRI signal and the high-level energy of lower component
Distribution.The EDV value of Healthy People, CHF patient and their alternate datas is calculated, as shown in table 2.
The EDV value of 2 Healthy People of table and CHF patient
* Healthy People and CHF patient T inspection result p < 0.0001 are represented.
* represents alternate data and its initial data T inspection result p < 0.0001.
From table 2 it can be seen that Healthy People and CHF patient EDV value have marked difference, CHF patient EDV value compares Healthy People
It is much higher.The EDV value of Healthy People illustrates there is higher energy component is high-level, implies that high-rise component of degree n n has more less than 0
High adjusting strength.The EDV value of two crowds and their alternate data EDV value have difference, the EDV of human body HRV
It is significantly less than random sequence.
Heart rate variability (HRV) refers to the time-variance number measured between continuous cardiac cycle, precisely, it should be to survey
Measure the variance of the difference between the normal P-P interphase continuously occurred.But since P wave is obvious not as good as R wave or there is time width at P wave crest end
It is blunt, so the research to heart rate variability signals is usually replaced with the R -- R interval signal (RRI) equal with P-P interphase.Research
Show that HRV can be used as the Non-invaive examination index of activity of vegetative nervous system, especially in the prognosis for judging certain cardiovascular diseases
Aspect important in inhibiting.
Clinical with practical application, the present invention proposes that the letter of heartbeat in short-term of a continuous RR interphase of 512 (n=8) can be used
Number it is effective for above-mentioned EED analysis, and has certain data surplus；Show that EED decomposition reaches 7 by the studies above
Component level (n=7) is existing well as a result, data volume needed for then at least can use N=2n+1=27+1=256.
Claims (6)
1. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method, which comprises the steps of:
(1), the ECG signal for obtaining given time and the unknown state under given sample frequency, then denoises ECG signal
Pretreatment therefrom extracts RRI signal, obtains the RRI signal x (t) of unknown state；
(2), it regard RRI signal x (t) as original signal, all Local Extremums of original signal is found out, then by original signal
All maximum points and all minimum points link up using spline curve and be respectively formed coenvelope line emaxWith lower envelope line
emin, obtain envelope mean value signal m (t)=(e of upper and lower envelopemax+emin)/2；
(3), original signal x (t) is subtracted into envelope mean value signal m (t), obtains h (t)=x (t)-m (t)；Then judging h (t) is
The no decision condition for meeting extreme value mode function, if conditions are not met, h (t) is back to step (3) as original signal, until
hk(t) decision condition for meeting extreme value mode function, then remember c1(t)=hk(t), as first extreme value mode function component；
(4), original signal x (t) is subtracted into first extreme value mode function component c1(t), surplus r is obtained1(t)=x (t)-c1
(t), then judge hk(t) whether meet stopping criterion, if conditions are not met, by r1(t) it as new original series x (t), returns
To step (2) and (3)；If hk(t) when meeting stopping criterion but n < 8, return step (1) reacquires original signal；If
hk(t) when meeting stopping criterion and n >=8, obtain the 2nd, 3 ..., n extreme value mode function component and surplus rnIt (t), then will be former
Beginning signal x (t) is decomposed into n extreme value mode function component and a surplus, i.e.,
(5), to extreme value mode function component ci(t), i=1,2 ..., n carry out spectrum analysis and obtain each extreme value mode function component
Centre frequency；
(6), n extreme value mode function component for decomposing original signal x (t) represents point of original signal different frequency range
Amount, then calculates the energy of its each component
Ei=∫ | ci(t)|2Dt, i=1,2 ..., n
Each energy value is normalized, normalized Energy distribution vector is obtained
pi=Ei/ E, i=1,2 ..., n
Wherein,One-component p1The energy for indicating highest frequency range, represents signal in highest band limits
The ratio of Energy distribution, the last one component pnIndicate signal in the ratio of peak low band range energy distribution；According to normalizing
The Energy distribution vector-drawn normalized energy distribution map of change, wherein abscissa indicates that component level, ordinate indicate normalization
Energy distribution vector value, curve indicate average value, error bar indicate standard deviation；
(7), second component p is chosen2To the 7th component p7, energy difference calculated different value EDV, EDV=(p2+p3+p4)-(p5+p6+
p7), when EDV≤first standard value, then determine that the RRI signal is normal heart rate variability signals, as the first standard value < EDV
When the second standard value of <, then determine that the RRI signal is doubtful abnormal cardiac rate Variability Signals, when EDV >=second standard value, then
Determine the RRI signal for abnormal cardiac rate Variability Signals.
2. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method according to claim 1, feature exist
In: minimal data amount N=2 needed for the original signal x (t)n+1, wherein n is the number of the extreme value mode function component decomposited
Amount.
3. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method according to claim 1, feature exist
In: noise suppression preprocessing in the step (1) method particularly includes: filter ECG signal by 40Hz zero phase FIR low pass filter
Wave eliminates high-frequency noise, then removes baseline drift by median filter.
4. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method according to claim 1, feature exist
In the decision condition of extreme value mode function in the step (3) are as follows: (a), in entire data sequence, the quantity of extreme point with
The quantity of zero crossing it is equal or difference one；(b), at any time, upper and lower envelope is for time axial symmetry.
5. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method according to claim 1, feature exist
In: h in the step (4)k(t) meet stopping criterion formula are as follows:
ε indicates screening thresholding, takes between 0.2~0.3.
6. a kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method according to claim 1, feature exist
In: the first standard value in the step (7) is -0.15, and the second standard value is 0.08.
Priority Applications (3)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910841605.6A CN110464337B (en) | 2019-09-06 | 2019-09-06 | Heart rate variability signal analysis method based on extreme value energy decomposition method |
| PCT/CN2019/121416 WO2021042591A1 (en) | 2019-09-06 | 2019-11-28 | Heart rate variability signal analysis method based on extremum energy decomposition method |
| AU2019101756A AU2019101756A4 (en) | 2019-09-06 | 2019-11-28 | Method for analyzing heart rate variability signal based on extremum energy decomposition method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910841605.6A CN110464337B (en) | 2019-09-06 | 2019-09-06 | Heart rate variability signal analysis method based on extreme value energy decomposition method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110464337A true CN110464337A (en) | 2019-11-19 |
| CN110464337B CN110464337B (en) | 2020-09-01 |
Family
ID=68515060
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910841605.6A Active CN110464337B (en) | 2019-09-06 | 2019-09-06 | Heart rate variability signal analysis method based on extreme value energy decomposition method |
Country Status (3)
| Country | Link |
|---|---|
| CN (1) | CN110464337B (en) |
| AU (1) | AU2019101756A4 (en) |
| WO (1) | WO2021042591A1 (en) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2021042591A1 (en) * | 2019-09-06 | 2021-03-11 | 江苏华康信息技术有限公司 | Heart rate variability signal analysis method based on extremum energy decomposition method |
| CN116548928A (en) * | 2023-07-11 | 2023-08-08 | 西安浩阳志德医疗科技有限公司 | Nursing service system based on internet |
| CN117503153A (en) * | 2024-01-05 | 2024-02-06 | 北华大学 | Patient postoperative rehabilitation evaluation method based on artificial intelligence |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116738153B (en) * | 2023-08-15 | 2023-10-17 | 中国科学院东北地理与农业生态研究所 | Organic fertilizer utilization effect evaluation method based on spectral analysis |
| CN117390380B (en) * | 2023-12-12 | 2024-02-13 | 泰安金冠宏食品科技有限公司 | Data analysis method in oil-residue separation system |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080269628A1 (en) * | 2007-04-25 | 2008-10-30 | Siemens Medical Solutions Usa, Inc. | Denoising and Artifact Rejection for Cardiac Signal in a Sensis System |
| CN101496716A (en) * | 2009-02-26 | 2009-08-05 | 周洪建 | Measurement method for detecting sleep apnoea with ECG signal |
| CN101642369A (en) * | 2008-08-04 | 2010-02-10 | 南京大学 | Autonomic nervous function biological feedback method and system |
| CN108814579A (en) * | 2018-04-16 | 2018-11-16 | 西安交通大学 | A method of electrocardio, breathing combined calculation heart rate variability based on EMD decomposition |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110464337B (en) * | 2019-09-06 | 2020-09-01 | 江苏华康信息技术有限公司 | Heart rate variability signal analysis method based on extreme value energy decomposition method |
-
2019
- 2019-09-06 CN CN201910841605.6A patent/CN110464337B/en active Active
- 2019-11-28 AU AU2019101756A patent/AU2019101756A4/en active Active
- 2019-11-28 WO PCT/CN2019/121416 patent/WO2021042591A1/en active Application Filing
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080269628A1 (en) * | 2007-04-25 | 2008-10-30 | Siemens Medical Solutions Usa, Inc. | Denoising and Artifact Rejection for Cardiac Signal in a Sensis System |
| CN101642369A (en) * | 2008-08-04 | 2010-02-10 | 南京大学 | Autonomic nervous function biological feedback method and system |
| CN101496716A (en) * | 2009-02-26 | 2009-08-05 | 周洪建 | Measurement method for detecting sleep apnoea with ECG signal |
| CN108814579A (en) * | 2018-04-16 | 2018-11-16 | 西安交通大学 | A method of electrocardio, breathing combined calculation heart rate variability based on EMD decomposition |
Non-Patent Citations (2)
| Title |
|---|
| 庞春颖等: "一种基于改进经验模态分解的癫痫脑电识别新方法", 《中国生物医学工程学报》 * |
| 董红生等: "基于Hilbert谱的心率变异信号时频分析方法", 《仪器仪表学报》 * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2021042591A1 (en) * | 2019-09-06 | 2021-03-11 | 江苏华康信息技术有限公司 | Heart rate variability signal analysis method based on extremum energy decomposition method |
| CN116548928A (en) * | 2023-07-11 | 2023-08-08 | 西安浩阳志德医疗科技有限公司 | Nursing service system based on internet |
| CN116548928B (en) * | 2023-07-11 | 2023-09-08 | 西安浩阳志德医疗科技有限公司 | Nursing service system based on internet |
| CN117503153A (en) * | 2024-01-05 | 2024-02-06 | 北华大学 | Patient postoperative rehabilitation evaluation method based on artificial intelligence |
| CN117503153B (en) * | 2024-01-05 | 2024-03-15 | 北华大学 | Patient postoperative rehabilitation evaluation method based on artificial intelligence |
Also Published As
| Publication number | Publication date |
|---|---|
| AU2019101756A4 (en) | 2020-11-12 |
| WO2021042591A1 (en) | 2021-03-11 |
| CN110464337B (en) | 2020-09-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110464337A (en) | A kind of heart rate variability signals analysis method based on extreme value Energy Decomposition method | |
| AU2019101755A4 (en) | Method for quantitatively analyzing electrocardiogram signal based on extremum energy decomposition method | |
| Acharya et al. | Comprehensive analysis of cardiac health using heart rate signals | |
| Li et al. | On an automatic delineator for arterial blood pressure waveforms | |
| CN105286815B (en) | A kind of pulse wave signal feature point detecting method based on waveform time domain feature | |
| Podziemski et al. | Fetal heart rate discovery: algorithm for detection of fetal heart rate from noisy, noninvasive fetal ECG recordings | |
| CN108471961A (en) | A kind of computational methods of physiological parameter and corresponding Medical Devices | |
| Singh et al. | Effects of RR segment duration on HRV spectrum estimation | |
| CN109938719A (en) | A kind of Driver Fatigue Detection based on physiological parameter | |
| CN110558959B (en) | HRV signal analysis method for meditation training based on extreme value energy decomposition method | |
| CN112656393A (en) | Heart rate variability detection method and system | |
| CA2888513A1 (en) | Device and method for detecting and reporting of a stress condition of a person | |
| Faust et al. | Heart rate variability analysis for different age and gender | |
| CN110558974B (en) | Electrocardiogram signal analysis method based on extreme value energy decomposition method | |
| Garcia-Gonzalez et al. | Errors in the estimation of approximate entropy and other recurrence-plot-derived indices due to the finite resolution of RR time series | |
| Nenova et al. | An automated algorithm for fast pulse wave detection | |
| Maciejewski et al. | ECG parameter extraction and classification in noisy signals | |
| Chen et al. | A fast ECG diagnosis using frequency-based compressive neural network | |
| Daqrouq et al. | Arrhythmia detection using wavelet transform | |
| Kovalchuk et al. | Multifractal analysis of cardiac series and predictors of sudden cardiac death | |
| Fahruzi et al. | Screening of non-overlapping apnea and non-apnea from single lead ECG-apnea recordings using time-frequency approach | |
| Daoud et al. | HRV spectral estimation based on constrained Gaussian modeling in the nonstationary case | |
| Schack et al. | Adaptive methods of trend detection and their application in analysing biosignals | |
| Patil et al. | Discrimination between atrial fibrillation (AF) & normal sinus rhythm (NSR) using linear parameters | |
| Castiglioni et al. | How to check steady-state condition from cardiovascular time series |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| CB03 | Change of inventor or designer information |
Inventor after: Ning Xinbao Inventor after: Zeng Peng Inventor after: Zhou Zuojian Inventor after: Jiang Xiaodong Inventor after: Liu Guo Inventor after: Wang Binbin Inventor after: Wang Hua Inventor before: Zhou Zuojian Inventor before: Ning Xinbao Inventor before: Wang Binbin Inventor before: Jiang Xiaodong Inventor before: Wang Hua |
|
| CB03 | Change of inventor or designer information | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |