WO2019140998A1 - 一种准确提取qrs内异常电位的方法 - Google Patents

一种准确提取qrs内异常电位的方法 Download PDF

Info

Publication number
WO2019140998A1
WO2019140998A1 PCT/CN2018/116214 CN2018116214W WO2019140998A1 WO 2019140998 A1 WO2019140998 A1 WO 2019140998A1 CN 2018116214 W CN2018116214 W CN 2018116214W WO 2019140998 A1 WO2019140998 A1 WO 2019140998A1
Authority
WO
WIPO (PCT)
Prior art keywords
qrs
signal
aiqp
ecg signal
msd
Prior art date
Application number
PCT/CN2018/116214
Other languages
English (en)
French (fr)
Inventor
闫相国
吴宁
郑崇勋
王刚
Original Assignee
西安交通大学
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by 西安交通大学 filed Critical 西安交通大学
Priority to US16/609,174 priority Critical patent/US10912479B2/en
Publication of WO2019140998A1 publication Critical patent/WO2019140998A1/zh

Links

Images

Classifications

    • 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]
    • A61B5/346—Analysis of electrocardiograms
    • A61B5/349—Detecting specific parameters of the electrocardiograph cycle
    • A61B5/366—Detecting abnormal QRS complex, e.g. widening
    • 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]
    • A61B5/346—Analysis of electrocardiograms
    • A61B5/347—Detecting the frequency distribution of 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
    • A61B5/318—Heart-related electrical modalities, e.g. electrocardiography [ECG]
    • A61B5/346—Analysis of electrocardiograms
    • A61B5/349—Detecting specific parameters of the electrocardiograph cycle
    • A61B5/35—Detecting specific parameters of the electrocardiograph cycle by template matching
    • 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]
    • A61B5/346—Analysis of electrocardiograms
    • A61B5/349—Detecting specific parameters of the electrocardiograph cycle
    • A61B5/352—Detecting R peaks, e.g. for synchronising diagnostic apparatus; Estimating R-R interval
    • 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]
    • A61B5/346—Analysis of electrocardiograms
    • A61B5/349—Detecting specific parameters of the electrocardiograph cycle
    • A61B5/364—Detecting abnormal ECG interval, e.g. extrasystoles, ectopic heartbeats
    • 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/7221—Determining signal validity, reliability or quality
    • 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/7253—Details of waveform analysis characterised by using transforms

