Electrocardiosignal denoising method, equipment and storage medium based on BP neural network and improved EMD method
阅读说明：本技术 一种基于bp神经网络及改进的emd方法的心电信号去噪方法、设备及存储介质 (Electrocardiosignal denoising method, equipment and storage medium based on BP neural network and improved EMD method ) 是由 李康 陈阳 于 2021-07-29 设计创作，主要内容包括：本发明涉及一种基于BP神经网络及改进的EMD方法的心电信号去噪方法、设备及存储介质,本发明通过向原始信号添加高斯白噪声的方法解决了EMD分解方法中重构信号含有残余噪声的问题,通过BP神经网络模型来评估IMF分量中信号与噪声的比重并将其分类,克服了传统IMF分量分类方法指标单一的缺点,能够更准确地区分哪些分量需要进一步去噪,哪些分量可以直接用于信号重构,哪些信号需要舍弃,保证了噪声去除的同时有用信号不丢失,进而使得重构后的信号更加准确。(The invention relates to an electrocardiosignal denoising method, equipment and a storage medium based on a BP neural network and an improved EMD method, and solves the problem that a reconstructed signal contains residual noise in an EMD decomposition method by adding white Gaussian noise to an original signal.)
1. An electrocardiosignal denoising method based on a BP neural network and an improved EMD method is characterized by comprising the following steps:
s1, selecting a plurality of electrocardiosignal samples from the MIT-BIH database;
s2, decomposing all electrocardiosignals by using an improved EMD method, decomposing each electrocardiosignal into a plurality of IMF components, marking the IMF components obtained by decomposing the same electrocardiosignal as one group, randomly dividing the finally obtained IMF components into two groups, taking the large group with more number as a training set, and taking the rest IMF components as a test set;
s3, calculating the correlation coefficient of each IMF component and the corresponding electrocardiosignal before decomposition, drawing a correlation coefficient distribution graph, and dividing the IMF components into three categories of noise leading, signal leading and invalid components according to the correlation coefficient;
s4, calculating five characteristic parameters of the IMF component, including margin factor, kurtosis coefficient, baseline energy ratio, QRS energy ratio and peak-to-average ratio;
s5, constructing and training a BP neural network model;
s6, classifying the IMF components by using the trained BP neural network model, reserving the IMF components leading the signals, directly discarding invalid components, and performing next threshold denoising on the components marked as leading noise;
s7, further denoising the noise-dominant component by adopting a wavelet threshold method;
and S8, adding the denoised noise dominant component and the preserved electrocardiosignal dominant component to reconstruct the electrocardiosignal.
2. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein in step S2, the improved EMD method is as follows:
s21: adding noise with amplitude distribution obeying Gaussian distribution and power spectral density distribution obeying uniform distribution into the electrocardiosignal obtained in the step S1;
s22: decomposing the electrocardiosignal subjected to the noise addition in the step S21 to obtain a first IMF component, adding m kinds of noises with different noise intensities in the step S21 according to the principle that the noise intensity is gradually decreased until m times of EMD decomposition are completed, wherein the formula of the EMD decomposition is shown as a formula (I):
in formula (I), d (n) is the original electrocardiosignal, IMFk (n) is the k-th IMF component after decomposition, and ak (n) is the remainder of the signal after k-th decomposition.
3. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein in step S3, the IMF component classification method is as follows:
s31: calculating the correlation coefficient rho of each IMF component and the original signal, namely the corresponding electrocardiosignal before decompositionxyThe calculation formula is shown as formula (II):
in formula (II), x (i) is the signal of the correlation coefficient to be calculated, and y (i) is the original signal;
s32: drawing a correlation coefficient distribution graph; in the correlation coefficient distribution diagram, the abscissa is an IMF component, and the ordinate is a correlation coefficient corresponding to the IMF component;
s33: the correlation coefficient distribution map has n points, which are denoted as 1-n, and the corresponding IMF components are denoted as IMF1～IMFnFinding the first turning point in the correlation coefficient distribution diagram and recording the turning point as P, and recording the corresponding IMF component as IMFP，IMF1To IMFPAll IMF components in the interval are marked as noise-dominant components, and then further denoising is needed; from IMFPTo IMFnAnd finding out IMF components with the correlation coefficient smaller than 0.15 from the IMF components in between, recording the IMF components as invalid components, removing the invalid components, and recording the rest components as signal dominant components.
4. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein in step S4, the margin factor is a statistic for reflecting the interference bearing capacity of the signal, the larger the margin factor is, the larger the proportion of the useful information is, the smaller the interference degree is, and t is1To express the margin factor, the calculation formula of the margin factor is shown as formula (iii):
in the formula (III), xpIs the peak value, xiIs the ith signal sampling value, and n is the number of sampling points;
in step S4, the kurtosis coefficient is a physical quantity used to represent the steepness of the top of the data curve; the crest factor of a clean and intact electrocardiosignal is more than 5, and if myoelectric interference, baseline drift, power frequency interference or Gaussian distribution random noise exists, the crest factor of the electrocardiosignal is lower than 5;
kurtosis coefficient t2The formula (iv) is shown as formula (iv):
in the formula (IV), xiFor the ith signal sample value,is the signal sample average;
in step S4, the baseline energy ratio is a characteristic parameter used to reflect the baseline interference degree, the characteristic is defined as a quotient between energy of a 1-40Hz frequency band and energy of a 0-40Hz frequency band, the frequency of baseline drift is 0.15-0.3Hz, the baseline energy ratio effectively represents that there is enough baseline interference with a large influence, and the smaller the baseline energy ratio, the more the baseline interference; base line energy ratio t3The formula (c) is shown in formula (v):
in the formula (V), P (f) is a function of the power spectral density of the electrocardiosignal, and f is frequency;
in step S4, the QRS energy ratio is a physical quantity for reflecting the degree of electromyographic interference; the characteristic is defined as the ratio of the energy of 5-15Hz frequency band to the energy of 5-40Hz frequency band; QRS energy ratio t4The formula (c) is shown in formula (vi):
in step S4, the peak-to-average ratio is a parameter reflecting the degree of difference between the data peak value and the data mean value; the calculation formula of the peak-to-average ratio is shown as the formula (VII):
5. the electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein the BP neural network model comprises an input layer, a hidden layer and an output layer;
the input layer has 5 neurons;
the output layer has 3 neurons; three categories of IMF components are used as output vectors of the BP neural network model, and the three categories of the IMF components comprise noise leading components, signal leading components and invalid components which are respectively marked as 0, 1 and 2;
number of neurons in hidden layerm is the number of neurons of the hidden layer, n is the number of neurons of the input layer, l is the number of neurons of the output layer, and alpha is 1-10;
in the BP neural network model, the activation functions of the hidden layer and the output layer are selected from an S function, and the formula of the S function f (x) is shown as the formula (VIII):
further preferably, the number of neurons of the hidden layer is 4.
6. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein in step S5, the training process of the BP neural network model is as follows:
(1) normalization treatment: normalizing the five characteristic parameters obtained in the step S4;
(2) inputting: inputting the characteristic parameters after normalization processing into a BP neural network model as input vectors;
(3) and a forward propagation stage: the neurons of the training set transmitted into the input layer interact with the weight threshold value and then are transmitted to the neurons of the hidden layer, the neurons are transmitted into the next layer through the calculation of the activation function, and the neurons are gradually transmitted to the output layer by layer in the same way;
(4) and (3) a back propagation stage: comparing the result of the output layer with the expected result to obtain the error of the output result, transmitting the error layer by layer in a reverse direction to obtain the deviation of the threshold value and the weight value of each layer, and adjusting the threshold value and the weight value through the deviation;
(5) classifying IMF components of a training set by using a trained BP neural network model, evaluating the BP neural network model, wherein the classification accuracy of the BP neural network is 95%, and if the classification accuracy of the BP neural network is not 95%, continuing training, specifically: and calculating characteristic parameters of the IMF components in the training set, classifying the characteristic parameters as the input of the BP neural network model, and calculating the accuracy of the BP neural network model by comparing the classification result of the BP neural network model with the class label of the IMF components.
7. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in claim 1, wherein in step (3), the forward propagation stage is realized by the following steps:
a. computing hidden layer output HjThe formula is shown as formula (IX):
in the formula (IX), HjFor the output of the hidden layer, f () is the activation function of the hidden layer, vijIs the connection weight, x, of the input layer and the hidden layeriFor input layer input, ajIs a threshold value; i is the number of neurons in the input layer, and j is the number of neurons in the hidden layer;
b. computing output layer output OlThe formula is shown in formula (X):
in the formula (X), f () is the output layer activation function, wklConnection weight for hidden layer and output layer, HkFor hidden layer output, blIs a threshold value; j is the number of neurons in the hidden layer, and k is the number of neurons in the output layer;
c. calculating error ElFormula (XI):
in the formula (XI), ElIs the mean square error, ylIs the output of the assumed BP neural network model, i.e. the expected output, OlIs the actual output of the BP neural network model, i.e. the output of the output layer;
in the step (3), the specific implementation steps of the back propagation stage include:
d. updating the weight vijAnd wkl，vijIncrease (delta)ij，wklIncrease (delta)kl；i＝1,2,...,n，j＝1,2,...,p；
e. Updating the threshold ajAnd bl，ajIncrease (delta)ij，blIncrease (delta)kl；α is a learning rate, and α is 0.05.
8. The electrocardiosignal denoising method based on the BP neural network and the improved EMD method as claimed in any one of claims 1-7, wherein in step S7, the noise dominant component is denoised by wavelet threshold method, the process is as follows:
selecting an IMF component dominated by noise as a wavelet basis function, and adopting a soft threshold function as a threshold function, wherein the formula is shown as formula (XII):
in formula (xii), w is a wavelet coefficient, sgn () is a sign function, λ is a threshold, p is 1 to N, q is 1 to N, and N is a signal length of the IMF component.
9. A computer device comprising a memory and a processor, wherein the memory stores a computer program, and wherein the processor when executing the computer program implements the steps of the method for denoising electrocardiographic signals based on the BP neural network and the improved EMD method according to any one of claims 1 to 8.
10. A computer-readable storage medium, on which a computer program is stored, wherein the computer program, when being executed by a processor, implements the steps of the method for denoising electrocardiographic signals based on the BP neural network and the improved EMD method according to any one of claims 1 to 8.
Technical Field
The invention belongs to the technical field of biological signal processing, and particularly relates to an electrocardiosignal processing method, electrocardiosignal processing equipment and a storage medium.
Background
Electrocardiosignals, one of the most important physiological signals of the human body, contain a large amount of information and are widely used in the fields of emotion recognition, biomedical treatment and the like. However, the electrocardiographic signal is a very weak signal, which is inevitably accompanied by the occurrence of noise in the measurement process, so that the originally weak signal is more difficult to measure.
The traditional methods for processing electrocardiosignals comprise an FIR filtering method, a self-adaptive filtering method, a wavelet filtering method and the like, however, the local filtering effect of the methods is poor, and sometimes, due to unreasonable setting, useful signals are inevitably removed.
An Empirical Mode Decomposition (EMD) method is a signal decomposition method proposed by Huang et al, and has been currently applied to various fields. The method can decompose a signal into a limited number of intrinsic mode components (IMF), the frequency is arranged from high to low, components which do not accord with the characteristics of the signal to be researched can be removed according to the characteristics and the requirements of the signal to be researched, other components which accord with the characteristics of the signal are processed, and finally the residual components are superposed and reconstructed, so that the denoised signal is obtained.
The most critical step of the EMD decomposition method is the selection of IMF components. There are many choices of IMF components in the conventional EMD decomposition method, and there is no unified standard at present. The traditional selection method mainly depends on subjective experience of people, according to the frequency of the IMF component, the high-frequency IMF component is directly considered to mainly contain a noise signal, the low-frequency IMF component mainly contains a useful signal, and the high-frequency component is directly abandoned, and the residual component is reconstructed. However, in general, the high-frequency IMF component still contains part of the useful signal, and directly discarding it will result in loss of the useful signal, so that the reconstructed signal is inaccurate. The subsequent improved method is to determine the complexity of the signal by calculating the sample entropy of the IMF component, the larger the sample entropy is, the more complex the signal is, the more noise is contained therein, and although the method for calculating the sample entropy has a clear evaluation index, the index is too single, and the final result may not be accurate. Yet another method is to calculate the correlation coefficient of the IMF component with the original signal, draw a correlation coefficient profile, and determine the class of the signal by observing the discontinuities in the correlation coefficient profile. The method has high accuracy, but the workload is too large and the efficiency is low, so that errors are inevitably generated, the classification of IMF components is mistaken, and finally the reconstructed signal is inaccurate.
Disclosure of Invention
Aiming at the defects of the prior art, the invention provides an improved electrocardiosignal processing method, which is a technology of fusing a BP neural network, an improved Empirical Mode Decomposition (EMD) method and a wavelet threshold method.
According to the method, the problem of residual noise in the reconstructed signal after EMD decomposition is solved by adding white Gaussian noise into the original signal, the EMD decomposition is adopted to complete the decomposition of the original signal, the characteristic parameters and the categories of all IMF components are calculated to train a BP neural network model, the IMF components to be classified are classified through the BP neural network model, the IMF components containing noise are denoised through a wavelet threshold method, and finally the effective IMF components are reconstructed to obtain the denoised electrocardiosignals. For selection of IMF components, the method references a plurality of indexes, has high accuracy, greatly improves the efficiency by neural network classification, and has good final denoising effect.
The invention also provides computer equipment and a storage medium.
Interpretation of terms:
the MIT-BIH database is a database provided by the Massachusetts institute of technology for studying arrhythmia, and is an internationally recognized ECG database which can be used as a standard.
The technical scheme of the invention is as follows:
an electrocardiosignal denoising method based on a BP neural network and an improved EMD method comprises the following steps:
s1, selecting a plurality of electrocardiosignal samples from the MIT-BIH database;
s2, decomposing all electrocardiosignals by using an improved EMD method, decomposing each electrocardiosignal into a plurality of IMF components, marking the IMF components obtained by decomposing the same electrocardiosignal as one group, randomly dividing the finally obtained IMF components into two groups, taking the large group with more number as a training set, and taking the rest IMF components as a test set;
s3, calculating the correlation coefficient of each IMF component and the corresponding electrocardiosignal before decomposition, drawing a correlation coefficient distribution graph, and dividing the IMF components into three categories of noise dominant, signal dominant and invalid components according to the correlation coefficient, wherein the signs are IMF respectivelyNoise(s)，IMFSignal，IMFInvalidation；
S4, calculating five characteristic parameters of the IMF component, including margin factor, kurtosis coefficient, baseline energy ratio, QRS energy ratio and peak-to-average ratio;
s5, constructing and training a BP neural network model;
s6, classifying the IMF components by using the trained BP neural network model, reserving the IMF components leading the signals, directly discarding invalid components, and performing next threshold denoising on the components marked as leading noise;
s7, further denoising the noise-dominant component by adopting a wavelet threshold method, and recording the denoised noise-dominant component as IMFDe-noising；
And S8, adding the denoised noise dominant component and the preserved electrocardiosignal dominant component to reconstruct the electrocardiosignal.
Preferably, in step S2, the EMD improvement method is as follows:
s21: adding noise with amplitude distribution obeying Gaussian distribution and power spectral density distribution obeying uniform distribution into the electrocardiosignal obtained in the step S1;
s22: decomposing the electrocardiosignal subjected to the noise addition in the step S21 to obtain a first IMF component, adding m kinds of noises with different noise intensities in the step S21 according to the principle that the noise intensity is gradually decreased until m times of EMD decomposition are completed, wherein the formula of the EMD decomposition is shown as a formula (I):
in formula (I), d (n) is the original ECG signal, IMFk(n) is the decomposed kth IMF component, akAnd (n) is the remainder of the signal after the k-th decomposition.
Preferably, in step S3, the IMF component classification method includes the following steps:
s31: calculating the correlation coefficient rho of each IMF component and the original signal, namely the corresponding electrocardiosignal before decompositionxyThe calculation formula is shown as formula (II):
in formula (II), x (i) is the signal of the correlation coefficient to be calculated, and y (i) is the original signal;
s32: drawing a correlation coefficient distribution graph; in the correlation coefficient distribution diagram, the abscissa is an IMF component, and the ordinate is a correlation coefficient corresponding to the IMF component;
s33: the correlation coefficient distribution map has n points, which are denoted as 1-n, and the corresponding IMF components are denoted as IMF1～IMFnFinding the first turning point in the correlation coefficient distribution diagram and recording the turning point as P, and recording the corresponding IMF component as IMFP，IMF1To IMFPAll IMF components in the interval are marked as noise-dominant components, and then further denoising is needed; from IMFPTo IMFnAnd finding out IMF components with the correlation coefficient smaller than 0.15 from the IMF components in between, recording the IMF components as invalid components, removing the invalid components, and recording the rest components as signal dominant components.
Preferably, in step S4, the margin factor is a statistic for reflecting the interference tolerance of the signal, and the larger the margin factor is, the larger the proportion of the useful information is, the smaller the interference tolerance is, and t is1To express the margin factor, the calculation formula of the margin factor is shown as formula (iii):
in the formula (III), xpIs the peak value, xiN is the number of sampling points for the ith signal sample.
Preferably, in step S4, the kurtosis coefficient is a physical quantity used to indicate the steepness of the top of the data curve; the crest factor of a clean and intact electrocardiosignal is more than 5, and if myoelectric interference, baseline drift, power frequency interference or Gaussian distribution random noise exists, the crest factor of the electrocardiosignal is lower than 5;
kurtosis coefficient t2The formula (iv) is shown as formula (iv):
in the formula (IV), xiFor the ith signal sample value,is the signal sample average.
Preferably, in step S4, the baseline energy ratio is a characteristic parameter used to reflect the baseline interference degree, the characteristic is defined as a quotient between energy in a frequency range of 1-40Hz and energy in a frequency range of 0-40Hz, the frequency of baseline drift is 0.15-0.3Hz, the baseline energy ratio effectively represents that there is enough baseline interference with a larger influence, and a smaller baseline energy ratio indicates more baseline interference; base line energy ratio t3The formula (c) is shown in formula (v):
in the formula (V), P (f) is a function of the power spectral density of the electrocardiosignal, and f is frequency.
Preferably, in step S4, the QRS energy ratio is a physical quantity reflecting the degree of electromyographic interference; the characteristic is defined as the ratio of the energy of 5-15Hz frequency band to the energy of 5-40Hz frequency band; QRS energy ratio t4The formula (c) is shown in formula (vi):
preferably, in step S4, the peak-to-average ratio is a parameter reflecting the difference between the data peak value and the data mean value; the calculation formula of the peak-to-average ratio is shown as the formula (VII):
preferably, according to the present invention, the BP neural network model includes an input layer, a hidden layer, and an output layer;
the input layer has 5 neurons;
the output layer has 3 neurons; three categories of IMF components are used as output vectors of the BP neural network model, and the three categories of the IMF components comprise noise leading components, signal leading components and invalid components which are respectively marked as 0, 1 and 2;
number of neurons in hidden layerm is the number of neurons in the hidden layer, n is the number of neurons in the input layer, l is the number of neurons in the output layer, and alpha is 1-10.
Further preferably, the number of neurons of the hidden layer is 4.
Preferably, in the BP neural network model, the activation functions of the hidden layer and the output layer are both selected from S functions, and the formula of S function f (x) is shown in formula (viii):
preferably, in step S5, the training process of the BP neural network model is as follows:
(1) normalization treatment: normalizing the five characteristic parameters obtained in the step S4;
(2) inputting: inputting the characteristic parameters after normalization processing into a BP neural network model as input vectors;
(3) and a forward propagation stage: the neurons of the training set transmitted into the input layer interact with the weight threshold value and then are transmitted to the neurons of the hidden layer, the neurons are transmitted into the next layer through the calculation of the activation function, and the neurons are gradually transmitted to the output layer by layer in the same way;
(4) and (3) a back propagation stage: and comparing the result of the output layer with the expected result to obtain the error of the output result, transmitting the error layer by layer in a reverse direction to obtain the deviation of the threshold and the weight of each layer, and adjusting the threshold and the weight through the deviation.
(5) Classifying IMF components of a training set by using a trained BP neural network model, evaluating the BP neural network model, wherein the classification accuracy of the BP neural network is 95%, and if the classification accuracy of the BP neural network is not 95%, continuing training, specifically: and calculating characteristic parameters of the IMF components in the training set, classifying the characteristic parameters as the input of the BP neural network model, and calculating the accuracy of the BP neural network model by comparing the classification result of the BP neural network model with the class label of the IMF components.
Preferably, in step (3), the forward propagation stage includes:
a. computing hidden layer output HjThe formula is shown as formula (IX):
in the formula (IX), HjFor the output of the hidden layer, f () is the activation function of the hidden layer, vijIs the connection weight, x, of the input layer and the hidden layeriFor input layer input, ajIs a threshold value; i is the number of neurons in the input layer, and j is the number of neurons in the hidden layer;
b. computing output layer output OlThe formula is shown in formula (X):
in the formula (X), f () is the output layer activation function, wklConnection weight for hidden layer and output layer, HkFor hidden layer output, blIs a threshold value; j is the number of neurons in the hidden layer, and k is the number of neurons in the output layer;
c. calculating error ElFormula (XI):
in the formula (XI), ElIs the mean square error, ylIs a presumed BPOutput of the neural network model, i.e. expected output, OlIs the actual output of the BP neural network model, namely the output layer output.
Preferably, in step (3), the implementation steps of the back propagation stage include:
d. updating the weight vijAnd wkl，vijIncrease (delta)ij，wklIncrease (delta)kl；i＝1,2,...,n，j＝1,2,...,p；
e. Updating the threshold ajAnd bl，ajIncrease (delta)ij，blIncrease (delta)kl；α is a learning rate, and α is 0.05.
Preferably, in step S7, denoising the noise-dominant component by wavelet thresholding is performed, which includes the following steps:
selecting an IMF component dominated by noise as a wavelet basis function, and adopting a soft threshold function as a threshold function, wherein the formula is shown as formula (XII):
in formula (xii), w is a wavelet coefficient, sgn () is a sign function, λ is a threshold, p is 1 to N, q is 1 to N, and N is a signal length of the IMF component.
A computer device comprises a memory and a processor, wherein the memory stores a computer program, and the processor realizes the steps of the electrocardiosignal denoising method based on a BP neural network and an improved EMD method when executing the computer program.
A computer-readable storage medium, on which a computer program is stored which, when being executed by a processor, carries out the steps of a method for denoising electrocardiographic signals based on a BP neural network and an improved EMD method.
The invention has the beneficial effects that:
the method solves the problem that the reconstructed signal in the EMD decomposition method contains residual noise by adding Gaussian white noise to the original signal, evaluates the proportion of the signal and the noise in the IMF component through a BP neural network model and classifies the signal and the noise, overcomes the defect of single index of the traditional IMF component classification method, can more accurately distinguish which components need to be further denoised, can be directly used for signal reconstruction, and can discard which signals, thereby ensuring that useful signals are not lost while the noise is removed, and further ensuring that the reconstructed signal is more accurate.
Drawings
FIG. 1 is a schematic flow chart of an electrocardiosignal denoising method based on a BP neural network and an improved EMD method according to the present invention;
FIG. 2 is a schematic diagram of waveforms of the electrocardiographic signals before denoising;
FIG. 3(a) shows the IMF component IMF after EMD decomposition1A schematic diagram of the waveform of (a);
FIG. 3(b) shows the IMF component IMF after EMD decomposition2A schematic diagram of the waveform of (a);
FIG. 3(c) shows the IMF component IMF after EMD decomposition3A schematic diagram of the waveform of (a);
FIG. 3(d) shows the IMF component IMF after EMD decomposition4A schematic diagram of the waveform of (a);
FIG. 3(e) shows the IMF component IMF after EMD decomposition5A schematic diagram of the waveform of (a);
FIG. 3(f) shows the IMF component IMF after EMD decomposition6A schematic diagram of the waveform of (a);
FIG. 3(g) shows IMF components IMF after EMD decomposition7A schematic diagram of the waveform of (a);
FIG. 3(h) shows the IMF component IMF after EMD decomposition8A schematic diagram of the waveform of (a);
FIG. 3(i) shows the IMF component IMF after EMD decomposition9A schematic diagram of the waveform of (a);
FIG. 3(j) is a waveform diagram of the IMF component res after EMD decomposition;
FIG. 4 is a schematic diagram of a BP neural network model structure;
fig. 5 is a schematic diagram of a reconstructed cardiac signal waveform.
Fig. 6 is a correlation coefficient distribution diagram.
Detailed Description
The invention is further defined in the following, but not limited to, the figures and examples in the description.
Example 1
An electrocardiosignal denoising method based on a BP neural network and an improved EMD method is shown in figure 1 and comprises the following steps:
s1, selecting 30 electrocardiosignal samples of 201 and 230 from the MIT-BIH database;
and S2, decomposing all the electrocardiosignals by using an improved EMD method, wherein each electrocardiosignal is decomposed into 9 IMF components and a remainder, and the IMF components obtained by decomposing the same electrocardiosignal are marked as a group for subsequent signal reconstruction. Randomly dividing the finally obtained IMF components into two groups, wherein the large group with a large number is used as a training set, and the rest IMF components are used as a test set; the training set accounts for 80% and the testing set accounts for 20%. The method comprises the following steps:
(1) taking the electrocardiograph signal 201 as an example, let the original signal be d (t), and calculate the first derivative to obtain all local maximum and minimum points of d (t).
(2) Fitting the envelope function by cubic spline interpolation using all local extreme points, and recording the upper envelope curve of the maximum fitting as amax(t) the lower envelope curve fitted to the minima is denoted amin(t)。
(3) And averaging the two envelope curves to obtain an average envelope curve e (t).
(4) The intermediate signal is m (t), and the original signal is used for subtracting the average envelope curve to obtain the intermediate signal, namely
m(t)＝d(t)-e(t)
If m (t) satisfies the condition that there are extreme points and zero-crossing points differing by at most one in number and the upper and lower envelopes are locally symmetrical about the time axis, then m (t) can be regarded as an IMF component, denoted as IMF1(t) of (d). If not, the steps (1) to (4) need to be repeated with m (t) as a new input signal until the condition is satisfied.
(5) Let remainder a (t), a (t) ═ d (t) -imf (t), first remainder a1(t) isd (t) subtracting IMF1(t) to obtain a new input signal, and repeating the steps (1) - (4) to obtain a second IMF component denoted as IMF2(t), and so on, completing EMD decomposition and finally obtaining IMF1(t)......IMF9(t) and remainder a9(t), EMD decomposition is completed. FIG. 3(a) shows the IMF component IMF after EMD decomposition1(i.e., IMF)1(t)) a waveform schematic diagram; FIG. 3(b) shows the IMF component IMF after EMD decomposition2(i.e., IMF)2(t)) a waveform schematic diagram; FIG. 3(c) shows the IMF component IMF after EMD decomposition3(i.e., IMF)3(t)) a waveform schematic diagram; FIG. 3(d) shows the IMF component IMF after EMD decomposition4(i.e., IMF)4(t)) a waveform schematic diagram; FIG. 3(e) shows the IMF component IMF after EMD decomposition5(i.e., IMF)5(t)) a waveform schematic diagram; FIG. 3(f) shows the IMF component IMF after EMD decomposition6(i.e., IMF)6(t)) a waveform schematic diagram; FIG. 3(g) shows IMF components IMF after EMD decomposition7(i.e., IMF)7(t)) a waveform schematic diagram; FIG. 3(h) shows the IMF component IMF after EMD decomposition8(i.e., IMF)8(t)) a waveform schematic diagram; FIG. 3(i) shows the IMF component IMF after EMD decomposition9(i.e., IMF)9(t)) a waveform schematic diagram; fig. 3(j) is a waveform diagram of the IMF component res (i.e., a9(t)) after EMD decomposition.
And IMF components and the remainder obtained by decomposing the same electrocardiosignal are marked, so that the subsequent electrocardiosignal reconstruction is facilitated.
S3, calculating the correlation coefficient of each IMF component and the corresponding electrocardiosignal before decomposition, drawing a correlation coefficient distribution graph, and dividing the IMF components into three categories of noise dominant, signal dominant and invalid components according to the correlation coefficient as shown in figure 6; the method comprises the following steps:
(1) taking the IMF component obtained by decomposing the electrocardiographic signal 201 as an example, the correlation coefficient between each IMF component and the original electrocardiographic signal is calculated, and the formula is as follows:
in the formula, IMF (t) represents an IMF component, and d (t) represents the electrocardiographic signal 201 for decomposition.
(2) Drawing the correlation coefficient distribution graph of IMF component and original signal, finding the first correlation coefficient turning point, and recording as IMFP，IMF1～IMFPDenoted as the noise-dominant component.
(3) From IMFPTo IMFnAnd finding out IMF components with the correlation coefficient smaller than 0.15 from the IMF components in between, recording the IMF components as invalid components, removing the invalid components, and recording the rest components as signal dominant components.
S4, calculating five characteristic parameters of the IMF component, including margin factor, kurtosis coefficient, baseline energy ratio, QRS energy ratio and peak-to-average ratio; taking the electrocardiosignal 201 as an example, the electrocardiosignal is decomposed into 9 IMF components by an EMD method, and then characteristic parameters of each IMF component, including a margin factor, a kurtosis coefficient, a baseline energy ratio, a QRS energy ratio and a peak-to-average ratio, are calculated by using a listed formula.
The margin factor is used for reflecting the statistic of the interference bearing capacity of the signal, the larger the margin factor is, the larger the proportion of the useful information is, the smaller the interference degree is, and the larger t is1To express the margin factor, the calculation formula of the margin factor is shown as formula (iii):
in the formula (III), xpIs the peak value, xiN is the number of sampling points for the ith signal sample.
The kurtosis coefficient is a physical quantity used for representing the steepness degree of the top of a data curve; the crest factor of a clean and intact electrocardiosignal is more than 5, and if myoelectric interference, baseline drift, power frequency interference or Gaussian distribution random noise exists, the crest factor of the electrocardiosignal is lower than 5;
kurtosis coefficient t2The formula (iv) is shown as formula (iv):
in the formula (IV), xiFor the ith signal sample value,is the signal sample average.
The baseline energy ratio is a characteristic parameter used for reflecting the baseline interference degree, the characteristic is defined as the quotient of the energy of a 1-40Hz frequency band and the energy of a 0-40Hz frequency band, the frequency of baseline drift is 0.15-0.3Hz, the baseline interference with larger influence can exist through the effective representation of the baseline energy ratio, and the smaller the baseline energy ratio is, the more the baseline interference is; base line energy ratio t3The formula (c) is shown in formula (v):
in the formula (V), P (f) is a function of the power spectral density of the electrocardiosignal, and f is frequency.
The QRS energy ratio is a physical quantity used for reflecting the electromyographic interference degree; the characteristic is defined as the ratio of the energy of 5-15Hz frequency band to the energy of 5-40Hz frequency band; the QRS wave energy is mainly concentrated in a frequency bandwidth of 10Hz, the center width is 10Hz, when myoelectric interference occurs, high-frequency components in signals can be increased, the energy ratio can be reduced, whether large myoelectric interference exists can be effectively represented through the QRS energy ratio, and the smaller the value is, the more myoelectric interference exists. QRS energy ratio t4The formula (c) is shown in formula (vi):
the peak-to-average ratio is a parameter reflecting the degree of difference between the data peak value and the data mean value; because Q, R, S, T waves of the electrocardiosignals have different peak values and the peak values of the electrocardiosignals are very high, the difference degree between the peak power and the average power can be effectively reflected by calculating the peak-to-average ratio, and the higher the difference degree is, the higher the peak value is. The calculation formula of the peak-to-average ratio is shown as the formula (VII):
s5, constructing and training a BP neural network model;
s6, classifying the IMF components by using the trained BP neural network model, reserving the IMF components leading the signals, directly discarding invalid components, and performing next threshold denoising on the components marked as leading noise;
s7, further denoising the noise-dominant component by adopting a wavelet threshold method;
and S8, adding the denoised noise dominant component and the preserved electrocardiosignal dominant component to reconstruct the electrocardiosignal.
Example 2
The electrocardiosignal denoising method according to the embodiment 1 is characterized in that:
as shown in fig. 4, the BP neural network model includes an input layer, a hidden layer, and an output layer; the input layer has 5 neurons; the output layer has 3 neurons; three categories of IMF components are used as output vectors of the BP neural network model, and the three categories of the IMF components comprise noise leading components, signal leading components and invalid components which are respectively marked as 0, 1 and 2;
the number of neurons in the hidden layer lacks a theoretical basis and is usually calculated by adopting an empirical formulam is the number of neurons in the hidden layer, n is the number of neurons in the input layer, l is the number of neurons in the output layer, and alpha is 1-10. The number of neurons of the hidden layer is 4.
In the BP neural network model, the activation functions of the hidden layer and the output layer are selected from an S function, and the formula of the S function f (x) is shown as the formula (VIII):
in step S5, the training process of the BP neural network model is as follows:
(1) normalization treatment: normalizing the five characteristic parameters obtained in the step S4;
(2) inputting: inputting the characteristic parameters after normalization processing into a BP neural network model as input vectors;
(3) and a forward propagation stage: the neurons of the training set transmitted into the input layer interact with the weight threshold value and then are transmitted to the neurons of the hidden layer, the neurons are transmitted into the next layer through the calculation of the activation function, and the neurons are gradually transmitted to the output layer by layer in the same way;
(4) and (3) a back propagation stage: and comparing the result of the output layer with the expected result to obtain the error of the output result, transmitting the error layer by layer in a reverse direction to obtain the deviation of the threshold and the weight of each layer, and adjusting the threshold and the weight through the deviation. The whole counter-propagation can be regarded as a negative feedback regulation means.
(5) Classifying IMF components of a training set by using a trained BP neural network model, evaluating the BP neural network model, wherein the classification accuracy of the BP neural network is 95%, and if the classification accuracy of the BP neural network is not 95%, continuing training, specifically: and calculating characteristic parameters of the IMF components in the training set, classifying the characteristic parameters as the input of the BP neural network model, and calculating the accuracy of the BP neural network model by comparing the classification result of the BP neural network model with the class label of the IMF components.
In step (3), the specific implementation steps of the forward propagation stage include:
a. computing hidden layer output HjThe formula is shown as formula (IX):
in the formula (IX), HjFor the output of the hidden layer, f () is the activation function of the hidden layer, vijIs the connection weight, x, of the input layer and the hidden layeriFor input layer input, ajIs a threshold value; i is the number of neurons in the input layer, j is the number of neurons in the hidden layer；
b. Computing output layer output OlThe formula is shown in formula (X):
in the formula (X), f () is the output layer activation function, wklConnection weight for hidden layer and output layer, HkFor hidden layer output, blIs a threshold value; j is the number of neurons in the hidden layer, and k is the number of neurons in the output layer;
c. calculating error ElFormula (XI):
in the formula (XI), ElIs the mean square error, ylIs the output of the assumed BP neural network model, i.e. the expected output, OlIs the actual output of the BP neural network model, namely the output layer output.
In the step (3), the weight values between the input layer and the hidden layer, and between the hidden layer and the output layer are updated by errors in the back propagation stage, which is essentially a negative feedback process, and the specific implementation steps of the back propagation stage include:
d. updating the weight vijAnd wkl，vijIncrease (delta)ij，wklIncrease (delta)kl；i＝1,2,...,n，j＝1,2,...,p；
e. Updating the threshold ajAnd bl，ajIncrease (delta)ij，blIncrease (delta)kl；α is a learning rate, and α is 0.05.
And calculating characteristic parameters of the IMF components in the training set, classifying the characteristic parameters as the input of the BP neural network, and calculating the accuracy of the BP neural network model by comparing the classification result of the BP neural network with the class label of the IMF components. In this example, the effect of neural network model classification is evaluated by calculating the ratio of the predicted correct samples to all the predicted samples, i.e. the accuracy (Acc), and the formula is as follows:
the method comprises the steps of predicting the number of components with signal dominance actually being signal dominance, predicting the number of components with signal dominance actually being invalid components, predicting the number of components with signal dominance actually being noise dominance, predicting the number of components with noise dominance actually being invalid components, predicting the number of components with noise dominance actually being signal dominance, g the number of components with invalid components predicted as being signal dominance actually being signal dominance, h the number of components with invalid components predicted as being invalid components actually being invalid components, and i the number of components with invalid components predicted as being noise dominance actually being noise dominance. The accuracy of the neural network model finally obtained in the embodiment is 95.5%, and the requirement of denoising is met.
Taking 109 electrocardiosignals as signals to be processed, decomposing 109 electrocardiosignals by using an EMD method to obtain 10 IMF components, calculating five characteristic parameters of each component, classifying the IMF components by using a trained neural network model, reserving the IMF components with the categories as the signal leading, discarding the IMF components with the categories as invalid components, and further performing wavelet threshold denoising on the IMF components with the categories as noise leading.
Selecting IMF components with noise dominance as wavelet basis functions, and denoising each noise dominance component by adopting a soft threshold function as a threshold function, wherein the formula is as follows:
wherein IMFq(p) is the p-th sample point value of the q-th noise dominant component, sgn () is the sign function, IMFq(p)' is IMFq(p) denoisingValue of after, thetaqThe threshold for the qth noise dominant component is p 1 to N, q 1 to N, and N is the signal length of the IMF component.
Threshold value thetaqSelecting a fixed threshold which is most widely applied, wherein the formula is as follows:
wherein σqIs the standard deviation of the q-th noise dominant component.
The IMF component with the noise-removed class as the noise dominance, the IMF component with the previously reserved class as the signal dominance and the remainder a9(t) carrying out signal reconstruction to obtain the final electrocardiosignal without noise, wherein the formula is as follows:
d(t)＝IMFsignal(t)+IMFDe-noising(t)+a9(t)
Fig. 2 is a schematic diagram of an electrocardiographic signal waveform before denoising, and fig. 5 is a schematic diagram of an electrocardiographic signal waveform after reconstruction.
In the embodiment, firstly, an EMD decomposition method is adopted to decompose a large number of electrocardiosignals in a sample to obtain a plurality of IMF components and calculate five characteristic parameters of the IMF components, the IMF components are classified according to traditional evaluation indexes and are labeled by categories, a neural network model is trained by taking the five characteristic parameters and the category labels as sample data, then the electrocardiosignals to be denoised are decomposed to obtain the IMF components, the trained neural network is used to classify the IMF components of the signals, the IMF components with the categories as the leading noise are denoised by wavelet threshold, the IMF components with the categories as the leading noise are retained, and finally the denoised IMF components with the leading noise, the IMF components with the leading signal and the remainder are reconstructed to obtain the electrocardiosignals after the noise is removed. The selection of IMF components is more accurate, the noise in the electrocardiosignal can be effectively removed, and the useful information in the signal is reserved.
Example 3
A computer device, comprising a memory and a processor, wherein the memory stores a computer program, and the processor implements the steps of the electrocardiosignal denoising method based on the BP neural network and the improved EMD method in embodiment 1 or 2 when executing the computer program.
Example 4
A computer-readable storage medium, on which a computer program is stored, which, when being executed by a processor, implements the steps of the electrocardiosignal denoising method of embodiment 1 or 2 based on the BP neural network and the improved EMD method.