# Hybrid Prediction Method for ECG Signals Based on VMD, PSR, and RBF Neural Network

Authors: Fuying Huang; Tuanfa Qin; Limei Wang; Haibin Wan
Journal: BioMed Research International
 (2021)

---

## Abstract

To explore a method to predict ECG signals in body area networks (BANs), we propose a hybrid prediction method for ECG signals in this paper. The proposed method combines variational mode decomposition (VMD), phase space reconstruction (PSR), and a radial basis function (RBF) neural network to predict an ECG signal. To reduce the nonstationarity and randomness of the ECG signal, we use VMD to decompose the ECG signal into several intrinsic mode functions (IMFs) with finite bandwidth, which is helpful to improve the prediction accuracy. The input parameters of the RBF neural network affect the prediction accuracy and computational burden. We employ PSR to optimize input parameters of the RBF neural network. To evaluate the prediction performance of the proposed method, we carry out many simulation experiments on ECG data from the MIT-BIH Arrhythmia Database. The experimental results show that the root mean square error (RMSE) and mean absolute error (MAE) of the proposed method are of 10-3 magnitude, while the RMSE and MAE of some competitive prediction methods are of 10-2 magnitude. Compared with other several prediction methods, our method obviously improves the prediction accuracy of ECG signals.


---

## Body

## 1. Introduction

ECG signals are very important for doctors to diagnose diverse kinds of heart diseases. It is of significance to predict ECG signals accurately. Accurate prediction of ECG signals can help doctors know the patient’s condition in advance, while it can also reduce the energy consumption of sensors in body area networks (BANs). In BANs, some sensors are placed under the skin or inside the body, and their batteries are inconvenient to replace. How to reduce the energy consumption and prolong the lifetime of such sensors is a challenge. If a prediction model is established in both a sensor node and the sink node, when the prediction error exceeds the specified threshold value, the sensor node will send the measured data; otherwise, it will not send the measured data, and the sink node will use the predicted data as the measured data [1]. In other words, prediction can reduce the volume of transmitted data [2], thus reducing the energy consumption of the sensor. In BANs, there are many physiological signals, such as the ECG, body temperature, and blood pressure. If these physiological signals can be accurately predicted from the past and current data, the amount of data transmission and the energy consumption of the sensor will be greatly reduced.There are many methods of predicting physiological signals. Among them, the artificial neural network (ANN) has become a widely used prediction method because of its nonlinearity, self-adaption, and fast calculation speed. As a type of ANN, a radial basis function (RBF) neural network has the advantages of fast training speed and ability to not easily fall into the local optimum; thus, RBF neural networks have attracted increasingly more attention in time series prediction [3–6]. The key problem of designing a RBF neural network is to determine the network structure parameters, including the center of the hidden layer, base width, and weight of the output layer. In recent years, several algorithms have been proposed to solve the problem of parameter optimization [7–10]. To solve global optimization problems, an enhanced MSIQDE (improved quantum-inspired differential evolution) algorithm based on mixing multiple strategies was proposed in [7]. This algorithm has better optimization ability than the DE (differential evolution), QDE (quantum-inspired differential evolution), QGA (quantum genetic algorithm), and MSIQDE algorithms. In [8], a novel improved DE algorithm with the wavelet basis function was proposed. This algorithm shows better optimization ability than other several DE algorithms and effectively solved a real-world airport gate assignment problem. To optimize the parameters of photovoltaic models and enhance the conversion efficiency of solar energy, Song et al. [9] proposed a new multipopulation parallel coevolutionary DE algorithm that showed higher accuracy and reliability than other several methods in extracting parameters of solar modules.When using a RBF neural network to predict time series signals, it is not enough to only optimize its internal structure parameters. One must also optimize the input parameters of its input layer. The input parameters of a RBF neural network include input dimension and delay time. Different input parameters produce different prediction results. It is difficult to determine the optimal input parameters. Phase space reconstruction (PSR) [11] is one of the effective methods used to solve the input parameter optimization problem of RBF neural networks. PSR reconstructs the multidimensional phase space structure of the original system through a one-dimensional time series signal. The embedding dimension and delay time of PSR are taken as the input dimension and delay time of the neural network, respectively. Several researchers have used PSR to optimize the input parameters of neural networks and obtained good prediction results [12–14].This paper focuses on the prediction of ECG signals. In BANs, the ECG signal has the largest amount of data among all physiological signals. The sensor will consume a significant amount energy to transmit vast quantities of raw ECG data. To save the energy of the sensor, it is necessary to reduce the transmission of ECG data without affecting the normal use of ECG signals. An ECG signal is a nonlinear nonstationary time series signal and has an inherent random feature. Although an ECG signal is considered a periodic signal, its period does not keep a fixed value. These features hinder the accurate prediction of ECG signals. Several papers have discussed the prediction of ECG signals. Wei et al. [15] developed a model for highly accurate prediction of ECGs and EEGs by combining CNN and bidirectional recurrent neural network (BRNN). In [16], Sun et al. proposed a prediction method of an ECG signal using an error backpropagation neural network (BPNN) and VMD technology. Sun et al. [17] proposed another prediction method of an ECG signal. The method was based on PSR theory and BPNN, with accuracy close to that of the previously mentioned method. To reduce the amount of data transmission, a waveform prediction lightweight algorithm was proposed in [1]. To improve the prediction accuracy, the algorithm used wavelet transform (WT) to preprocess the noise. Liu [18] proposed a data fusion algorithm based on WT and a least-squares support vector machine (LS-SVM). The algorithm used the LS-SVM model to predict an ECG signal. In [19], Heurtefeux et al. investigated the tradeoff between accuracy and complexity to predict ECG values in BANs. They suggested using an ARMA model to predict ECG values, but they could not build an ARMA model to predict ECG signals experimentally. Sun et al. [20] used an ARIMA model to analyze ECG data streams, but they did not use one to predict an ECG signal. An ECG signal prediction method based on an ARIMA model and a discrete wavelet transform (DWT) was proposed in [21]. This method obtained a good prediction effect, but it needed a significant amount of calculation because it used many high-order autoregressive (AR) models. In [22], an ECG signal prediction method was proposed in which the PSR theory and TS fuzzy model were used to predict an ECG signal. The accuracy of the prediction is close to that of [16, 17].To improve the prediction accuracy of ECG signals, we propose a novel hybrid prediction method for ECG signals using variational mode decomposition (VMD), PSR, and a RBF neural network. In this paper, we make the following primary contributions.(a)
We use VMD to decompose a nonstationary ECG signal into several relatively stationary intrinsic mode functions (IMFs), which is helpful to improve the prediction accuracy. The parameterK represents the number of IMF decomposed by VMD. We determine K by the ratio of the residual energy to original signal energy.(b)
We use PSR to optimize the input parameters of the RBF neural network. We determine the embedding dimension by false nearest neighbors (FNN) method. The embedding dimension is the input dimension of the RBF neural network. We determine the delay time by comparing the mutual information (MI) method with the experimental method.(c)
We propose a novel prediction method for ECG signals based on VMD, PSR, and a RBF neural network. Compared with other prediction methods, the proposed method improves the prediction accuracy of ECG signals.The rest of this paper is organized as follows. In Section2, we describe the theories of VMD, PSR, and the RBF neural network and propose a novel hybrid prediction method for ECG signals. The simulation experiments and the analysis of its results are presented in Section 3. In Section 4, we provide concluding remarks.