Definitions

  • the invention belongs to the field of medical signal processing, and particularly relates to a method for accurately extracting an abnormal potential in a QRS from an electrocardiographic signal.
  • SCD sudden cardiac deaths
  • Sudden cardiac death has a clear trend of rejuvenation. In some patients, sudden deaths often have no obvious aura symptoms, and the body looks healthy. Therefore, early warning of sudden cardiac death is particularly important, but there is currently no effective means of detection.
  • the non-invasive detection of abnormalities in the ventricular depolarization process mainly includes: (1) ventricular late potential (VLP) based on signal average electrocardiogram (SAECG), and (2) abnormal potentials within QRS (AIQPs), ( 3) Fragmentation of QRS waves (fQRS).
  • VLP ventricular late potential
  • SAECG signal average electrocardiogram
  • AIQPs abnormal potentials within QRS
  • fQRS Fragmentation of QRS waves
  • Ventricular late potential examination refers to the high-frequency low-amplitude fragmentation potential of the body surface information superimposed on the end of the QRS complex and extending into the ST segment. It reflects the delayed electrical activity of the ischemic region and is due to the ischemic area. Caused by slow and irregular re-entry activities in the myocardium. Late potential examination is of great value in the prediction of sudden death after acute myocardial infarction.
  • the most widely used clinically is the time domain VLP detection method. Although this analysis has a high negative predictive value, the positive predictive effect is not good.
  • the main reason for the low positive predictive rate of the VLP detection method is that only the QRS complex terminal is detected and extended to the high frequency low amplitude fragmentation potential in the ST segment.
  • a large number of animal myocardial infarction models and human body mapping basic research results clearly indicate that these high-frequency low-amplitude fragmentation potentials not only exist in the terminal region of the QRS complex, but may also be hidden in the QRS complex. In some myocardial infarction sites, the high-frequency low-amplitude fragmentation potential may only exist in the QRS complex and is not reflected in the end of the QRS complex. If these high-frequency low-amplitudes in the QRS complex are accurately extracted Fragmentation potential will significantly improve the reliability of early warning of sudden cardiac death.
  • AIQPs QRS internal abnormal potentials
  • AIQPs parameters can significantly improve the diagnostic accuracy for patients with high risk of ventricular arrhythmia.
  • the extraction accuracy of AIQPs extraction technology still cannot meet the requirements, and the robustness of extraction technology needs to be further improved.
  • the original ECG signal sampled during the ECG measurement is represented by x(i), which includes the power frequency interference p(i), the baseline drift b(i), the measurement noise n(i), and the abnormal potential in the QRS to be extracted.
  • AIQP(i) of (AIQPs) an ideal ECG signal not containing the above components is represented by x p (i);
  • x(i) x p (i)+AIQP(i)+p(i)+b(i)+n(i), (1)
  • AIQP(i) is the abnormal potential in the QRS to be extracted in the effective interval of AIQPs, and the other parts take the value of 0.
  • the signal after x(i) removal of power frequency interference and baseline drift is represented by x 2 (i).
  • the signal y(i) containing the abnormal potential and noise interference in the QRS to be extracted can be obtained:
  • the object of the present invention is to provide a method for accurately extracting an abnormal potential in a QRS, using a nonlinear transform prediction technique, combined with a spline method, to obtain an abnormal potential and other interference components not included in the QRS.
  • the ideal ECG signal x p (i) and finally extracts the abnormal potential within the QRS.
  • a method for accurately extracting an abnormal potential in a QRS comprising the following steps:
  • Step 1 pre-processing the original ECG signal x 1 (i); to obtain a pre-processed ECG signal x 2 (i); when the original ECG signal is a measured single-heart beat ECG signal, The low-pass filter and the power frequency trap are used to eliminate the influence of baseline drift and power frequency interference on the subsequent process; when the original ECG signal is the measured ECG signal containing multiple heart beats, the signal is utilized The average technique is processed to eliminate the effects of baseline drift, power frequency interference, and measurement noise on subsequent processes.
  • Step 2 Perform feature point detection on the pre-processed ECG signal x 2 (i), determine the feature point position and QRS range, and obtain the estimated ideal ECG signal by nonlinear transformation; first, perform feature points on the pre-processed ECG signal Detecting, determining the position of the feature point and the range of the QRS; secondly, filtering the pre-processed ECG signal obtained by the step 1 using two low-pass filters of different filter frequencies; and then, the two filters obtained The result is subtracted, and the difference signal is obtained, and the first zero-crossing position before and after each feature point position of the difference signal is searched; and then, the first zero-crossing position before and after each feature point position is included.
  • the time range is replaced by the low-pass filter filtering result of the higher filter frequency of the above two different frequencies, and the other part is replaced by the low-pass filter filtering result of the lower frequency of the above two different frequencies to obtain a composite signal; finally, Performing low-pass filtering on the obtained composite signal to obtain an estimated ideal ECG signal;
  • Step 3 According to the pre-processed ECG signal, the feature point position and the estimated ideal ECG signal, the spline method is used to obtain an accurate estimation of the ideal ECG signal.
  • the pre-processed ECG signal obtained by the step 1 and the estimated ideal ECG signal obtained in the second step are subtracted to obtain an error signal, and the error signal is searched for the zero-crossing position; then, the error signal obtained in the search is over.
  • the sample strip weight is 1 and the others are 0.
  • the cubic smooth spline is used to obtain accurate estimation. Ideal ECG signal;
  • Step 4 Subtracting the pre-processed ECG signal obtained by the step 1 and accurately estimating the ideal ECG signal obtained by the third step, and the subtraction result is filtered by a band pass filter to obtain a filtering result; according to the filtering The result and the QRS range obtained in step 2, the abnormal potential in the QRS is obtained by moving the standard deviation analysis technique;
  • Step 5 Perform credibility evaluation on the abnormal potential in the obtained QRS. Using the standard deviation analysis method, the credibility of the abnormal potential in the QRS obtained in the fourth step is evaluated, and it is confirmed whether the abnormal potential in the QRS obtained in the fourth step is credible, and the evaluation result is output.
  • the second step is specifically as follows:
  • the ECG feature points should at least include a QRS starting point, a QRS ending point, and a Q, R, S waveform peak point;
  • f h is the low pass filter filter frequency, 100Hz ⁇ f h ⁇ 200Hz;
  • x h (i) is the filtering result of the higher frequency low pass filter pair x 2 (i)
  • x l (i) is the filtering result of the lower frequency low pass filter pair x 2 (i).
  • the signal x s (i) is synthesized by the formula (6),
  • the estimated ideal ECG signal x 3 (i) is obtained, and f 3 is the low-pass filter filtering frequency, 100 Hz ⁇ f 3 ⁇ 200 Hz.
  • the third step is specifically as follows:
  • x 3 (i) is the estimated ideal ECG signal
  • the smoothed spline is used to obtain an accurately estimated ideal electrocardiographic signal x 4 (i).
  • step four described is specifically as follows:
  • x 4 (i) is an ideal ECG signal that is accurately estimated
  • the reference MSD value ref msd is calculated, and the interval from the first 100 ms to the QRS b of the QRS start position QRS b is defined as a reference interval, and the mean ref_msd mean and the standard deviation ref_msd std of the msd(i) in the reference interval are first calculated respectively.
  • ref msd using equation (11):
  • Ref msd ref_msd mean + ⁇ *ref_msd std , (11)
  • ⁇ is generally selected to be greater than 2;
  • AIQP e 0 and the search is stopped to determine whether AIQP e is searched. If AIQP e is equal to 0, the exit flag is returned and the failure flag is returned, otherwise continue.
  • AIQP b is the starting position of the searched AIQPs
  • AIQP e is the end position of the searched AIQPs.
  • step five described is specifically as follows:
  • ref std is the standard deviation of the reference interval y(i)
  • QRS std is the interval from the QRS starting position QRS b to the QRS ending position QRS e (i) Standard deviation.
  • the failure flag is returned, otherwise the success flag is returned and the extracted abnormal potential AIQP(i) in the QRS is also returned.
  • the invention provides a method for accurately extracting an abnormal potential in a QRS by using an ideal electrocardiographic signal secondary estimation technique.
  • nonlinear transformation technology can be used to effectively track the trend of non-ECG feature point regions, and it can effectively eliminate the possible impact of ECG feature points on extracting ideal ECG signals.
  • the spline method is used to further estimate the ideal ECG signal, and the ideal ECG signal can be accurately estimated.
  • the method of the invention requires fewer parameters to be selected, and the result is more reliable.
  • the present invention also evaluates the reliability of the extracted abnormal potential in the QRS, thereby ensuring the reliability of the results of application development using the method of the present invention.
  • a prominent feature of the method of the present invention is that the single-pulse electrocardiographic signal can be extracted in the QRS, which can greatly expand the application range of the AIQPs analysis technology.
  • the scenarios and scopes applicable to the method of the present invention include: 1) integrating the method of the present invention in an electrocardiograph to evaluate the risk of sudden cardiac death in patients undergoing conventional electrocardiographic measurements; 2) integrating the present in a conventional multi-parameter monitor Inventive method for real-time dynamic tracking and monitoring of patients with myocardial infarction; 3) development of portable or wearable devices that are easy to use, to achieve early warning of sudden cardiac death in a family environment; 4) in mobile devices (such as mobile phones) Integrating the method of the invention can provide a convenient and efficient warning method for sudden death risk of the device user.
  • Figure 1 is a flow chart of the present invention.
  • FIG. 3 is a flow chart of preprocessing according to an embodiment of the present invention.
  • FIG. 4 is a schematic diagram of a baseline drift cancellation and power frequency interference removal process according to an embodiment of the present invention.
  • FIG. 5 is a flow chart of estimating an ideal ECG signal according to an embodiment of the present invention.
  • FIG. 6 is a schematic diagram of a process for estimating an ideal ECG signal according to an embodiment of the present invention.
  • FIG. 7 is a flow chart of accurately estimating an ideal ECG signal in accordance with an embodiment of the present invention.
  • FIG. 8 is a schematic diagram of a process for accurately estimating an ideal ECG signal using a cubic smooth spline according to an embodiment of the present invention.
  • Fig. 9 is a flow chart showing the evaluation of the abnormal potential extraction and the reliability of the result in the QRS according to the embodiment of the present invention.
  • FIG. 10 is a schematic diagram showing the process of extracting abnormal potentials and evaluating credibility of results in the QRS according to an embodiment of the present invention.
  • Fig. 11 is a graph showing the results of the abnormal potential extraction in the QRS by the method of the present invention for the single beat ECG and the multiple beat ECG superimposition average.
  • Figure 12 is a graph showing the results of an abnormal potential extraction within a QRS of a single beat ECG of a patient with myocardial infarction by the method of the present invention.
  • Figure 13 is a graph showing the results of an abnormal potential extraction within QRS of a single beat ECG of a healthy person by the method of the present invention.
  • FIG. 1 there is shown a flow chart of the present invention, a method for accurately extracting an abnormal potential in a QRS, comprising the following steps:
  • Step 1 pre-processing the original ECG signal x 1 (i); to obtain a pre-processed ECG signal x 2 (i); when the original ECG signal is a measured single-heart beat ECG signal, The low-pass filter and the power frequency trap are used to eliminate the influence of baseline drift and power frequency interference on the subsequent process; when the original ECG signal is the measured ECG signal containing multiple heart beats, the signal is utilized The average technique is processed to eliminate the effects of baseline drift, power frequency interference, and measurement noise on subsequent processes.
  • Step 2 Perform feature point detection on the pre-processed ECG signal x 2 (i), determine the feature point position and QRS range, and obtain the estimated ideal ECG signal by nonlinear transformation; first, perform feature points on the pre-processed ECG signal Detecting, determining the position of the feature point and the range of the QRS; secondly, filtering the pre-processed ECG signal obtained by the step 1 using two low-pass filters of different filter frequencies; and then, the two filters obtained The result is subtracted, and the difference signal is obtained, and the first zero-crossing position before and after each feature point position of the difference signal is searched; and then, the first zero-crossing position before and after each feature point position is included.
  • the time range is replaced by the low-pass filter filtering result of the higher filter frequency of the above two different frequencies, and the other part is replaced by the low-pass filter filtering result of the lower frequency of the above two different frequencies to obtain a composite signal; finally, Performing low-pass filtering on the obtained composite signal to obtain an estimated ideal ECG signal;
  • Step 3 According to the pre-processed ECG signal, the feature point position and the estimated ideal ECG signal, the spline method is used to obtain an accurate estimation of the ideal ECG signal.
  • the pre-processed ECG signal obtained by the step 1 and the estimated ideal ECG signal obtained in the second step are subtracted to obtain an error signal, and the error signal is searched for the zero-crossing position; then, the error signal obtained in the search is over.
  • the sample strip weight is 1 and the others are 0.
  • the cubic smooth spline is used to obtain accurate estimation. Ideal ECG signal;
  • Step 4 Subtracting the pre-processed ECG signal obtained by the step 1 and accurately estimating the ideal ECG signal obtained by the third step, and the subtraction result is filtered by a band pass filter to obtain a filtering result; according to the filtering The result and the QRS range obtained in step 2, the abnormal potential in the QRS is obtained by moving the standard deviation analysis technique;
  • Step 5 Perform credibility evaluation on the abnormal potential in the obtained QRS. Using the standard deviation analysis method, the credibility of the abnormal potential in the QRS obtained in the fourth step is evaluated, and it is confirmed whether the abnormal potential in the QRS obtained in the fourth step is credible, and the evaluation result is output.
  • Step 1 is specifically as follows: Participating in FIG. 3 and FIG. 4, in FIG. 3, 101 is to remove the baseline drift, and the baseline drift cancellation method is used to remove the baseline drift, and the removed signal x 1 (i) is obtained.
  • x 1 (i) in Figure 4 is designed with a third-order, 1 Hz Butterworth digital filter to design low-pass filter parameters.
  • the x(i) is obtained by bidirectional zero phase filtering.
  • 102 is to remove power frequency interference, and according to the obtained signal x 1 (i), the digital power frequency trap is used to obtain the signal x 2 (i) after removing the power frequency interference and the baseline drift.
  • x 2 (i) in Fig. 4 is a result obtained by using a digital power frequency trap.
  • Step 2 is specifically as follows: Participating in FIG. 5 and FIG. 6, wherein FIG. 6 is a partial view in the vicinity of the QRS section.
  • 201 is a feature point detection algorithm, and the ECG feature point detection is performed by using x 2 (i), and the QRS range (starting position QRS b , ending position QRS e ) and ECG feature point position p(j) are obtained, and the number of feature points is obtained.
  • Is M, j 1, 2, ..., M.
  • the ECG feature points should include at least a QRS start point, a QRS end point, and Q, R, S waveform peak points.
  • 202 is a low-pass filter that uses a higher-frequency low-pass filter to filter x 2 (i) to obtain x h (i), f h is the 202 low-pass filter filter frequency, 100 Hz ⁇ f h ⁇ 200Hz.
  • 203 is a low-pass filter that uses a lower-frequency low-pass filter to filter x 2 (i) to obtain x l (i), f l is 203 low-pass filter filter frequency, 40 Hz ⁇ f l ⁇ 80Hz.
  • x d (i) in Fig. 6 is the result of signal subtraction.
  • 205 is the zero-crossing detection and the set set h calculation.
  • the point set set(j) is constructed according to p b (j) and p f (j):
  • 206 is a signal synthesis, and according to the set set h , the signal x s (i) is synthesized by the formula (6),
  • x s (i) in Fig. 6 is the resulting composite signal.
  • 207 is a low-pass filter. After filtering x s (i), the estimated ideal ECG signal x 3 (i) is obtained, and f 3 is the 207 low-pass filter filter frequency, 100 Hz ⁇ f 3 ⁇ 200 Hz. .
  • Step 3 is specifically as follows: Participating in FIG. 7 and FIG. 8, wherein FIG. 8 is a partial view near the QRS section.
  • 301 is a signal subtraction operation, and the error signal x e (i) is calculated using equation (7):
  • x e (i) in Fig. 8 is the result of signal subtraction.
  • w(i) in Fig. 8 is a spline weight obtained according to the formula (7).
  • 303 is a cubic smoothing spline operation, and based on the estimated ideal electrocardiographic signal x 3 (i) and the spline weight w(i), the cubic smooth spline is used to obtain an accurately estimated ideal electrocardiographic signal x 4 (i ),
  • Step 4 is specifically as follows: Participate in FIG. 7, FIG. 9, and FIGS.
  • 304 is a signal subtraction operation, and the difference signal e(i) is calculated by using equation (9):
  • 305 is a band pass filter, and by performing band pass filtering on e(i), a signal y(i) containing an abnormal potential in the QRS to be extracted is obtained, and f 1 and f 2 are respectively low frequencies of the band pass filter. For frequency and high frequency low frequencies, f 1 and f 2 can be selected for specific subsequent applications.
  • 401 is to calculate the moving window variance msd(i) for the signal y(i), and calculate msd(i) using the formula (10):
  • the window length is 2k+1, and k is generally in the range of 2ms to 5ms.
  • QRS QRS start position of the front section B of about 100ms to define B as a reference QRS interval, calculating the first reference interval msd (i) the mean and standard deviation values ref_msd mean ref_msd std, respectively, using equation (11) is calculated ref msd:
  • Ref msd ref_msd mean + ⁇ *ref_msd std , (11)
  • ⁇ is generally selected to be greater than 2, and the specific selection may be determined according to actual conditions.
  • Fig. 9 is to determine the start position AIQP b of the AIQPs based on ref msd .
  • the specific method is to start the forward search from the QRS starting position QRS b . If the duration of msd(i)>ref msd is greater than or equal to m, the search is stopped, and the position at this time is set to t b , and the starting position of the AIQPs is AIQP b. Equation (12) is calculated:
  • 404 is to determine whether AIQP b is searched. If AIQP b is equal to 0, it exits and returns a failure flag, otherwise it continues.
  • Fig. 9 is to determine the end position AIQP e of the AIQPs based on ref msd .
  • the specific method is to start the reverse search from about 50 ms after the QRS termination position QRS e . If the duration of msd(i)>ref msd is greater than or equal to m, the search is stopped, and the position at this time is set to t e , the end of the AIQPs
  • the position AIQP e is calculated according to formula (13):
  • 406 is to determine whether AIQP e is searched, and if AIQP e is equal to 0, exit and return a failure flag, otherwise continue.
  • AIQP(i) is the extracted abnormal potential AIQPs within QRS, where two vertical dashed lines represent AIQP b and AIQP e , respectively .
  • Step 5 is specifically as follows: Participate in FIG. 9 and FIG. In FIG. 9, 408 is to calculate the reference interval standard deviation ref std and the QRS regional standard deviation QRS std , ref std is the standard deviation of the reference interval y(i), and the QRS std is the QRS start position QRS b to the QRS end position QRS e The standard deviation of y(i) in the interval.
  • the failure flag is returned, otherwise the success flag is returned and the extracted abnormal potential AIQP(i) in the QRS is also returned.
  • AIQP(i) was extracted by the method of the present invention, respectively.
  • the method of the present invention can extract AIQP(i) only by a single heart beat ECG, and the AIQP (i) of single heart beat ECG extraction is small when the measurement noise is small compared with the multiple superimposed average. More accurate. It can also be seen from Fig. 11 that as the average number of superpositions increases, the interference of the reference interval in step 4 of the present invention becomes smaller and smaller. Therefore, if the ECG signal obtained by the measurement is largely interfered with, the extraction from the single heart beat ECG fails, and the multiple heart beat ECGs can be superimposed and averaged.
  • Figure 12 is a result of AIQP(i) extraction of a single heart beat ECG of a patient with myocardial infarction by the method of the present invention
  • Figure 13 is a single heart beat ECG of a heart healthy person by the method of the present invention.
  • the result of AIQP(i) extraction Comparing Fig. 12 and Fig. 13, it can be found that the amplitude of AIQP(i) in patients with myocardial infarction is significantly greater than that in healthy subjects, and the morphology of both AIQP(i) is also significantly different. These features can be used for sudden cardiac death. Early warning.

