Proceedings of the

International Conference on Mechanical Engineering 2011

(ICME2011) 18-20 December 2011, Dhaka, Bangladesh

ICME11-RT-001

SIMPLIFIED MATHEMATICAL MODEL for GENERATING ECG

SIGNAL and FITTING THE MODEL USING NONLINEAR LEAST
SQUARE TECHNIQUE
Md. Abdul Awal1, Sheikh Shanawaz Mostafa1 and Mohiuddin Ahmad2
1
Department of Biomedical Engineering, Khulna University of Engineering & Technology, Khulna
Bangladesh
2
Dept. of Electrical and Electronic Engineering, Khulna University of Engineering & Technology, Khulna,
Bangladesh

ABSTRACT
ECG models are complex and their computational time is high. In this paper, we propose a Gaussian
wave-based model which can simulate ECG wave as well as P, Q, R, S and T waves individually. In
addition, this model is capable of simulating various kinds of practical phenomena. The coefficient of the
model was calculated by nonlinear least square technique using Gauss-Newton algorithm. In order to
evaluate the effectiveness of the model, different kind time domain and frequency domain techniques such
as PSD and MSC were used. The goodness of fitting was calculated using MSE, NMSE, RMSE, NRMSE
and PRD and compared with real and model ECG signal. The lower value of these error and higher
cross-correlation coefficient of 0.9208 between model and real ECG indicates the outstanding performance
of the model. The model is also successful in generating noisy ECG signal.

Keywords: ECG Signal, Gaussian Wave, Nonlinear Least Square Technique, Goodness Of Fitting

1. INTRODUCTION used seven Gaussian functions for modeling of ECG by

The heart is a hollow muscular organ that beats in means of 3D state-space model which require numerical
rhythm to generate the force for pumping blood through integration using a fourth-order Runge-Kutta method. S.
the whole body. Each beats triggered by a bioelectric Paravena et al. [5] used a large number of Gaussian (4 to
signals originated at sinoatrial (SA) node and spreads 133) with no base line drift factor based on minimum
throughout the body. And the electrocardiogram (ECG or bank method and zero crossing method. But fitting this
EKG) is an investigative tool for a wide range of heart model to the real ECG signal, starting and end point of
conditions, from minor to life threatening that measure any interval using zero crossing method is not efficient.
and records the electrical activity of the heart in superb In addition, growing number of Gaussian functions
detail. So ECG signal modeling and processing is the involve much time to run the program. [6] proposed a
most significant topics in biomedical engineering. model using Gaussian function. However they cannot
For modeling of ECG, different techniques have been represent QRS wave individually as well as it is unable to
developed in the past. A pole-zero models of the ECG fit with the real ECG at a significant level. They used
was represented by [1] for feature extraction and data double differentiation of the Gaussian function which is
compression. Another research [2] reported that, the time consuming and need complex mathematical
poles and zeros form clusters and the clusters can be operation. The fitting techniques were incompetent
related to the constituent waves of the ECG models. because they were not capable to fit any negative values
Transform-type methods like nonlinear transform using in their model which was quite common in real data.
multiplication backward difference for detecting QRS This paper propose a Gaussian wave base model
proposed by different researchers. However, this types which can simulate ECG wave as well as P, Q, R, S and T
of modeling cannot provide a direct representation of wave individually and is very simple as compared to
the constituent waves in the ECG as cardiac specialist earlier mentioned model. Nevertheless, ECG signals are
are needed for making diagnoses. corrupted by various kinds of noises like other electrical
Chip Away Decomposition (ChAD) algorithm which signal, , such as (i) power lines interference [7], (ii)
is an iterative method for Gaussian parameter high-frequency electromyography (EMG) noise, (iii)
determination was used for decomposing and motion artifacts, (iv) impedance changes at the
representing the ECG model by [3]. Clifford et al. [4] skin/electrode, (v) baseline drifts [8], (vi) electrosurgical
© ICME2011 1 RT-011
noise [9], and (vii) white noise. Therefore, noisy ECG is each other wave. To solve this problem, ECG
important as normal ECG to generate realistic ECG components are fitted to their right position using
signal. In addition, noisy ECG can be generated by the shifting and zero padding method.
model. The coefficient of the model is calculated section N j ,SNR(t ) is the noise parameter in the model. j
by section by nonlinear least square technique using
gauss-Newton algorithm. The performance of the indicates various types of ECG noise, such as white and
proposed model and fitting algorithm are evaluated by color noise, muscle artifact (MA), electrode movement
using ECG from recorded data. Thereby, for illustrating (EM), and baseline wander (BW). SNR indicates the
model performance, goodness of fitting are calculated. input signal-to-noise ratio.
Real data are pre-processed for better fitting by using As real ECG is contaminated by these noises, so noise
Butterworth filtering .Finally; the frequency-domain should be taken into account for more realistic modeling.
analysis of ECG signals is demonstrated. Various noises are modeled in the following session.
The rest of this paper is organized as follows. Section 2
introduces the proposed model ECG model. Details (a)Synthetic noises
methods of the proposed model are presented in Section White noise is a random signal with a flat power spectral
3. Section 4 provides experimental results and features of density. It has all frequency components. Flicker noise or
the model and conclusion and future work come in color noise is a type of low-frequency electronic noise
Section 5. with an inversely proportional power spectral density
compared with the frequency. Resistance fluctuation is
2. ECG MODEL the main reason for flicker noise generation and that’s
ECG signal is resembling of the combination of bell why all resistors has flicker noise. For the current study,
curve like P, Q, R, S and T waves; it falls toward both we have modeled the noise color by a single parameter