## 2. Methods

### 2.1. Variational Mode Decomposition

Variational mode decomposition (VMD) can decompose a nonstationary signal into several discrete band-limited IMFs around the center frequency, meeting the condition that the sum of the estimated bandwidths of each mode is the smallest [23]. The decomposition steps can be summarized as follows. The parameter K represents the number of modal components decomposed by VMD and it must be predetermined before VMD.Step 1.
Calculate the analytic signal of each modal functionukt by Hilbert transform
(1)δt+jπt∗ukt.Step 2.
Multiply the analytical signal by the estimated center frequencye−jwkt, and move it to the base frequency spectrum
(2)δt+jπt∗ukte−jwkt.Step 3.
Estimate the bandwidth of each mode by Gaussian smoothing of the demodulated signal, i.e., theL2 norm of the gradient. The constrained variational model is
(3)minukwk∑kδt+jπt∗ukte−jwkt22s.t.∑kuk=x,where x is the input signal and ·2 is the Euclidian distance.Step 4.
Turn the constrained variational model into an unconstrained variational model by introducing the quadratic penalty factorα and Lagrange multiplication operator λt. The extended Lagrange expression is
(4)Luk,wk,λt=α∑kδt+jπt∗ukte−jwkt22+x−∑kuk22+<λt,x−∑kuk>.
Find the saddle point of the extended Lagrange expression using the alternating direction multiplier method (optimization algorithm) to solve the minimization problem of (3). The saddle point is the optimal solution.Step 5.
Obtain the saddle point of the extended Lagrange expression by alternately updatingukn+1, wkn+1, and λkn+1. The equations for this are
(5)u⌣kn+1w=x⌣w−∑i<ku⌣in+1w−∑i>ku⌣inw+λ⌣w/21+2αw−wkn2,(6)wkn+1=∫0∞wu⌣kn+1w2dw∫0∞u⌣kn+1w2dw,(7)λ⌣n+1u⌣kn+1w=λ⌣nu⌣kn+1w+τx⌣u⌣kn+1w−∑ku⌣kn+1u⌣kn+1w,where u⌣w, λ⌣w, and x⌣w are the Fourier transforms of the signals ut, λt, andxt, respectively.Step 6.
Repeat Step 5 until the convergence condition is reached(8)∑ku⌣kn+1−u⌣kn22/u⌣kn22<ε.


### 2.2. Phase Space Reconstruction

Dutch mathematician Floris Takens [24] proved that as long as the embedding dimension is large enough, the reconstructed attractor retains the topological properties of the original attractor. In other words, in the trajectory of the reconstructed space, the reconstructed phase space and differential homeomorphism of the original power system are maintained, and the phase space of the delay time reconstruction keeps the geometric structure of the original system, along with its dynamic characteristics.For a time seriesx1,x2,⋯,xN, if the delay time τ and embedding dimension m can be properly determined, then the phase space can be reconstructed as
(9)X=X1X2⋮XM=x1x1+τ⋯x1+m−1τx2x2+τ⋯x2+m−1τ⋯⋯⋯⋯xMxM+τ⋯xN, where M=N−m−1τ is the length of the vector sequence and Xii=1,2,⋯,M is a point of the phase space.WhenX is used as the input vector , the one-step prediction output vector Y shows as follows:
(10)Y=Y1Y2⋮YM=x2+m−1τx3+m−1τ⋮xN+1.The key to PSR is to correctly determine the embedding dimensionm and delay time τ. The delay time can be determined by the average displacement, autocorrelation function, complex correlation, and mutual information (MI) methods. Methods to determine the embedding dimension include false nearest neighbors (FNN) method, Cao method, and G-P algorithm. Common methods to determine the embedding dimension and delay time at the same time include the window embedding and C-C methods.


#### 2.2.1. Mutual Information Method