Abstract

一种准确提取QRS内异常电位方法，在理想心电信号预估阶段，采用非线性变换技术，预估理想心电信号；根据预估理想心电信号，利用样条方法对理想心电信号进一步估计，精确估计理想心电信号；根据精确估计出的理想心电信号，通过移动标准差分析技术，可精确提取QRS内异常电位；该方法不但可用于多次叠加平均心电信号，也可用于单次搏动心电信号。

Description

一种准确提取QRS内异常电位的方法 技术领域
本发明属于医学信号处理领域，特别涉及一种从心电信号中准确提取QRS内异常电位的方法。
背景技术
据统计，我国每年心脏猝死(SCD)的总人数高达50多万，平均每分钟就有3人因心脏原因在发病1小时内死亡，而抢救成功率却不到1％。心脏猝死出现明显的年轻化趋势，部分患者发生猝死事件往往无明显先兆症状，身体貌似健康。所以，对心脏猝死进行早期预警显得尤为重要，但目前缺乏有效的检测手段。
大量研究表明，心脏内局部区域出现的电传导延迟是心脏猝死的一个重要诱因，这种传导延迟会引起信号折返，进而会导致严重的室性心率失常。区域传导延迟可以在体表心电图(ECG)的QRS复合波中出现切迹或不易察觉的小幅波动。这些特征往往在常规心电图中无法明确反映，需要使用采样频率更高的高频心电图。目前，临床上无创检测这种心室去极化过程出现的异常特征主要包括：(1)基于信号平均心电图(SAECG)的心室晚电位(VLP)，(2)QRS内异常电位(AIQPs)，(3)碎裂QRS波(fQRS)。
心室晚电位检查是指体表信息叠加心电图于QRS复合波终末并延伸到ST段内的高频低振幅碎裂电位，它反映了缺血区心肌迟发的电活动，是由于缺血区心肌内缓慢而不规则的拆返活动引起。在急性心肌梗死后猝死的预测中，晚电位检查具有重要价值。
临床上使用较为广泛的是时域VLP检测方法，虽然这种分析具有高的阴性预测值，但阳性预测效果不佳。VLP检测方法阳性预测率不高的主要原因是仅检测QRS复合波终末并延伸到ST段内的高频低振幅碎裂电位。大量的动物心肌梗死模型和人体标测基础研究结果明确表明，这些高频低振幅碎裂电位不但存在于QRS复合波终末区域，而且也可能隐藏在QRS复合波内。在某些心肌梗死部位，产生的高频低振幅碎裂电位可能仅存在于QRS复合波内，并不反映在QRS复合波终末处，如果能准确提取QRS复合波内的这些高频低振幅碎裂电位,将可显著提高心脏猝死早期预警的可靠性。
QRS复合波内的这些高频低振幅碎裂电位又称为QRS内异常电位(AIQPs),提取AIQPs是一个极具挑战性的任务，因为AIQPs是嵌入在QRS波内的、非常微弱、且快速变化不可预测信号。目前，提出了多种方法来应对这一挑战。离散余弦变换(DCT)域的自回归移动平均(ARMA)模型已被用于估计AIQPs，其基本思想是用一个低阶ARMA模型来模拟正常QRS波，从而可提取出不可预测的AIQPs。也有学者提出利用小波变换分析QRS波内的高频成分，并以此为基础诊断恶性室性心律失常。很多提取方法使用的是线性模型或线性变换技术，由于人体心脏工作机制十分复杂和精细，将QRS波群建模为非线性信号可能更接近实际情况。有学者提出一种用具有相同平滑度的径向基函数非线性神经网络逼近QRS波的方法，取得了较好的效果，但主要缺点是需要调整的参数过多，神经元最佳数量因人而异，参数设置不当会高估或低估RBF神经网络逼近误差，影响提取精度。
大量研究表明，对室性心律失常高风险患者，AIQPs参数可显著提高诊断准确率，但目前AIQPs提取技术的提取精度仍不能满足要求，而且提取技术的鲁棒性也需要进一步提高。
在心电测量过程中采样得到的原始心电信号用x(i)表示，其中包含工频干扰p(i)、基线漂移b(i)、测量噪声n(i)和含有待提取QRS内异常电位(AIQPs)的AIQP(i)；不包含上述成分的理想心电信号用x p(i)表示；
x(i)＝x p(i)+AIQP(i)+p(i)+b(i)+n(i)，  (1)
其中，AIQP(i)在AIQPs有效区间内为待提取QRS内异常电位，其他部分取值为0。
x(i)去除工频干扰和基线漂移后的信号用x 2(i)表示。
利用公式(2)可得到包含待提取QRS内异常电位和噪声干扰的信号y(i)：
y(i)＝x 2(i)-x p(i)＝AIQP(i)+n(i)，  (2)
在AIQPs有效区间内，当n(i)的标准差与y(i)的标准差相比较小时，可以认为y(i)≈AIQP(i)，即为待提取QRS内异常电位。问题的关键是如何尽可能准确得到理想心电信号x p(i)、如何确定AIQPs有效区间以及如何定量评价提取结果是否可信。
发明内容
为了克服上述现有技术的缺陷，本发明的目的在于提供一种准确提取QRS内异常电位的方法，采用非线性变换预估技术，结合样条方法，获得不包含QRS内异常电位和其他干扰成分的理想心电信号x p(i)，并最终提取QRS内异常电位。
为了达到上述目的，本发明的技术方案为：
一种准确提取QRS内异常电位的方法，包括以下步骤：
步骤一、对原始心电信号x 1(i)；进行预处理，得到预处理心电信号x 2(i)；当原始心电信号是测量得到的单次心搏心电信号时，对其利用低通滤波器和工频陷波器进行处理，消除基线漂移和工频干扰对后续过程的影响；当原始心电信号是测量得到的包含多个心搏的心电信号，对其利用信号平均技术进行处理，消除基线漂移、工频干扰和测量噪声对后续过程的影响。
步骤二、对预处理心电信号x 2(i)进行特征点检测，确定特征点位置和QRS范围，利用非线性变换得到预估理想心电信号；首先，对预处理心电信号进行特征点检测，确定特征点位置和QRS范围；其次，分别用两个不同滤波频率的低通滤波器，对经步骤一处理后得到的预处理心电信号进行滤波；然后，对得到的这两个滤波结果相减，得到差值信号，搜索该差值信号每个特征点位置前、后的第一个过零点位置；再后，每个特征点位置前、后的第一个过零点位置包含的时间范围，用上述两个不同频率中较高滤波频率的低通滤波器滤波结果代替，其他部分用上述两个不同频率中较低频率的低通滤波器滤波结果代替，得到合成信号；最后，对得到的合成信号进行低通滤波，得到预估理想心电信号；
步骤三、根据预处理心电信号、特征点位置和预估理想心电信号，利用样条方法得到精确估计理想心电信号。首先对经步骤一处理后得到的预处理心电信号和经步骤二得到的预估理想心电信号相减，得到误差信号，搜索该误差信号过零点位置；然后，在搜索得到的误差信号过零点位置和步骤二得到的特征点位置处，取样条权重为1，其他为0；最后，根据步骤一得到的预估理想心电信号和得到的样条权重，利用三次平滑样条得到精确估计理想心电信号；
步骤四、对经步骤一处理后得到的预处理心电信号和经步骤三得到的精确估计理想心电信号相减，该相减结果经一个带通滤波器滤波后得到滤波结果；根据该滤波结果和步骤二得到的QRS范围，通过移动标准差分析技术，获取QRS内异常电位；
步骤五、对获取的QRS内异常电位进行可信度评价。利用标准差分析方法，对经步骤四获取的QRS内异常电位可信性进行评价，确认步骤四获取的QRS内异常电位是否可信，并输出评价结果。
所述的步骤二具体为：
一、利用x 2(i)进行ECG特征点检测，得到QRS范围，起始位置QRS b、终止位置QRS e和ECG特征点位置p(j)，特征点数量为M，j＝1,2,…,M，所述ECG特征点至少应包含QRS起始点、QRS终止点，以及Q、R、S波形峰值点；
二、使用较高频率的低通滤波器对x 2(i)进行滤波，得到x hi，f h为低通滤波器滤波频率，100Hz≤f h≤200Hz；
三、使用较低频率的低通滤波器对x 2(i)进行滤波，得到x l(i)，f l为低通滤波器滤波频率，40Hz≤f l≤80Hz；
四、利用公式(3)计算差值信号x d(i)：
x d(i)＝x h(i)-x l(i)，  (3)
其中，x h(i)是较高频率低通滤波器对x 2(i)的滤波结果，x l(i)是较低频率低通滤波器对x 2(i)的滤波结果。
五、根据信号x d(i)，对每个ECG特征点时间位置p(j)，j＝1,2,…,M，反向、正方向分别搜索差值信号x d(i)，得到前、后第一个过零点，分别得到对应的时间位置p b(j)和p f(j)；
六、根据p b(j)和p f(j)构成该点集合set(j)：
set(j)＝{p b(j),p b(j)+1,…,p f(j)-1,p f(j)}，  (4)
以此为基础，构成集合set h：
set h＝{set(1),set(2),…,set(M)}，  (5)
根据集合set h，通过公式(6)合成信号x s(i)，
Figure PCTCN2018116214-appb-000001
七、对x s(i)进行低通滤波后得到预估理想心电信号x 3(i)，f 3为低通滤波器滤波频率，100Hz≤f 3≤200Hz。
所述的步骤三具体为：
一、利用公式(7)计算误差信号x e(i)：
x e(i)＝x 2(i)-x 3(i)，  (7)
x 3(i)是预估理想心电信号；
二、利用公式(8)计算样条权重w(i)：
Figure PCTCN2018116214-appb-000002
三、根据预估理想心电信号x 3(i)和样条权重w(i)，利用三次平滑样条得到精确估计的理想心电信号的x 4(i)。
所述的步骤四具体为：
一、利用公式(9)计算差值信号e(i)：
e(i)＝x 2(i)-x 4(i)  (9)
x 4(i)是精确估计的理想心电信号；
二、对e(i)进行带通滤波，得到包含待提取QRS内异常电位的信号y(i)，带通滤波带宽根据具体需要选择；
三、对信号y(i)计算移动窗口方差msd(i)，利用公式(10)计算msd(i)：
Figure PCTCN2018116214-appb-000003
其中窗口长度为2k+1，k取值范围为2ms～5ms，msd(i)计算结果是k＝2ms。
四、是计算参考MSD值ref msd，把QRS起始位置QRS b前100ms到QRS b的区间定义为参考区间，先分别计算参考区间内msd(i)的均值ref_msd mean和标准差值ref_msd std，利用公式(11)计算ref msd：
ref msd＝ref_msd mean+α*ref_msd std， (11)
其中，α一般选择大于2；
五、根据ref msd确定AIQPs的开始位置AIQP b，具体方法是，从QRS起始位置QRS b开始正向搜索，如果msd(i)>ref msd持续时间大于等于预先设定常数m，则停止搜索，此时的位置设为t b，AIQPs的开始位置AIQP b按公式(12)计算：
AIQP b＝t b-m-k， (12)
其中，m一般取值为5ms；如果搜索到QRS终止位置QRS e，则AIQP b＝0并停止搜索。
六、判断是否搜索到AIQP b，如果AIQP b等于0则退出并返回失败标志，否则继续。
七、根据ref msd确定AIQPs的结束位置AIQP e，具体方法是，从QRS终止位置QRS e后的大约50ms处开始反向搜索，如果msd(i)>ref msd持续时间大于等于m，则停止搜索，此时的位置设为t e，AIQPs的结束位置AIQP e按公式(13)计算：
AIQP e＝t e+m+k  (13)；
如果搜索到AIQPs的开始位置AIQP b，则AIQP e＝0并停止搜索，判断是否搜索到AIQP e，如果AIQP e等于0则退出并返回失败标志，否则继续。
八、提取QRS内异常电位AIQP(i)，按公式(14)计算：
Figure PCTCN2018116214-appb-000004
其中，AIQP b是搜索得到的AIQPs开始位置，AIQP e是搜索得到的AIQPs结束位置。
所述的步骤五具体为：
一、计算参考区间标准差差ref std和QRS区域标准差QRS std，ref std为参考区间y(i)的标准差，QRS std为QRS起始位置QRS b到QRS终止位置QRS e的区间内y(i)的标准差。
提取结果的可信度判断按公式(15)计算：
Figure PCTCN2018116214-appb-000005
其中，β>1，具体选择可根据实际情况确定。
如果可信度等于0则返回失败标志，否则返回成功标志并同时返回提取的QRS内异常电位AIQP(i)。
有益效果
本发明提出一种利用理想心电信号二次估计技术，精确提取QRS内异常电位的方法。在理想心电信号预估阶段，采用非线性变换技术，即可对非ECG特征点区域变化趋势进行有效跟踪，又可有效消除ECG特征点对提取理想心电信号可能造成的影响。根据预估理想心电信号，利用样条方法对理想心电信号进一步估计，可精确估计出理想心电信号。与现有方法相比，本发明方法需要选择的参数少，结果更可靠。本发明还对提取的QRS内异常电位进行可信度评价，从而保证利用本发明方法进行应用开发的结果可靠性。与传统的多次叠加平均方法相比，本发明方法的一个突出特点是可对单次搏动心电信号进行QRS内异常电位提取，可大大扩展AIQPs分析技术的应用范围。
本发明方法可应用的场景和范围包括：1)在心电图机中集成本发明方法，可在进行常规心电图测量时，对患者进行心脏猝死危险性评价；2)在常规多参数监护仪中集成本发明方法，可对心肌梗死患者的病情变化进行实时动态跟踪监测；3)开发便于使用的便携式或可穿戴装置，实现在家庭环境下的心脏猝死危险性预警；4)在移动设备(如手机)中集成本发明方法，可为设备使用者提供一种便捷、高效的心脏猝死危险性预警手段。
附图说明
图1是本发明的流程图。
图2是本发明的原理描述使用的模拟信号。
图3是本发明实施例的预处理流程图。
图4是本发明实施例的基线漂移消除和工频干扰去除过程示意图。
图5是本发明实施例的预估理想心电信号流程图。
图6是本发明实施例的预估理想心电信号获取过程示意图。
图7是本发明实施例的精确估计理想心电信号流程图。
图8是本发明实施例的利用三次平滑样条得到精确估计理想心电信号过程示意图。
图9是本发明实施例的QRS内异常电位提取与结果可信度评价流程图。
图10是本发明实施例的QRS内异常电位提取与结果可信度评价的过程示意图。
图11是对单次搏动ECG和多次搏动ECG叠加平均，分别用本发明方法进行QRS内异常电位提取的结果图。
图12是用本发明所述方法对一例心肌梗死患者的单次搏动ECG进行QRS内异常电位提取的结果。
图13是用本发明所述方法对一例心脏健康者的单次搏动ECG进行QRS内异常电位提取的结果。
具体实施方式
下面结合附图对本发明的原理作详细描述。
参照图1，是本发明的流程图，一种准确提取QRS内异常电位的方法，包括以下步骤：
步骤一、对原始心电信号x 1(i)；进行预处理，得到预处理心电信号x 2(i)；当原始心电信号是测量得到的单次心搏心电信号时，对其利用低通滤波器和工频陷波器进行处理，消除基线漂移和工频干扰对后续过 程的影响；当原始心电信号是测量得到的包含多个心搏的心电信号，对其利用信号平均技术进行处理，消除基线漂移、工频干扰和测量噪声对后续过程的影响。
步骤二、对预处理心电信号x 2(i)进行特征点检测，确定特征点位置和QRS范围，利用非线性变换得到预估理想心电信号；首先，对预处理心电信号进行特征点检测，确定特征点位置和QRS范围；其次，分别用两个不同滤波频率的低通滤波器，对经步骤一处理后得到的预处理心电信号进行滤波；然后，对得到的这两个滤波结果相减，得到差值信号，搜索该差值信号每个特征点位置前、后的第一个过零点位置；再后，每个特征点位置前、后的第一个过零点位置包含的时间范围，用上述两个不同频率中较高滤波频率的低通滤波器滤波结果代替，其他部分用上述两个不同频率中较低频率的低通滤波器滤波结果代替，得到合成信号；最后，对得到的合成信号进行低通滤波，得到预估理想心电信号；
步骤三、根据预处理心电信号、特征点位置和预估理想心电信号，利用样条方法得到精确估计理想心电信号。首先对经步骤一处理后得到的预处理心电信号和经步骤二得到的预估理想心电信号相减，得到误差信号，搜索该误差信号过零点位置；然后，在搜索得到的误差信号过零点位置和步骤二得到的特征点位置处，取样条权重为1，其他为0；最后，根据步骤一得到的预估理想心电信号和得到的样条权重，利用三次平滑样条得到精确估计理想心电信号；
步骤四、对经步骤一处理后得到的预处理心电信号和经步骤三得到的精确估计理想心电信号相减，该相减结果经一个带通滤波器滤波后得到滤波结果；根据该滤波结果和步骤二得到的QRS范围，通过移动标准差分析技术，获取QRS内异常电位；
步骤五、对获取的QRS内异常电位进行可信度评价。利用标准差分析方法，对经步骤四获取的QRS内异常电位可信性进行评价，确认步骤四获取的QRS内异常电位是否可信，并输出评价结果。
图2是用于本发明原理描述的模拟心电信号x(i)，包含理想心电信号用x p(i)和模拟包含待提取QRS内异常电位AIQP s(i)，以及工频干扰p(i)、基线漂移b(i)和测量噪声n(i)，采样率1000Hz，数据长度N＝800。
所述的步骤一具体为：参加图3和图4，图3中101是去除基线漂移，利用基线漂移消除方法去除基线漂移，得到去除后的信号x 1(i)。
目前，基线漂移消除方法有多种方法，由于基线漂移对最终提取结果没有明显影响，图4中的x 1(i)是用三阶、1Hz的Butterworth数字滤波器设计低通滤波器参数，然后对x(i)用双向零相位滤波得到。
图3中102是去除工频干扰，根据得到的信号x 1(i)，利用数字工频陷波器得到去除工频干扰和基线漂移后的信号x 2(i)。图4中的x 2(i)是利用数字工频陷波器得到的结果。
所述的步骤二具体为：参加图5和图6，其中图6是QRS区间附近的局部图。图5中201是特征点检测算法，利用x 2(i)进行ECG特征点检测，得到QRS范围(起始位置QRS b、终止位置QRS e)和ECG特征点位置p(j)，特征点数量为M，j＝1,2,…,M。所述ECG特征点至少应包含QRS起始点、QRS终止点，以及Q、R、S波形峰值点。
图5中202是一个低通滤波器，使用较高频率的低通滤波器对x 2(i)进行滤波，得到x h(i)，f h为202低通滤波器滤波频率，100Hz≤f h≤200Hz。图6中的x h(i)是用三阶、f h＝150Hz的Butterworth数字滤波器设计低通滤波器参数，然后对x 2(i)用双向零相位滤波得到。
图5中203是一个低通滤波器，使用较低频率的低通滤波器对x 2(i)进行滤波，得到x l(i)，f l为203低通滤波器滤波频率，40Hz≤f l≤80Hz。图6中的x l(i)是用三阶、f l＝60Hz的Butterworth数字滤波器设计低通滤波器参数，然后对x 2(i)用双向零相位滤波得到。
图5中204是信号相减运算，利用公式(3)计算差值信号x d(i)：
x d(i)＝x h(i)-x l(i)，  (3)
图6中的x d(i)是信号相减后的结果。
图5中205是过零点检测与集合set h计算，根据信号x d(i)，对每个ECG特征点时间位置p(j)，j＝1,2,…,M，反向、正方向分别搜索差值信号x d(i)，得到前、后第一个过零点，分别得到对应的时间位置p b(j)和p f(j)。
图6中的x d(i)中，虚线为零值线，符号“O”代表时间位置p(j)处对应x d(i)上的点，符号“Δ”和
Figure PCTCN2018116214-appb-000006
分别代表时间位置p b(j)和p f(j)对应x d(i)上的点，M＝4。由于x d(i)为时间离散信号，x d(i)过零点处对应的实际值并不一定为0。
对每个p(j)，根据p b(j)和p f(j)构成该点集合set(j)：
set(j)＝{p b(j),p b(j)+1,…,p f(j)-1,p f(j)}，  (4)
以此为基础，构成集合set h：
set h＝{set(1),set(2),…,set(M)}， (5)
图5中206是信号合成，根据集合set h，通过公式(6)合成信号x s(i)，
Figure PCTCN2018116214-appb-000007
图6中的x s(i)是得到的合成信号。
图5中207是一个低通滤波器，对x s(i)进行滤波后得到预估理想心电信号x 3(i)，f 3为207低通滤波器滤波频率，100Hz≤f 3≤200Hz。图6中的x 3(i)是是用三阶、f 3＝150Hz的Butterworth数字滤波器设计低通滤波器参数，然后对x s(i)用双向零相位滤波得到的预估理想心电信号
所述的步骤三具体为：参加图7和图8，其中图8是QRS区间附近的局部图。图7中301是信号相减运算，利用公式(7)计算误差信号x e(i)：
x e(i)＝x 2(i)-x 3(i)，  (7)
图8中的x e(i)是信号相减后的结果。
图7中302是计算样条权重，利用公式(8)计算样条权重w(i)：
Figure PCTCN2018116214-appb-000008
图8中的w(i)是根据公式(7)得到的样条权重。
图7中303是三次平滑样条运算，根据预估理想心电信号x 3(i)和样条权重w(i)，利用三次平滑样条得到精确估计的理想心电信号的x 4(i)，
图8中的x 4(i)是得到的精确估计的理想心电信号。
所述的步骤四具体为：参加图7、图9和图8、10。
图7中304是信号相减运算，利用公式(9)计算差值信号e(i)：
e(i)＝x 2(i)-x 4(i)。  (9)
图7中305是一个带通滤波器，通过对e(i)进行带通滤波，得到包含待提取QRS内异常电位的信号y(i)，f 1和f 2分别为带通滤波器的低频频率和高频低频，f 1和f 2可根据具体后续应用选择。图8中的y(i)是用五阶、f 1＝70Hz、f 2＝300Hz的Butterworth数字滤波器设计带通滤波器参数，然后对e(i)用双向零相位滤波。
图9中401是对信号y(i)计算移动窗口方差msd(i)，利用公式(10)计算msd(i)：
Figure PCTCN2018116214-appb-000009
其中窗口长度为2k+1，k一般取值范围为2ms～5ms。图10中msd(i)计算结果是k＝2ms。
图9中402是计算参考MSD值ref msd。把QRS起始位置QRS b前大约100ms到QRS b的区间定义为参考区间，先分别计算参考区间内msd(i)的均值ref_msd mean和标准差值ref_msd std，利用公式(11)计算ref msd：
ref msd＝ref_msd mean+α*ref_msd std， (11)
其中，α一般选择大于2，具体选择可根据实际情况确定。
图10msd(i)中的水平虚线的幅度值代表ref msd值，α＝3。
图9中403是根据ref msd确定AIQPs的开始位置AIQP b。具体方法是，从QRS起始位置QRS b开始正向搜索，如果msd(i)>ref msd持续时间大于等于m，则停止搜索，此时的位置设为t b，AIQPs的开始位置AIQP b按公式(12)计算：
AIQP b＝t b-m-k， (12)
其中，m一般取值为5ms；如果搜索到QRS终止位置QRS e，则AIQP b＝0并停止搜索。
图9中404是判断是否搜索到AIQP b，如果AIQP b等于0则退出并返回失败标志，否则继续。
图9中405是根据ref msd确定AIQPs的结束位置AIQP e。具体方法是，从QRS终止位置QRS e后的大约50ms处开始反向搜索，如果msd(i)>ref msd持续时间大于等于m，则停止搜索，此时的位置设为t e，AIQPs的结束位置AIQP e按公式(13)计算：
AIQP e＝t e+m+k； (13)
如果搜索到AIQPs的开始位置AIQP b，则AIQP e＝0并停止搜索。
图9中406是判断是否搜索到AIQP e，如果AIQP e等于0则退出并返回失败标志，否则继续。
图9中407是提取QRS内异常电位AIQP(i)，按公式(14)计算：。
Figure PCTCN2018116214-appb-000010
图10AIQP(i)为提取的QRS内异常电位AIQPs，其中两条垂直虚线分别代表AIQP b和AIQP e。
所述的步骤五具体为：参加图9和图10。图9中408是计算参考区间标准差差ref std和QRS区域标准差QRS std，ref std为参考区间y(i)的标准差，QRS std为QRS起始位置QRS b到QRS终止位置QRS e的区间内y(i)的标准差。
图9中409是提取结果的可信度判断，可信度按公式(15)计算：
Figure PCTCN2018116214-appb-000011
其中，β>1，具体选择可根据实际情况确定。
如果可信度等于0则返回失败标志，否则返回成功标志并同时返回提取的QRS内异常电位AIQP(i)。
图11是使用本发明描述中的模拟信号产生方法，分别产生1次、10次、和100次模拟心搏ECG，然后分别对多次模拟心跳ECG进行叠加平均，对得到的1次心搏ECG、以及10次和100次模拟心搏平均ECG，分别用本发明所述方法提取AIQP(i)。1次、以及10次和100次叠加平均提取的AIQP(i)，与模拟待提取AIQP s(i)相比，相关系数和均方误差分别为：0.95、0.1；0.87、0.32；0.87、0.29。结果表明，本发明方法仅需要单次心搏ECG就可提取AIQP(i)，与多次叠加平均相比，在测量噪声足有小的情况下，单次心搏ECG提取的AIQP(i)准确度更高。从图11中也可看出，随着叠加平均次数增加，本发明步骤四中参考区间的干扰越来越小。所以，如果由于测量得到的ECG信号干扰较大，导致从单次心搏ECG中提取失败，可对多个心搏ECG进行叠加平均提取。
图12是用本发明所述方法对一例心肌梗死患者的单次心搏ECG进行AIQP(i)提取的结果，图13是用本发明所述方法对一例心脏健康者的单次心搏ECG进行AIQP(i)提取的结果。对比图12和图13可以发现，心肌梗死患者的AIQP(i)幅度明显大于心脏健康者的AIQP(i)幅度，且两者AIQP(i)的形态也有明显差异，这些特征可用于心脏猝死的早期预警。