wave. If i  P , Q, R, S, T  then Gaussian wave for each

sites which are one of the characteristic of Gaussian representing the slope of a spectral density function that
decreases monotonically with frequency by following Eq.

2
[10]

S( f ) 
component of ECG wave have following parameter: Mi

f
is height of curves peak, ti is the center position of the (7)
peak and Wi controls the width then ECG components
Where f is the frequency,  is the variance of the
can be written as 2
 t  P 
 
2

 2W 
original signal and β is the slope; a measure of noise color.
  White noise (β=0), pink or flicker noise (β=1), and brown
P wave : M P e P
(1)
 t  
noise or the random walk process (β=2), are three of the
 t  Q 1   
 
2 2

 2W  
 MQ 2 e
most commonly referenced noises .White noise and color
 2WQ
Q2

 Q1    noise are simulated having 3dB input SNR using

Q wave : MQ1 e 2
(2) MATLAB.
 t  R 
d   
2

2WR  (b) Real Noise

R wave : M R e (3) Real noises are extracted from the noise stress test
dt
 t  S 
database in MIT-BIH [11]. Low-amplitude muscle noise
  
2

S wave :  M S e

is common in ECG. Muscle noise is, in contrast to
 2WS  baseline wander and 50/60 Hz interference.
(4)
 t  T 
Mechanical movement of recording electrodes or skin
 
2

 2W 
stretching results in changes in potential. Due to
  alteration in the physical dimensions of the
T wave : MT e T
(5)
electrode-skin half cell thus modifying cell potential and
Simply the general equation (1-5) can be written as skin-electrode impedance [12]. The proposed model is

M e
 t tij 
tested on the ECG data recorded from BIOPAC data
 t ti 
 d

 
2
  j 2
2

   Mi e  
2 wi   wij 
fi 
i acquisition system [13]. The noisy signal was generated
 

 dt 
by adding BW, EM, and MA artifacts (Noise Stress Test
iQ , j 1
ij
iP , R,S ,T
(6) Database of MIT-BIH) and white and color noise were
added to the clean ECG signals [12].
 N j ,SNR(t ) 3. MATERIALS AND METHODS
 d