Given that there are two systemsS=s1,s2,⋯,sn1 and Q=q1,q2,⋯,qn2, the information entropy of S and Q is
(11)HS=−∑i=1n1Pssilog2Pssi,(12)HQ=−∑j=1n2Pqqjlog2Pqqj,respectively, wherePssi and Pqqj are probabilities of si and qj, respectively.The joint entropy ofS and Q is
(13)HS,Q=−∑i=1n1∑j=1n2Psqsi,qjlog2Psqsi,qj,where Psqsi,qj is the joint probability of si and qj.Given thatS has been measured, the conditional entropy of Q is
(14)HQS=∑iPssiHQsi=−∑i,jPssiPqjsilog2Pqjsi=−∑i,jPsqsi,qjlog2Pqjsi=−∑i,jPsqsi,qjlog2Psqsi,qjPssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPsqsi,qjlog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPssiPqjsilog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑iPssilog2Pssi∑jPqjsi=HS,Q−HS.The MI ofQ and S is
(15)IQ,S=∑j∑iPsqqi,sjIqi,sj=∑j∑iPsqqi,sjIqi−Isjqi=∑j∑iPsqqi,sjlog21Pqqi−log21Pqi/sj=∑j∑iPsqqi,sjlog2Pqi/sjPqqi=∑j∑iPsqqi,sjlog2Pqisj−∑j∑iPsqqi,sjlog2Pqqi=−HQS−∑j∑iPqqiPsjqilog2Pqqi=−HQS−∑iPqqilog2Pqqi∑jPsjqi=−HQS+HQ=HQ−HQS.According to (14) and (15), the following results are obtained.
(16)IQ,S=HQ−HQS=HQ−HS,Q−HS=HQ+HS−HS,Q=IS,Q.Given thatQ is the delayed sequence of S, that is, Q=s1+τ,s2+τ,⋯,sn1+τ, the MI of S and Q is
(17)Iτ=IS,Q=IQ,S.Theτ corresponding to the first local minimum of Iτ is the optimal delay time τ.


#### 2.2.2. False Nearest Neighbors Method

False nearest neighbors (FNN) are phase points that are adjacent in low-dimensional space, but not adjacent after mapping to a certain high-dimensional space. SupposeXmi=xi,xi+τ,⋯,xi+m−1τ is a phase vector in m-dimensional reconstruction space, and the nearest neighbor point of Xmi is XmNNi. If XmNNi satisfies the following inequality when the dimension of reconstruction space increases from m dimensions to m+1 dimensions, XmNNi is the FNN of Xmi.
(18)Xm+1i−Xm+1NNi2−Xmi−XmNNi2Xmi−XmNNi2≥RT,where ·2 is the Euclidian distance and RT is the threshold value.As the embedding dimensionm increases from small to large, the proportion of FNN points is calculated. When the proportion of FNN points is very small or the number of FNNs no longer changes as m increases, then m is the best embedding dimension.


### 2.3. RBF Neural Network

Unlike a BP neural network, which is a global approximation network, a RBF neural network is a kind of local approximation network. As long as there are enough hidden neurons, a RBF neural network can approximate any continuous nonlinear function with any accuracy. Compared with a BP neural network, a RBF neural network has the advantage of fast training speed, and it does not easily fall into local minima.The basic idea of a RBF neural network is that the RBF of the hidden layer node transforms the input vector, maps the low-dimensional input data to the high-dimensional space, weights the sums of the output of the node, and maps the results from the high-dimensional space to the low-dimensional space for output [25].The network structure of a RBF neural network generally consists of an input layer, a hidden layer, and an output layer, as shown in Figure1.Figure 1
Structure of RBF neural network.In Figure1, a Gaussian function is usually used as the RBF of a hidden layer as follows.
(19)hj=exp−X−Cj22σj2,j=1,2,⋯,m,where X=x1,x2,⋯,xn, Cj is the center of the Gaussian function, σj is the variance of the Gaussian function, and m is the number of hidden layer nodes.The output of the RBF neural network is(20)yk=∑j=1mwjkhj=∑j=1mwjkexp−X−Cj22σj2,k=1,2,⋯,r,where wjk is the connection weight from the hidden layer to the output layer.The learning algorithm of the RBF neural network solves for three parameters, which are the center of the RBF, the variance of the RBF, and the weight of the connection between the hidden layer node and output layer node. Common learning algorithms of RBF neural networks includek-means, the gradient training method, and the orthogonal least square product (OLS) algorithm.


### 2.4. Proposed ECG Signal Prediction Method

Based on the study of ECG prediction, this paper proposes a hybrid method of ECG signal prediction based on VMD, PSR, and a RBF neural network. Its flowchart is shown in Figure2.Figure 2
Flowchart of the proposed method.The prediction steps of the proposed method are as follows.Step 1. Decompose an ECG signal intoK IMFs by VMD. In the experiment, we take K=10 according to the ratio of residual energy to the original signal energy.Step 2. Determine the embedding dimensionm and delay time τ of each IMF. We determine the embedding dimension by the FNN method and determine the delay time by comparing the MI method with the experimental method.Step 3. Reconstruct the phase space of each IMF according to the embedding dimension and delay time.Step 4. According to the trained set of each IMF, establish a RBF neural network and use it to predict the test set of each IMF. The embedding dimension is the input dimension of the RBF neural network.Step 5. Add the prediction results of the RBF neural network to obtain the final ECG signal prediction result.Step 6. Analyze the prediction error and compare it with some traditional prediction methods.


## 2.1. Variational Mode Decomposition

Variational mode decomposition (VMD) can decompose a nonstationary signal into several discrete band-limited IMFs around the center frequency, meeting the condition that the sum of the estimated bandwidths of each mode is the smallest [23]. The decomposition steps can be summarized as follows. The parameter K represents the number of modal components decomposed by VMD and it must be predetermined before VMD.Step 1.
Calculate the analytic signal of each modal functionukt by Hilbert transform
(1)δt+jπt∗ukt.Step 2.
Multiply the analytical signal by the estimated center frequencye−jwkt, and move it to the base frequency spectrum
(2)δt+jπt∗ukte−jwkt.Step 3.
Estimate the bandwidth of each mode by Gaussian smoothing of the demodulated signal, i.e., theL2 norm of the gradient. The constrained variational model is
(3)minukwk∑kδt+jπt∗ukte−jwkt22s.t.∑kuk=x,where x is the input signal and ·2 is the Euclidian distance.Step 4.
Turn the constrained variational model into an unconstrained variational model by introducing the quadratic penalty factorα and Lagrange multiplication operator λt. The extended Lagrange expression is
(4)Luk,wk,λt=α∑kδt+jπt∗ukte−jwkt22+x−∑kuk22+<λt,x−∑kuk>.
Find the saddle point of the extended Lagrange expression using the alternating direction multiplier method (optimization algorithm) to solve the minimization problem of (3). The saddle point is the optimal solution.Step 5.
Obtain the saddle point of the extended Lagrange expression by alternately updatingukn+1, wkn+1, and λkn+1. The equations for this are
(5)u⌣kn+1w=x⌣w−∑i<ku⌣in+1w−∑i>ku⌣inw+λ⌣w/21+2αw−wkn2,(6)wkn+1=∫0∞wu⌣kn+1w2dw∫0∞u⌣kn+1w2dw,(7)λ⌣n+1u⌣kn+1w=λ⌣nu⌣kn+1w+τx⌣u⌣kn+1w−∑ku⌣kn+1u⌣kn+1w,where u⌣w, λ⌣w, and x⌣w are the Fourier transforms of the signals ut, λt, andxt, respectively.Step 6.
Repeat Step 5 until the convergence condition is reached(8)∑ku⌣kn+1−u⌣kn22/u⌣kn22<ε.