Claims (5)

  1. 一种准确提取QRS内异常电位的方法，其特征在于，包括以下步骤：
    步骤一、对原始心电信号x 1(i)；进行预处理，得到预处理心电信号x 2(i)；当原始心电信号是测量得到的单次心搏心电信号时，对其利用低通滤波器和工频陷波器进行处理，消除基线漂移和工频干扰对后续过程的影响；当原始心电信号是测量得到的包含多个心搏的心电信号，对其利用信号平均技术进行处理，消除基线漂移、工频干扰和测量噪声对后续过程的影响；
    步骤二、对预处理心电信号x 2(i)进行特征点检测，确定特征点位置和QRS范围，利用非线性变换得到预估理想心电信号；首先，对预处理心电信号进行特征点检测，确定特征点位置和QRS范围；其次，分别用两个不同滤波频率的低通滤波器，对经步骤一处理后得到的预处理心电信号进行滤波；然后，对得到的这两个滤波结果相减，得到差值信号，搜索该差值信号每个特征点位置前、后的第一个过零点位置；再后，每个特征点位置前、后的第一个过零点位置包含的时间范围，用上述两个不同频率中较高滤波频率的低通滤波器滤波结果代替，其他部分用上述两个不同频率中较低频率的低通滤波器滤波结果代替，得到合成信号；最后，对得到的合成信号进行低通滤波，得到预估理想心电信号；
    步骤三、根据预处理心电信号、特征点位置和预估理想心电信号，利用样条方法得到精确估计理想心电信号。首先对经步骤一处理后得到的预处理心电信号和经步骤二得到的预估理想心电信号相减，得到误差信号，搜索该误差信号过零点位置；然后，在搜索得到的误差信号过零点位置和步骤二得到的特征点位置处，取样条权重为1，其他为0；最后，根据步骤一得到的预估理想心电信号和得到的样条权重，利用三次平滑样条得到精确估计理想心电信号；
    步骤四、对经步骤一处理后得到的预处理心电信号和经步骤三得到的精确估计理想心电信号相减，该相减结果经一个带通滤波器滤波后得到滤波结果；根据该滤波结果和步骤二得到的QRS范围，通过移动标准差分析技术，获取QRS内异常电位；
    步骤五、对获取的QRS内异常电位进行可信度评价。利用标准差分析方法，对经步骤四获取的QRS内异常电位可信性进行评价，确认步骤四获取的QRS内异常电位是否可信，并输出评价结果。
  2. 根据权利要求1所述的一种准确提取QRS内异常电位的方法，其特征在于，
    所述的步骤二具体为：
    一、利用x 2(i)进行ECG特征点检测，得到QRS范围，起始位置QRS b、终止位置QRS e和ECG特征点位置p(j)，特征点数量为M，j＝1,2,…,M，所述ECG特征点至少应包含QRS起始点、QRS终止点，以及Q、R、S波形峰值点；
    二、使用较高频率的低通滤波器对x 2(i)进行滤波，得到x hi，f h为低通滤波器滤波频率，100Hz≤f h≤200Hz；
    三、使用较低频率的低通滤波器对x 2(i)进行滤波，得到x l(i)，f l为低通滤波器滤波频率，40Hz≤f l≤80Hz；
    四、利用公式(3)计算差值信号x d(i)：
    x d(i)＝x h(i)-x l(i)，   (3)
    其中，x h(i)是较高频率低通滤波器对x 2(i)的滤波结果，x l(i)是较低频率低通滤波器对x 2(i)的滤波结果。
    五、根据信号x d(i)，对每个ECG特征点时间位置p(j)，j＝1,2,…,M，反向、正方向分别搜索差值信号x d(i)，得到前、后第一个过零点，分别得到对应的时间位置p b(j)和p f(j)；
    六、根据p b(j)和p f(j)构成该点集合set(j)：
    set(j)＝{p b(j),p b(j)+1,…,p f(j)-1,p f(j)}，   (4)
    以此为基础，构成集合set h：
    set h＝{set(1),set(2),…,set(M)}，   (5)
    根据集合set h，通过公式(6)合成信号x s(i)，
    Figure PCTCN2018116214-appb-100001
    七、对x s(i)进行低通滤波后得到预估理想心电信号x 3(i)，f 3为低通滤波器滤波频率，100Hz≤f 3≤200Hz。
  3. 根据权利要求1所述的一种准确提取QRS内异常电位的方法，其特征在于，所述的步骤三具体为：
    一、利用公式(7)计算误差信号x e(i)：
    x e(i)＝x 2(i)-x 3(i)，    (7)
    x 3(i)是预估理想心电信号；
    二、利用公式(8)计算样条权重w(i)：
    Figure PCTCN2018116214-appb-100002
    三、根据预估理想心电信号x 3(i)和样条权重w(i)，利用三次平滑样条得到精确估计的理想心电信号的x 4(i)。
  4. 根据权利要求1所述的一种准确提取QRS内异常电位的方法，其特征在于，所述的步骤四具体为：
    一、利用公式(9)计算差值信号e(i)：
    e(i)＝x 2(i)-x 4(i)   (9)
    x 4(i)是精确估计的理想心电信号；
    二、对e(i)进行带通滤波，得到包含待提取QRS内异常电位的信号y(i)，带通滤波带宽根据具体需要选择；
    三、对信号y(i)计算移动窗口方差msd(i)，利用公式(10)计算msd(i)：
    Figure PCTCN2018116214-appb-100003
    其中窗口长度为2k+1，k取值范围为2ms～5ms，msd(i)计算结果是k＝2ms。
    四、是计算参考MSD值ref msd，把QRS起始位置QRS b前100ms到QRS b的区间定义为参考区间，先分别计算参考区间内msd(i)的均值ref_msd mean和标准差值ref_msd std，利用公式(11)计算ref msd：
    ref msd＝ref_msd mean+α*ref_msd std，  (11)
    其中，α一般选择大于2；
    五、根据ref msd确定AIQPs的开始位置AIQP b，具体方法是，从QRS起始位置QRS b开始正向搜索，如果msd(i)>ref msd持续时间大于等于预先设定常数m，则停止搜索，此时的位置设为t b，AIQPs的开始位置AIQP b按公式(12)计算：
    AIQP b＝t b-m-k，  (12)
    其中，m一般取值为5ms；如果搜索到QRS终止位置QRS e，则AIQP b＝0并停止搜索。
    六、判断是否搜索到AIQP b，如果AIQP b等于0则退出并返回失败标志，否则继续。
    七、根据ref msd确定AIQPs的结束位置AIQP e，具体方法是，从QRS终止位置QRS e后的大约50ms处开始反向搜索，如果msd(i)>ref msd持续时间大于等于m，则停止搜索，此时的位置设为t e，AIQPs的结束位置AIQP e按公式(13)计算：
    AIQP e＝t e+m+k     (13)；
    如果搜索到AIQPs的开始位置AIQP b，则AIQP e＝0并停止搜索，判断是否搜索到AIQP e，如果AIQP e等于0则退出并返回失败标志，否则继续。
    八、提取QRS内异常电位AIQP(i)，按公式(14)计算：
    Figure PCTCN2018116214-appb-100004
    其中，AIQP b是搜索得到的AIQPs开始位置，AIQP e是搜索得到的AIQPs结束位置。
  5. 根据权利要求1所述的一种准确提取QRS内异常电位的方法，其特征在于，所述的步骤五具体为：
    一、计算参考区间标准差差ref std和QRS区域标准差QRS std，ref std为参考区间y(i)的标准差，QRS std为QRS起始位置QRS b到QRS终止位置QRS e的区间内y(i)的标准差。
    提取结果的可信度判断按公式(15)计算：
    Figure PCTCN2018116214-appb-100005
    其中，β>1，具体选择可根据实际情况确定。
    如果可信度等于0则返回失败标志，否则返回成功标志并同时返回提取的QRS内异常电位AIQP(i)。