In Eq. (6),    depend on i.
i 3.1 ECG

 dt 
The experiment is performed to collect the data for
this research work. The subject is a male of 26 years old

      iff i  R and      1 iff i  S .

 d d  d
i i with no known cardiovascular disorder. For ECG and
 dt   dt   dt 
measurement, required equipments are BIOPAC
electrode lead set (SS2L), BIOPAC disposable vinyl
electrodes (EL503), BIOPAC data acquisition unit
Otherwise the term does not exist. (MP36) [12] with cable and power. For ECG
If we perform normal sum of Eq. (6), the individual P, measurement, white lead was placed on right forearm,
Q, R, S, T component of ECG will be overlapped with

© ICME2011 2 RT-011
red lead on the left leg and the black lead was placed on required to solve this problem.
right leg as shown in Fig. 2. Subjects was seated in a chair Consider that an ECG function y= f(x) of a variable of
relaxing and asked to be as still as possible to ensure x tabulated at i values, where y1=f(x1),
lower motion artifact and EMG signal on the data. After y2= f(x2)……………yi = (xi ). Moreover, assume that the

parameters f ( x; 1 , 2 ,............., j ) and the set of i

running calibration sequence ECG data was recorded. known analysis form the function depending on j

equation will be

y1  f ( x1 ;1 ,2 ,.........

...., j ) 
...., j )
y2  f ( x2 ;1 ,2 ,......... 


(8)

yi  f ( xi ;1 ,2 ,.........
...., j ) 


1 , 2 ,............., j which satisfies our model properly. At

We have to solve the equation to obtain the value of

first an initial value is picked for  k and defined

(a)

d k  yk  f ( xk;1 , 2 ,.........., j ) (9)

And then the estimation for the change d k needed to

reduce d k to 0

d k   d dl
j
df
(10)
l 1 l
xk

For k=1,2,….,i where   (1 , 2 ,..., j )

Fig 1. (a)MP36 Biopac (b)Placement of lead in ECG This element can be written as a i  j matrix of partial
measureemnt. derivatives of

d1 d j
df df
.... .... .....
3.2 Preprocessing
x1 x1
The acquired ECG data are preprocessed to remove
d1 d j
noise, artifacts, and baseline wander using df df
.... .... ....
Akl 
(11)
Savitzky-Golay Filtering [13]. To this end, two x2 x2
frequency-selective fourth-order Butterwort filters [14]  
are used: one high-pass filter with cutoff frequency at  
0.5Hz and one low-pass filter with cutoff frequency at 90
d1 d j
df df
.... .... ....
xi
Hz. To suppress the interferences from the power line
xi
grid, a notch filter centered around 50 was used. Again,
this filter is implemented as a fourth-order Butterworth Then,

d k  Akl dl
filter like [15].
As from the properties of Gaussian wave it is known (12)
that, the wave do not cross the zero. But P, Q, R, S and T
waves cross the zero line in the real ECG signal. To solve And the brief equation is

d  Ad
this zero crossing problem 0.05 was added to real ECG
signal before calculating the coefficients, so that the (13)
baseline of real ECG lies in zero line.
If by defining a= ATA and b= AT dβ
3.3 Nonlinear Fitting
a d  b
We find,
For fitting the mathematical model with the real world

Then Eq. (14) is solved for d using Gaussian

(14)
data statistical hypothesis testing like: test of normality
of residuals, chi square test, analysis of variance, least
elimination techniques. This offset is applied to α and a
new d is calculated. By interactively applying this
square test etc [16] is needed. As it observed that, ECG
model is nonlinear in the coefficients. So the nonlinear
procedure until the elements of d become smaller than
least square techniques can be the best choice for this
fitting. Nonlinear models are more difficult to fit than

residuals is calculated by R 2  d .d after the final

linear models because the coefficient cannot be estimated desired limit, a solution is obtained. The sum of square
using simple techniques instead an iterative approaches

© ICME2011 3 RT-011
 [ x(n)  x
iteration. N 1
2
m ( n)]
NMSE  n 0

[ x(n)]
(18)
N 1
3.5 Performance Evaluation Parameters
3.5.1 Coherence 2

n 0
The magnitude squared coherence (MSC) estimate is a
function of frequency with values between 0 and 1 that Another measurement is Root Mean Square Error, which
indicates how well the model ECG corresponds to real

[ x(n)  x
is
N 1
ECG at each frequency. The MSC estimate Cxy of the
RMSE 
input signals (x and y) using Welch's averaged, modified 1 2
m (n)] (19)
n 0
periodogram method [17]. The MSC is nothing but a N
function of the power spectral densities (Pxx(f) and Pyy(f))
of x and y and the cross power spectral density (Pxy(f)) of

[ x(n)  x
x and y. The Normalized version of RMSE is
N 1

Cxy 
2
2
Pxy ( f ) m ( n)]
NRMSE  n 0

[ x(n)]
(15) (20)
Pxx ( f ) Pyy ( f ) N 1
2

n 0
In this paper; x, y represents the model and real ECG Percent Root Mean square Difference (PRD) can be
signal respectively and x and y must be the same length.

 [ x(n)  x
determined by Eq. (26).
N 1
3.5.2 Power spectral density (PSD) m ( n)]
2

PRD  n 0
 100%
 [ x(n)]
PSD represents the strength of the energy as a (21)
function of frequency [17].In other words, it shows at N 1
2
which frequencies variations are strong and at which
n 0
frequencies variations are weak. Energy can be obtained
within a specific frequency range by integrating PSD
within that frequency range. Computation of PSD is done 4. RESULT AND DISCUSSION
directly by the method called FFT or computing The coefficient of the model is calculated by
autocorrelation function and then transforming it. nonlinear least square technique using Gauss-Newton
algorithm having 95% confidence level and the
3.5.3 Cross-correlation coefficient coefficient are shown in Table 1.
If x(n) be the recorded or collected ECG signal and
xm (n) be the ECG signal generated by the mathematical Table 1: Coefficient to fit the model with BIOPAC
recorded ECG
model, then cross-correlation coefficient ρ between
ECG Coefficient
x(n) and xm (n) is given by [18] Component Ai Bi ti

[ x(n )   x ][ xm (n )  m ]
P 0.185 17.8 236.9


-0.1103 30.64 26.8
 x m
j=1
(16) Q
j=2 -0.1075 5.705 76.23


R 50 29.87 55.71
Where denotes the average calculated by S 0.509 9.09 11.9
summing over the observed time series, indexed by n. T 0.3255 29.78 154.3

x(n), and m and  m are the mean and standard

where µ x and σx are the mean and standard deviation of
Using the value of the model coefficient, the real and
model ECG are quite similar. This is represented in fig 2.
deviation of xm (n ) . A value of ρ ∼ 1 reflects a strong
correlation, ρ ∼−1 implies a strong anticorrelation, and ρ 1.2
Real ECG
∼ 0 indicates that x(n) and xm (n ) are uncorrelated. This 1 Model ECG

means that a value of ρ = 1 suggests that model and real 0.8

ECG are identical. 0.6

Amplitude

0.4
3.5.4 Error Assessment:

[ x(n)  x
0.2
The Mean Square Error (MSE) is defined as [19].
N 1
MSE 
0
1 2 (17)
m (n)]
n 0
N -0.2

The normalized form of MSE is -0.4

-0.6
0 100 200 300 400 500 600 700 800 900
Samples

Fig 2. Comparison between real and model ECG signal

© ICME2011 4 RT-011
To evaluate the performance of the model, goodness of (Slow heart rate), Tachycardia (Fast heart rate) etc. From
fitting was calculated using Eq. (17) ~ Eq. (21) and the Fig. 5, brachycardia was simulated for 5 sec possessing 4
values are given in Table 2. It is observed that, the values beat. So the beat per minute (BPM) was 48. Sinus rhythm
of different errors were quite low and the value of and tachycardia were simulated in the same way.
cross-correlation coefficient was 0.9205 which indicated
the highly correlation between the real and model ECG.

Table 2: Goodness of fitting

Goodness Value
of fitting (BIOPAC)
MSE 0.00779
NMSE 0.172477
RMSE 0.0882615
NRMSE 0.029748
PRD 2.974839

The model was not only evaluated in the time domain but Fig 5. ECG simulation by model for 5 Sec. For
also in the frequency domain. Figure 3 and Fig. 4 showed brachycardia with BPM 48, sinus rhythm with BPM 72,
the power spectral density (PSD) and magnitude square tachycardia with BPM 108.
coherence (MSC) respectively. In Fig 3, the model and
real ECG are reasonably similar and the MSC in fig 4, Another feature of the model is that, it can generate
was exists 1 in most of the frequencies which is a better realistic and simulated noisy ECG signal. N j ,SNR(t ) of
indication of the proposed model.
Eq. (6) indicates the noise parameter of the model.
Power spectral density Different noisy ECG signal (different j) were simulated
0.25
Model ECG for the input SNR of 3dB in fig 6.
Real ECG
0.2
1.5
White noise Colored noise
1
0.15 0.5
PSD

0.1 -0.5

-1
1.5
Real muscle artifacts Real electrode movements
0.05 1

0.5

0
0
0 10 20 30 40 50 60 -0.5
Frequency (Hz)
-1
1.5
Real baseline wander Composite noise
1
Fig 3. Comparison of power spectral density between 0.5
real and model ECG signal. 0

-0.5

1 -1
0 2000 4000 6000
0 2000 4000 6000
Sample Sample
0.9
Fig 6. Different types of noisy ECG signal.
0.8

0.7 5. CONCLUSION
0.6 The proposed model is capable of replicating many
important features of the human ECG wave. A number of
MSC

0.5

 Instead of using too many parameters like other,

0.4
features and applications of the model are:
0.3
the proposed system can generate ECG and
0.2
capable of simulating various kind of practical

 It does not need three dimensional state spaces

0.1 phenomena such as brachycardia, tachycardia.
0
0 10 20 30 40 50 60
Frequency(Hz) which are difficult for realization and

 Noisy ECG signal can be modeled by simply

Fig 4. Magnitude squared coherence (MSC) between real simulation.
and model ECG signal.

 In frequency domain the real and modeled ECG

adding a noise parameter with the model.
The model not only simulate ECG beat but can

 A realistic ECG database can be created by

simulate different heart condition like Brachycardia showed almost same properties.

© ICME2011 5 RT-011
 fitting the model with individual subject’s ECG. of fractal signals from noisy measurement using
This can be used for further analysis and for wavelets” ,IEEE Trans. Signal Process., vol. 40, no.

 By only saving the coefficient of the model,

education purpose. 3, pp. 611–623, Mar. 1992.
11. www. physionet.org
ECG data compression can be possible. 12. Biopac Student Lab 3.7, BIOPAC System, Inc., 2008.
However, the lower value of error parameters like 13. Awal, M.A., Mostafa, S.S., and Ahmad, M.,
MSE 0.00779 indicates the effectiveness of this “Performance Analysis of Savitzky-Golay
simplified model and the proposed model enables us to Smoothing Filter Using ECG Signal”, International
investigate wave morphology variation. It can be Journal of Computer and Information Technology,
improved for better fitting, model-based denoising, vol. 01, no. 02, January 2011.
compression and neural network based classification etc 14. Butterworth, S., “On the theory of filter amplifiers,”
which is under investigation. Finally, it is hoped that this Exp. Wireless Eng., vol. 7, pp. 536–541, Oct. 1930.
simple model will provide an efficient tool for testing and 15. Vullings, R., Vries, B.D., and Bergmans, J.W.M.,
processing of the ECG signals with different level of “An Adaptive Kalman Filter for ECG Signal
noise and/or motion artifact. Enhancement”, IEEE Transactions on Biomedical
Engineering, Vol. 58, No. 4,pp 1094-1103, April
7. REFERENCES 2011.
1. Murthy, Ivaturi S. N., Rangaraj, M. R. , Udupa, K.J., 16. Press, W.H., Teukolsky, S.A., Vetterling, W. T.,
Goyal, A. K., “Homomorphic Analysis and Flannery, B. P., Numerical recipes in C+ +, 2nd ed,
Modeling of ECG Signals,” IEEE Trans Biomed Cambridge University Press 2002.
Eng., vol. BME-26 Issue:6 pp. 330 - 344,June 1979. 17. Welch, P. D. ,“The use of fast Fourier transform for
2. Murthy, I.S.N. , Prasad, G.S.S.D. , “Analysis of the estimation of power spectra: a method based on
ECG from pole-zero models,” IEEE Trans Biomed time averaging over short, modified periodograms,”
Eng., vol: 39, issue:7 pp. 741 - 751,July 1992. IEEE Trans. Audio Electroacoustic, vol. AU-15, pp.
3. Suppappola, S., Sun, Y., Chiaramida. S.A., 70–73, June 1967.
“Gaussian pulse decomposition: An intuitive model 18. Clifford, G. D., Azuaje F.,. McSharry, P. E.,
of electrocardiogram waveforms”, Annals of Advanced Methods and Tools for ECG Data
Biomedical Engineering. vol. 25, no. 2 pp. 252-260, Analysis, Artech house, inc. 2006.
Mar-Apr, 1997. 19. Manikandan, M.S., Dandapat, S., “Wavelet energy
4. Clifford, G.D., Villarroel, M., “Model-Based based diagnostic distortion measure for ECG”,
Determination of QT Intervals”, IEEE Computers in Biomedical Signal Processing and Control, vol. 2,
Cardiology, vol.33 (k53) pp. 357-360, Sep 2006. issue 2, pp 80-96, Apr 2007.
5. Parvaneh, S., Pashna, M., “Electrocardiogram
Synthesis Using a Gaussian Combination Model 8. MAILING ADDRESS
(GCM)”, Computers in Cardiology, vol. 34, pp. Md. Abdul Awal
621−624, 2007. Department of Biomedical Engineering, Khulna
6. Richardson, J., Haywood, L.J., Murthy, V.K., and University of Engineering & Technology, Khulna-9203,
Harvey,G., “A mathematical model for ECG wave Bangladesh
forms and power spectra”, Mathematical E-mail: [email protected]
Biosciences, vol. 12, issues 3-4, pp. 321-328, Dec
1971. Sheikh Shanawaz Mostafa
7. Sornmo, L., and Laguna,P., Bio-Electrical Signal Department of Biomedical Engineering, Khulna
Processing in Cardiac And Neurological University of Engineering & Technology, Khulna-9203,
Applications, Elsevier Academic Press ed.1, June 15, Bangladesh
2005 . E-mail: [email protected]
8. Tam, H., and Webster, J. G., "Minimizing electrode
motion artifact by skin abrasion", IEEE Trans. Mohiuddin Ahmad
Biomed. Eng., vol. 24, pp. 134-139, 1977. Dept. of Electrical and Electronic Engineering, Khulna
9. Begg, R., Lai, D.T.H., Palaniswami, M., University of Engineering & Technology, Khulna-9203,
Computational Intelligence in Biomedical Bangladesh
Engineering, CRC Press, 2008. Email: [email protected]
10. Wornell, W. G. and Oppenheim, A. V., “Estimation

© ICME2011 6 RT-011