## 2.2. Phase Space Reconstruction

Dutch mathematician Floris Takens [24] proved that as long as the embedding dimension is large enough, the reconstructed attractor retains the topological properties of the original attractor. In other words, in the trajectory of the reconstructed space, the reconstructed phase space and differential homeomorphism of the original power system are maintained, and the phase space of the delay time reconstruction keeps the geometric structure of the original system, along with its dynamic characteristics.For a time seriesx1,x2,⋯,xN, if the delay time τ and embedding dimension m can be properly determined, then the phase space can be reconstructed as
(9)X=X1X2⋮XM=x1x1+τ⋯x1+m−1τx2x2+τ⋯x2+m−1τ⋯⋯⋯⋯xMxM+τ⋯xN, where M=N−m−1τ is the length of the vector sequence and Xii=1,2,⋯,M is a point of the phase space.WhenX is used as the input vector , the one-step prediction output vector Y shows as follows:
(10)Y=Y1Y2⋮YM=x2+m−1τx3+m−1τ⋮xN+1.The key to PSR is to correctly determine the embedding dimensionm and delay time τ. The delay time can be determined by the average displacement, autocorrelation function, complex correlation, and mutual information (MI) methods. Methods to determine the embedding dimension include false nearest neighbors (FNN) method, Cao method, and G-P algorithm. Common methods to determine the embedding dimension and delay time at the same time include the window embedding and C-C methods.


### 2.2.1. Mutual Information Method

Given that there are two systemsS=s1,s2,⋯,sn1 and Q=q1,q2,⋯,qn2, the information entropy of S and Q is
(11)HS=−∑i=1n1Pssilog2Pssi,(12)HQ=−∑j=1n2Pqqjlog2Pqqj,respectively, wherePssi and Pqqj are probabilities of si and qj, respectively.The joint entropy ofS and Q is
(13)HS,Q=−∑i=1n1∑j=1n2Psqsi,qjlog2Psqsi,qj,where Psqsi,qj is the joint probability of si and qj.Given thatS has been measured, the conditional entropy of Q is
(14)HQS=∑iPssiHQsi=−∑i,jPssiPqjsilog2Pqjsi=−∑i,jPsqsi,qjlog2Pqjsi=−∑i,jPsqsi,qjlog2Psqsi,qjPssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPsqsi,qjlog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPssiPqjsilog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑iPssilog2Pssi∑jPqjsi=HS,Q−HS.The MI ofQ and S is
(15)IQ,S=∑j∑iPsqqi,sjIqi,sj=∑j∑iPsqqi,sjIqi−Isjqi=∑j∑iPsqqi,sjlog21Pqqi−log21Pqi/sj=∑j∑iPsqqi,sjlog2Pqi/sjPqqi=∑j∑iPsqqi,sjlog2Pqisj−∑j∑iPsqqi,sjlog2Pqqi=−HQS−∑j∑iPqqiPsjqilog2Pqqi=−HQS−∑iPqqilog2Pqqi∑jPsjqi=−HQS+HQ=HQ−HQS.According to (14) and (15), the following results are obtained.
(16)IQ,S=HQ−HQS=HQ−HS,Q−HS=HQ+HS−HS,Q=IS,Q.Given thatQ is the delayed sequence of S, that is, Q=s1+τ,s2+τ,⋯,sn1+τ, the MI of S and Q is
(17)Iτ=IS,Q=IQ,S.Theτ corresponding to the first local minimum of Iτ is the optimal delay time τ.


### 2.2.2. False Nearest Neighbors Method

False nearest neighbors (FNN) are phase points that are adjacent in low-dimensional space, but not adjacent after mapping to a certain high-dimensional space. SupposeXmi=xi,xi+τ,⋯,xi+m−1τ is a phase vector in m-dimensional reconstruction space, and the nearest neighbor point of Xmi is XmNNi. If XmNNi satisfies the following inequality when the dimension of reconstruction space increases from m dimensions to m+1 dimensions, XmNNi is the FNN of Xmi.
(18)Xm+1i−Xm+1NNi2−Xmi−XmNNi2Xmi−XmNNi2≥RT,where ·2 is the Euclidian distance and RT is the threshold value.As the embedding dimensionm increases from small to large, the proportion of FNN points is calculated. When the proportion of FNN points is very small or the number of FNNs no longer changes as m increases, then m is the best embedding dimension.


## 2.2.1. Mutual Information Method