PCT/CN2018/116214 2018-01-22 2018-11-19 一种准确提取qrs内异常电位的方法 WO2019140998A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US16/609,174 US10912479B2 (en) 2018-01-22 2018-11-19 Method for accurately extracting abnormal potential within QRS

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201810061071.0A CN108523881B (zh) 2018-01-22 2018-01-22 一种准确提取qrs内异常电位的方法
CN201810061071.0 2018-01-22

Publications (1)

Publication Number Publication Date
WO2019140998A1 true WO2019140998A1 (zh) 2019-07-25

Family

ID=63485750

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2018/116214 WO2019140998A1 (zh) 2018-01-22 2018-11-19 一种准确提取qrs内异常电位的方法

Country Status (3)

Country Link
US (1) US10912479B2 (zh)
CN (1) CN108523881B (zh)
WO (1) WO2019140998A1 (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113705354A (zh) * 2021-07-30 2021-11-26 西安交通大学 一种混合式心电信号运动伪迹去除方法、系统、设备及可读存储介质

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108523881B (zh) 2018-01-22 2021-09-03 西安交通大学 一种准确提取qrs内异常电位的方法
CN111466905B (zh) * 2020-04-10 2021-01-22 西安交通大学 一种基于双向连通的心电波形提取方法
CN113509187B (zh) * 2021-07-12 2023-10-24 广州市康源图像智能研究院 一种实时的单导联心电r波检测方法
CN114742114B (zh) * 2022-06-09 2022-10-21 毕胜普生物科技有限公司 高频qrs波形曲线分析方法、装置、计算机设备与存储介质
CN115422976B (zh) * 2022-09-14 2023-11-21 湖南万脉医疗科技有限公司 一种基于人工网络的心肺耦合关系分析方法及监测系统

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1951320A (zh) * 2006-11-16 2007-04-25 上海交通大学 基于超完备特征的异常心电识别方法
US20090076403A1 (en) * 2007-09-14 2009-03-19 Bruce Hopenfeld Waveform feature value averaging system and methods for the detection of cardiac events
CN103705234A (zh) * 2013-12-05 2014-04-09 深圳先进技术研究院 动态心电信号数据中的检波方法和装置
WO2016183683A1 (en) * 2015-05-20 2016-11-24 University Health Network Method and system for assessing qrs components and the risk of ventricular arrhythmias
WO2017213533A1 (en) * 2016-06-07 2017-12-14 Institutul National De Cercetare-Dezvoltare Pentru Fizica Tehnica - Ift Iasi Apparatus for determining abnormal electric potentials in ventricular myocardium
CN108523881A (zh) * 2018-01-22 2018-09-14 西安交通大学 一种准确提取qrs内异常电位的方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6132381A (en) * 1997-08-14 2000-10-17 Agilent Technologies, Inc. Intramyocardial anomalous activity detection by subtracting modeled respiratory effect
WO2006123334A2 (en) * 2005-05-17 2006-11-23 Bio Signal Analysis Ltd. Analysis of electrocardiogram signals
US9408549B2 (en) * 2009-11-03 2016-08-09 Vivaquant Llc Detecting fiducial points in physiological signals
WO2012106729A1 (en) * 2011-02-04 2012-08-09 Phase Space Systems Corporation System and method for evaluating an electrophysiological signal
FR2988499B1 (fr) * 2012-03-23 2014-12-19 Univ Bretagne Occidentale Procede de detection d'au moins une anomalie dans un signal observe, produit programme d'ordinateur et dispositif correspondants

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1951320A (zh) * 2006-11-16 2007-04-25 上海交通大学 基于超完备特征的异常心电识别方法
US20090076403A1 (en) * 2007-09-14 2009-03-19 Bruce Hopenfeld Waveform feature value averaging system and methods for the detection of cardiac events
CN103705234A (zh) * 2013-12-05 2014-04-09 深圳先进技术研究院 动态心电信号数据中的检波方法和装置
WO2016183683A1 (en) * 2015-05-20 2016-11-24 University Health Network Method and system for assessing qrs components and the risk of ventricular arrhythmias
WO2017213533A1 (en) * 2016-06-07 2017-12-14 Institutul National De Cercetare-Dezvoltare Pentru Fizica Tehnica - Ift Iasi Apparatus for determining abnormal electric potentials in ventricular myocardium
CN108523881A (zh) * 2018-01-22 2018-09-14 西安交通大学 一种准确提取qrs内异常电位的方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113705354A (zh) * 2021-07-30 2021-11-26 西安交通大学 一种混合式心电信号运动伪迹去除方法、系统、设备及可读存储介质
CN113705354B (zh) * 2021-07-30 2024-04-05 西安交通大学 一种混合式心电信号运动伪迹去除方法、系统、设备及可读存储介质

Also Published As

Publication number Publication date
CN108523881B (zh) 2021-09-03
US10912479B2 (en) 2021-02-09
CN108523881A (zh) 2018-09-14
US20200337582A1 (en) 2020-10-29

Similar Documents

Publication Publication Date Title
WO2019140998A1 (zh) 一种准确提取qrs内异常电位的方法
Sathyapriya et al. Analysis and detection R-peak detection using Modified Pan-Tompkins algorithm
Stridh et al. Sequential characterization of atrial tachyarrhythmias based on ECG time-frequency analysis
Choudhary et al. Automatic detection of aortic valve opening using seismocardiography in healthy individuals
Pal et al. Empirical mode decomposition based ECG enhancement and QRS detection
Campolo et al. ECG-derived respiratory signal using empirical mode decomposition
Varanini et al. A multi-step approach for non-invasive fetal ECG analysis
Sasikala et al. Extraction of P wave and T wave in Electrocardiogram using Wavelet Transform
JP7232194B2 (ja) Ｅｍｇ信号を処理するための方法及びシステム
Rahul et al. Dynamic thresholding based efficient QRS complex detection with low computational overhead
TW201717845A (zh) 從具雜訊之心電圖資料決定心跳率的設備、系統和方法
Bae et al. ECG PQRST complex detector and heart rate variability analysis using temporal characteristics of fiducial points
Choudhary et al. Orthogonal subspace projection based framework to extract heart cycles from SCG signal
Song et al. New real-time heartbeat detection method using the angle of a single-lead electrocardiogram
Lin et al. Investigation on pulse wave forward peak detection and its applications in cardiovascular health
Chen et al. An ECG R-wave detection algorithm based on adaptive threshold
Murthy et al. Ecg signal denoising and ischemic event feature extraction using daubechies wavelets
KR20140114181A (ko) 심전도 신호에 기반하여 스트레스를 분석하고 추정하는 방법 및 장치
Elmansouri et al. New Electrocardiogram signal analysis in a research laboratory using LabVIEW
Zhao et al. An early warning of atrial fibrillation based on short-time ECG signals
Jindal et al. MATLAB based GUI for ECG arrhythmia detection using Pan-Tompkin algorithm
Illanes-Manriquez An automatic multi-lead electrocardiogram segmentation algorithm based on abrupt change detection
CN115024716A (zh) 基于心率标签生成的心冲击图信号重建方法
Zhuang et al. Improved ECG-derived respiration using empirical wavelet transform and Kernel principal component analysis
TWI696191B (zh) 降低雜訊影響之qrs波偵測演算法

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 18901316

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 18901316

Country of ref document: EP

Kind code of ref document: A1