Given that there are two systemsS=s1,s2,⋯,sn1 and Q=q1,q2,⋯,qn2, the information entropy of S and Q is
(11)HS=−∑i=1n1Pssilog2Pssi,(12)HQ=−∑j=1n2Pqqjlog2Pqqj,respectively, wherePssi and Pqqj are probabilities of si and qj, respectively.The joint entropy ofS and Q is
(13)HS,Q=−∑i=1n1∑j=1n2Psqsi,qjlog2Psqsi,qj,where Psqsi,qj is the joint probability of si and qj.Given thatS has been measured, the conditional entropy of Q is
(14)HQS=∑iPssiHQsi=−∑i,jPssiPqjsilog2Pqjsi=−∑i,jPsqsi,qjlog2Pqjsi=−∑i,jPsqsi,qjlog2Psqsi,qjPssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPsqsi,qjlog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑i,jPssiPqjsilog2Pssi=−∑i,jPsqsi,qjlog2Psqsi,qj−−∑iPssilog2Pssi∑jPqjsi=HS,Q−HS.The MI ofQ and S is
(15)IQ,S=∑j∑iPsqqi,sjIqi,sj=∑j∑iPsqqi,sjIqi−Isjqi=∑j∑iPsqqi,sjlog21Pqqi−log21Pqi/sj=∑j∑iPsqqi,sjlog2Pqi/sjPqqi=∑j∑iPsqqi,sjlog2Pqisj−∑j∑iPsqqi,sjlog2Pqqi=−HQS−∑j∑iPqqiPsjqilog2Pqqi=−HQS−∑iPqqilog2Pqqi∑jPsjqi=−HQS+HQ=HQ−HQS.According to (14) and (15), the following results are obtained.
(16)IQ,S=HQ−HQS=HQ−HS,Q−HS=HQ+HS−HS,Q=IS,Q.Given thatQ is the delayed sequence of S, that is, Q=s1+τ,s2+τ,⋯,sn1+τ, the MI of S and Q is
(17)Iτ=IS,Q=IQ,S.Theτ corresponding to the first local minimum of Iτ is the optimal delay time τ.


## 2.2.2. False Nearest Neighbors Method

False nearest neighbors (FNN) are phase points that are adjacent in low-dimensional space, but not adjacent after mapping to a certain high-dimensional space. SupposeXmi=xi,xi+τ,⋯,xi+m−1τ is a phase vector in m-dimensional reconstruction space, and the nearest neighbor point of Xmi is XmNNi. If XmNNi satisfies the following inequality when the dimension of reconstruction space increases from m dimensions to m+1 dimensions, XmNNi is the FNN of Xmi.
(18)Xm+1i−Xm+1NNi2−Xmi−XmNNi2Xmi−XmNNi2≥RT,where ·2 is the Euclidian distance and RT is the threshold value.As the embedding dimensionm increases from small to large, the proportion of FNN points is calculated. When the proportion of FNN points is very small or the number of FNNs no longer changes as m increases, then m is the best embedding dimension.


## 2.3. RBF Neural Network

Unlike a BP neural network, which is a global approximation network, a RBF neural network is a kind of local approximation network. As long as there are enough hidden neurons, a RBF neural network can approximate any continuous nonlinear function with any accuracy. Compared with a BP neural network, a RBF neural network has the advantage of fast training speed, and it does not easily fall into local minima.The basic idea of a RBF neural network is that the RBF of the hidden layer node transforms the input vector, maps the low-dimensional input data to the high-dimensional space, weights the sums of the output of the node, and maps the results from the high-dimensional space to the low-dimensional space for output [25].The network structure of a RBF neural network generally consists of an input layer, a hidden layer, and an output layer, as shown in Figure1.Figure 1
Structure of RBF neural network.In Figure1, a Gaussian function is usually used as the RBF of a hidden layer as follows.
(19)hj=exp−X−Cj22σj2,j=1,2,⋯,m,where X=x1,x2,⋯,xn, Cj is the center of the Gaussian function, σj is the variance of the Gaussian function, and m is the number of hidden layer nodes.The output of the RBF neural network is(20)yk=∑j=1mwjkhj=∑j=1mwjkexp−X−Cj22σj2,k=1,2,⋯,r,where wjk is the connection weight from the hidden layer to the output layer.The learning algorithm of the RBF neural network solves for three parameters, which are the center of the RBF, the variance of the RBF, and the weight of the connection between the hidden layer node and output layer node. Common learning algorithms of RBF neural networks includek-means, the gradient training method, and the orthogonal least square product (OLS) algorithm.


## 2.4. Proposed ECG Signal Prediction Method

Based on the study of ECG prediction, this paper proposes a hybrid method of ECG signal prediction based on VMD, PSR, and a RBF neural network. Its flowchart is shown in Figure2.Figure 2
Flowchart of the proposed method.The prediction steps of the proposed method are as follows.Step 1. Decompose an ECG signal intoK IMFs by VMD. In the experiment, we take K=10 according to the ratio of residual energy to the original signal energy.Step 2. Determine the embedding dimensionm and delay time τ of each IMF. We determine the embedding dimension by the FNN method and determine the delay time by comparing the MI method with the experimental method.Step 3. Reconstruct the phase space of each IMF according to the embedding dimension and delay time.Step 4. According to the trained set of each IMF, establish a RBF neural network and use it to predict the test set of each IMF. The embedding dimension is the input dimension of the RBF neural network.Step 5. Add the prediction results of the RBF neural network to obtain the final ECG signal prediction result.Step 6. Analyze the prediction error and compare it with some traditional prediction methods.


## 3. Results and Discussion

All ECG data in the simulation experiment are from the MIT-BIH Arrhythmia Database [26]. We randomly selected No. 100 ECG data, which consists of 2768 data points, for the experiment. We used two-thirds of the source data as the trained set (1845 data) and one-third of the source data as the test set (923 data).


### 3.1. Performance Measures

The common performance measures of prediction methods are root mean square error (RMSE), mean absolute error (MAE), and mean square error (MSE), defined as(21)RMSE=1N∑n=1NXn−X⌣n2,(22)MAE=1N∑n=1NXn−X⌣n,(23)MSE=1N∑n=1NXn−X⌣n2,where X⌣n is the predicted value of Xn and N is the length of Xn.


### 3.2. VMD of ECG Signal

A key problem of VMD is to set the numberK of modal components. K affects the accuracy of the final prediction. If K is too small, the original series may not be decomposed sufficiently. If the signal is too large, the difference between each component becomes very small, and it gives rise to unnecessary computing overhead. K can be determined according to the ratio Rres of residual energy to the original signal energy. Rres is defined as follows:
(24)Rres=1N∑n=1NXn−∑k=1KuknXn,where Xn is the original signal, ukn is the IMF, and N is the sample number. The decision rule of K is that when Rres is less than 1% and there is no significant downward trend, the number K can be determined [27]. For the No. 100 ECG data, Rres of VMD with different K are shown in Table 1.Table 1
Rres of VMD with different K.
K1234567891011Rres0.06580.02580.01160.00640.00350.00210.00180.00150.00180.00170.0017Table1 shows that Rres has no obvious downward trend when K=10. Therefore, we set K=10 in the experiment.The No. 100 ECG was decomposed into 10 IMFs by VMD as shown in Figure3. IMF1 is the residual, and IMF2-IMF10 is the component, sorted from low to high frequency.Figure 3
VMD of ECG signal.


### 3.3. Determination of Delay Time and Embedding Dimension

We determined the embedding dimensionm and delay time τ of each IMF. The embedding dimension m of IMF3 was determined by the FNN method, as shown in Figure 4.Figure 4
Determining embedding dimension by FNN method.From Figure4, we obtained m=5. The delay time τ of IMF3 was determined by the MI method, as shown in Figure 5.Figure 5
Determining delay time by the MI method.In Figure5, we obtained the delay time τ=7 according to the first local minimum of MI. The embedding dimension and delay time of each IMF are shown in Table 2.Table 2
The parameters of each IMF.
ParametersIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10τ311875432213m4455555444In the experiment, we obtained better prediction results when the delay time of each IMF was 1 (i.e.,τ=1). The comparison results are shown in Table 3.Table 3
The comparison results of different delay times.
Delay timeIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10RMSEMAEτ3118754322130.00540.0041τ∗11111111110.00350.0027According to the comparison in Table3, we finally took the delay of each IMF as 1 (i.e., τ=1).


### 3.4. Optimization of RBF Neural Network

In the simulation experiment, we used MATLAB functionnewrbe(P,T,spread) to design a RBF neural network. The parameter P is the input vectors of the RBF network, the parameter T is the output vectors of the RBF network, and the parameter spread is the spread of the RBF network. We used PSR to optimize P and T, that is, we used equations (9) and (10) to determine P and T, respectively. When the parameters P and T had been determined, we determined the parameter spread by experiments. Some experimental results are shown in Table 4.Table 4
Some experimental results of the parameterspread.
Spread0.10.20.30.40.50.60.70.80.91RMSE0.11150.00550.00370.00350.00350.00350.00360.00360.00360.0036MAE0.03460.00350.00280.00270.00270.00270.00270.00280.00280.0028We carried out the simulation experiment in the range of [0.1, 100] of the parameterspread. The best prediction results in the simulation experiment are RMSE=0.0035 and MAE=0.0027 when spread=0.4, 0.5, and 0.6.


### 3.5. Prediction Results of the Test Set

We used the proposed method to predict the test set of the No. 100 ECG signal. The prediction waveform is shown in Figure6.Figure 6
Prediction waveform: (a) prediction waveform of test set; (b) a part of prediction waveform of (a).
(a)(b)Figure6(b) is a partial ECG signal amplification from [180 : 210] in Figure 6(a). Figure 6 shows that the predicted waveform fits well with its original ECG waveform.The prediction errors are shown in Figure7.Figure 7
The prediction errors of test set.In Figure7, in addition to the large error at the peak point, the error in other places is relatively small.To evaluate the generality of the proposed method, we used ECG data other than the No.100 ECG data to carry out experiments. The experimental results are shown in Table5.Table 5
Prediction performance of different ECG numbers.
ECG numbersRMSEMSEMAE1000.00351.2323e−050.00271010.00644.1198e−050.00331030.00482.3258e−050.00361050.00266.8964e−060.00201150.00908.0978e−050.00391180.00826.6581e−050.00371230.00411.6513e−050.00302010.00256.0916e−060.00182070.00214.2296e−060.00152120.00431.8733e−050.00332190.00482.2814e−050.00372230.00361.3294e−050.00282340.00381.4448e−050.0030Some prediction waveforms are shown in Figure8.Figure 8
Prediction results of different ECG numbers: (a) No. 105; (b) No. 115; (c) No. 207; (d) No. 219.
(a)(b)(c)(d)From Table5 and Figure 8, the prediction performances of the proposed method on other ECG data are very close to that of the No. 100 ECG data. This illustrates that the generality of the proposed method is good and it is very suitable for ECG signal prediction.


### 3.6. Comparison with Other Prediction Methods

To illustrate the advantages of the proposed prediction method, we compared it with the prediction methods of Sun et al. [16], Sun et al. [17], and Su et al. [22]. We experimented on the same data source, the No. 100 ECG data, and the experimental results are shown in Table 6 and Figure 9.Table 6
Comparison with common prediction methods.
MethodsRMSEMAEThis paper0.00350.0027Sun et al. [16]0.02330.0157Sun et al. [17]0.04230.0240Su et al. [22]0.01460.0106Figure 9
Comparison with common prediction methods.As shown in Table6, the RMSE and MAE of the proposed method are of 10-3 magnitude, while the RMSE and MAE of [16, 17, 22] are of 10-2 magnitude. Table 6 and Figure 9 show that the RMSE and MAE of the proposed method are significantly smaller than those of [16], [17], and [22]. This illustrates that the prediction accuracy of this paper is much higher than that of [16, 17, 22].We also compared the proposed method with some traditional prediction models, such as ARMA, LS-SVM, SVM, and Kalman. The experimental results on the data source of No. 100 ECG are shown in Table7.Table 7
Comparison with traditional prediction models.
Prediction modelsRMSEThis paper0.0035LS-SVM [1]0.1754SVM [18]0.0669ARMA [18]0.4630Kalman [18]0.3967We can see from Table7 that the RMSE of the proposed method is much smaller than that of LS-SVM, SVM, ARMA, and Kalman.In addition, the proposed method (VMD-PSR-RBF for short) was compared with some approximate hybrid prediction methods, such as the method based on WT, PSR, and a RBF neural network (WT-PSR-RBF for short), the method based on EMD, PSR, and a RBF neural network (EMD-PSR-RBF for short), the method based on VMD, PSR, and a BP neural network (VMD-PSR-BP for short), and the method based on VMD, PSR, and a GRNN (VMD-PSR-GRNN for short). The experimental data are the No. 100 ECG data, and the experimental results are shown in Table8 and Figure 10.Table 8
Comparison with other hybrid prediction methods.
Prediction methodsRMSEMSEMAEVMD-PSR-RBF (this paper)0.00351.2323e−050.0027WT-PSR-RBF0.00522.6576e−050.0031EMD-PSR-RBF0.01281.6496e−040.0094VMD-PSR-BP0.01552.3991e−040.0125VMD-PSR-GRNN0.01291.6752e−040.0089Figure 10
Comparison with some approximate hybrid methods.Table8 and Figure 10 show that the prediction accuracy is better than that of several other approximate hybrid prediction methods.The comparison prediction waveforms are shown in Figures11 and 12.Figure 11
Comparison with VMD-PSR-BP and VMD-PSR-GRNN: (a) prediction waveform of test set; (b) a part of prediction waveform of (a).
(a)(b)Figure 12
Comparison with EMD-PSR-RBF and WT-PSR-RBF.
(a)(b)Figure11(b) is a partial ECG signal amplification from [192 : 198] in Figure 11(a). Figure 11(b) shows that the prediction curve of the proposed method (VMD-PSR-RBF) is closer to the original ECG curve than that of VMD-PSR-BP and VMD-PSR-GRNN.As shown in Figure12(b), the prediction curve of the proposed method (VMD-PSR-RBF) is also closer to the original ECG curve than that of EMD-PSR-BP and WT-PSR-RBF.


## 3.1. Performance Measures

The common performance measures of prediction methods are root mean square error (RMSE), mean absolute error (MAE), and mean square error (MSE), defined as(21)RMSE=1N∑n=1NXn−X⌣n2,(22)MAE=1N∑n=1NXn−X⌣n,(23)MSE=1N∑n=1NXn−X⌣n2,where X⌣n is the predicted value of Xn and N is the length of Xn.


## 3.2. VMD of ECG Signal

A key problem of VMD is to set the numberK of modal components. K affects the accuracy of the final prediction. If K is too small, the original series may not be decomposed sufficiently. If the signal is too large, the difference between each component becomes very small, and it gives rise to unnecessary computing overhead. K can be determined according to the ratio Rres of residual energy to the original signal energy. Rres is defined as follows:
(24)Rres=1N∑n=1NXn−∑k=1KuknXn,where Xn is the original signal, ukn is the IMF, and N is the sample number. The decision rule of K is that when Rres is less than 1% and there is no significant downward trend, the number K can be determined [27]. For the No. 100 ECG data, Rres of VMD with different K are shown in Table 1.Table 1
Rres of VMD with different K.
K1234567891011Rres0.06580.02580.01160.00640.00350.00210.00180.00150.00180.00170.0017Table1 shows that Rres has no obvious downward trend when K=10. Therefore, we set K=10 in the experiment.The No. 100 ECG was decomposed into 10 IMFs by VMD as shown in Figure3. IMF1 is the residual, and IMF2-IMF10 is the component, sorted from low to high frequency.Figure 3
VMD of ECG signal.


## 3.3. Determination of Delay Time and Embedding Dimension

We determined the embedding dimensionm and delay time τ of each IMF. The embedding dimension m of IMF3 was determined by the FNN method, as shown in Figure 4.Figure 4
Determining embedding dimension by FNN method.From Figure4, we obtained m=5. The delay time τ of IMF3 was determined by the MI method, as shown in Figure 5.Figure 5
Determining delay time by the MI method.In Figure5, we obtained the delay time τ=7 according to the first local minimum of MI. The embedding dimension and delay time of each IMF are shown in Table 2.Table 2
The parameters of each IMF.
ParametersIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10τ311875432213m4455555444In the experiment, we obtained better prediction results when the delay time of each IMF was 1 (i.e.,τ=1). The comparison results are shown in Table 3.Table 3
The comparison results of different delay times.
Delay timeIMF1IMF2IMF3IMF4IMF5IMF6IMF7IMF8IMF9IMF10RMSEMAEτ3118754322130.00540.0041τ∗11111111110.00350.0027According to the comparison in Table3, we finally took the delay of each IMF as 1 (i.e., τ=1).


## 3.4. Optimization of RBF Neural Network

In the simulation experiment, we used MATLAB functionnewrbe(P,T,spread) to design a RBF neural network. The parameter P is the input vectors of the RBF network, the parameter T is the output vectors of the RBF network, and the parameter spread is the spread of the RBF network. We used PSR to optimize P and T, that is, we used equations (9) and (10) to determine P and T, respectively. When the parameters P and T had been determined, we determined the parameter spread by experiments. Some experimental results are shown in Table 4.Table 4
Some experimental results of the parameterspread.
Spread0.10.20.30.40.50.60.70.80.91RMSE0.11150.00550.00370.00350.00350.00350.00360.00360.00360.0036MAE0.03460.00350.00280.00270.00270.00270.00270.00280.00280.0028We carried out the simulation experiment in the range of [0.1, 100] of the parameterspread. The best prediction results in the simulation experiment are RMSE=0.0035 and MAE=0.0027 when spread=0.4, 0.5, and 0.6.


## 3.5. Prediction Results of the Test Set

We used the proposed method to predict the test set of the No. 100 ECG signal. The prediction waveform is shown in Figure6.Figure 6
Prediction waveform: (a) prediction waveform of test set; (b) a part of prediction waveform of (a).
(a)(b)Figure6(b) is a partial ECG signal amplification from [180 : 210] in Figure 6(a). Figure 6 shows that the predicted waveform fits well with its original ECG waveform.The prediction errors are shown in Figure7.Figure 7
The prediction errors of test set.In Figure7, in addition to the large error at the peak point, the error in other places is relatively small.To evaluate the generality of the proposed method, we used ECG data other than the No.100 ECG data to carry out experiments. The experimental results are shown in Table5.Table 5
Prediction performance of different ECG numbers.
ECG numbersRMSEMSEMAE1000.00351.2323e−050.00271010.00644.1198e−050.00331030.00482.3258e−050.00361050.00266.8964e−060.00201150.00908.0978e−050.00391180.00826.6581e−050.00371230.00411.6513e−050.00302010.00256.0916e−060.00182070.00214.2296e−060.00152120.00431.8733e−050.00332190.00482.2814e−050.00372230.00361.3294e−050.00282340.00381.4448e−050.0030Some prediction waveforms are shown in Figure8.Figure 8
Prediction results of different ECG numbers: (a) No. 105; (b) No. 115; (c) No. 207; (d) No. 219.
(a)(b)(c)(d)From Table5 and Figure 8, the prediction performances of the proposed method on other ECG data are very close to that of the No. 100 ECG data. This illustrates that the generality of the proposed method is good and it is very suitable for ECG signal prediction.


## 3.6. Comparison with Other Prediction Methods

To illustrate the advantages of the proposed prediction method, we compared it with the prediction methods of Sun et al. [16], Sun et al. [17], and Su et al. [22]. We experimented on the same data source, the No. 100 ECG data, and the experimental results are shown in Table 6 and Figure 9.Table 6
Comparison with common prediction methods.
MethodsRMSEMAEThis paper0.00350.0027Sun et al. [16]0.02330.0157Sun et al. [17]0.04230.0240Su et al. [22]0.01460.0106Figure 9
Comparison with common prediction methods.As shown in Table6, the RMSE and MAE of the proposed method are of 10-3 magnitude, while the RMSE and MAE of [16, 17, 22] are of 10-2 magnitude. Table 6 and Figure 9 show that the RMSE and MAE of the proposed method are significantly smaller than those of [16], [17], and [22]. This illustrates that the prediction accuracy of this paper is much higher than that of [16, 17, 22].We also compared the proposed method with some traditional prediction models, such as ARMA, LS-SVM, SVM, and Kalman. The experimental results on the data source of No. 100 ECG are shown in Table7.Table 7
Comparison with traditional prediction models.
Prediction modelsRMSEThis paper0.0035LS-SVM [1]0.1754SVM [18]0.0669ARMA [18]0.4630Kalman [18]0.3967We can see from Table7 that the RMSE of the proposed method is much smaller than that of LS-SVM, SVM, ARMA, and Kalman.In addition, the proposed method (VMD-PSR-RBF for short) was compared with some approximate hybrid prediction methods, such as the method based on WT, PSR, and a RBF neural network (WT-PSR-RBF for short), the method based on EMD, PSR, and a RBF neural network (EMD-PSR-RBF for short), the method based on VMD, PSR, and a BP neural network (VMD-PSR-BP for short), and the method based on VMD, PSR, and a GRNN (VMD-PSR-GRNN for short). The experimental data are the No. 100 ECG data, and the experimental results are shown in Table8 and Figure 10.Table 8
Comparison with other hybrid prediction methods.
Prediction methodsRMSEMSEMAEVMD-PSR-RBF (this paper)0.00351.2323e−050.0027WT-PSR-RBF0.00522.6576e−050.0031EMD-PSR-RBF0.01281.6496e−040.0094VMD-PSR-BP0.01552.3991e−040.0125VMD-PSR-GRNN0.01291.6752e−040.0089Figure 10
Comparison with some approximate hybrid methods.Table8 and Figure 10 show that the prediction accuracy is better than that of several other approximate hybrid prediction methods.The comparison prediction waveforms are shown in Figures11 and 12.Figure 11
Comparison with VMD-PSR-BP and VMD-PSR-GRNN: (a) prediction waveform of test set; (b) a part of prediction waveform of (a).
(a)(b)Figure 12
Comparison with EMD-PSR-RBF and WT-PSR-RBF.
(a)(b)Figure11(b) is a partial ECG signal amplification from [192 : 198] in Figure 11(a). Figure 11(b) shows that the prediction curve of the proposed method (VMD-PSR-RBF) is closer to the original ECG curve than that of VMD-PSR-BP and VMD-PSR-GRNN.As shown in Figure12(b), the prediction curve of the proposed method (VMD-PSR-RBF) is also closer to the original ECG curve than that of EMD-PSR-BP and WT-PSR-RBF.


## 4. Conclusions

In this paper, we propose a hybrid prediction method for ECG signals using VMD, PSR, and a RBF neural network. Before prediction, we use VMD technology to preprocess the ECG signals in order to remove their nonstationarity and randomness. We determine the parameterK of VMD by the ratio of residual energy to original signal energy. To obtain the optimal input parameters of the RBF neural network, we employ PSR to determine the input dimension and the delay time of the RBF neural network. Using the ECG data from the MIT-BIH Arrhythmia Database as the data source, we evaluate the prediction performance of the proposed method. The RMSE and MAE of the proposed method are of 10-3 magnitude, while the RMSE and MAE of several common prediction methods are of 10-2 magnitude. The former is one order of magnitude smaller than the latter. The experimental results indicate that the proposed method is not only suitable for normal ECG signal prediction, such as No. 100 ECG, but also for abnormal ECG signal prediction, such as No. 207 ECG. Compared with other prediction methods, the proposed method is not only superior to some competitive prediction methods, such as [16, 17, 22], but also superior to some single prediction models, such as LS-SVM, SVM, ARMA, and Kalman. The proposed prediction method also outperforms some hybrid prediction methods, such as EMD-PSR-RBF and VMD-PSR-BP. WT-PSR-RBF is a competitive prediction method, its RMSE and MAE are 0.0052 and 0.0031, respectively, but the RMSE and MAE in this paper are smaller, only 0.0035 and 0.0027. As future work, we consider employing the prediction method proposed in this paper to reduce the energy consumption of the sensor in BANs.



---

*Source: 6624298-2021-03-16.xml*