can anyone tell me how to preprocess the ECG signal? I have used this code for Arrhythmia data from physionet.

if true
D = load('ECG.mat');
EKG = D.val;
Fs = 1000;                                                          % Sampling Frequency (Guess)
Fn = Fs/2;                                                          % Nyquist Frequency
Ts = 1/Fs;                                                          % Sampling Interval
t = linspace(0, 1, size(EKG,2))*Ts;                                 % Time Vector
Wp = [2  40]/Fn;                                                    % Passband
Ws = [1  60]/Fn;                                                    % Stopband
Rp =  5;                                                            % Passband Ripple
Rs = 20;                                                            % Stopband Ripple
[n,Wn] = buttord(Wp,Ws,Rp,Rs);                                      % Butterworth Filter Order
[b,a] = butter(n,Wn);                                               % Butterworth Filter Transfer Function Coefficients
[SOS,G] = tf2sos(b,a);                                              % Convert to Second-Order-Section For Stability
figure(1)
freqz(SOS, 4096, Fs)                                                % Assess Filter
EKGf = filtfilt(SOS,G,EKG');                                         % Filter EKGs
plot(t,EKGf)
  % code
end

[niyatha m]

hello i would like to know how you loaded ecg.mat and what is val in the dataset. Thank you.

Answer 1

[Star Strider]

That is a filter design I recognise (since I wrote it). It will remove baseline drift at the low end, and noise at the high end. I wrote it for a normal EKG, so it may not be appropriate for an arrhythmia EKG.

To filter an arrhythmia EKG (that requires a higher passband frequency), I would change that code to:

Fs = 1000;                                                          % Sampling Frequency (Guess)
Fn = Fs/2;                                                          % Nyquist Frequency
Ts = 1/Fs;                                                          % Sampling Interval
Wp = [2  97]/Fn;                                                    % Passband
Ws = [1  99]/Fn;                                                    % Stopband
Rp =  5;                                                            % Passband Ripple
Rs = 20;                                                            % Stopband Ripple
[n,Ws] = cheb2ord(Wp,Ws,Rp,Rs);                                     % Chebyshev Type II Filter Order
[z,p,k] = cheby2(n,Rs,Ws);                                          % Chebyshev Filter Transfer Function Coefficients
[SOS,G] = zp2sos(z,p,k);                                            % Convert to Second-Order-Section For Stability
figure(1)
freqz(SOS, 2^16, Fs) 

t = linspace(0, 1, size(EKG,2))*Ts;                                 % Time Vector
EKGf = filtfilt(SOS,G,EKG');                                        % Filter EKGs

This gives a much better filter (I’ve learned more since I wrote that code), and is compatible with an EKG with an arrhythmia.

If your EKG has 50 Hz or 60 Hz mains frequency noise, see the documentation tutorial on Remove the 60 Hz Hum from a Signal. It also presents a different way to design a filter.

[Neal]

Thanks for reply strider. Do you tell me how pre-processing of Sleep apnea EKG?

[Star Strider]

I would pre-process it as with any other EKG. Use a 100 Hz upper cutoff frequency (as I designed my filter) to allow for arrhythmias (assuming at least a 250 Hz sampling frequency).

Analysing it with the respiratory and other data afterwards will be the challenge.

[Neal]

why we are using 250 Hz sampling frequency for sleep apnea? previously assumed sampling frequency was 1000 Hz for NORMAl & Arrhythmiya. why this sampling frequency you have taken?

[Star Strider]

We aren’t using 250 Hz. What I thought I wrote was that my filter will work with any sampling frequency that is at least 250 Hz. (Many people read these posts, so I add details like that to my Answers and Comments. Some PhysioBank EKGs are sampled at 200 Hz or less, and would not work with this filter design.)

The filter I designed in my original Answer will work with your signal. It removes low-frequency baseline variations and high-frequency noise. It will not remove power-line frequency noise, so I included a link to a filter design that will do that in my original Answer.

[Neal]

Thanks,I have removed the baseline of Arrhythmia EKG data.I want to remove power-line frequency noise, I have used your given link ,it has not working properly. please design the filter for power-line noise removal. I have attached the EKG sample data with baseline removal of Arrhythmia EKG and raw EKG of sleep-apnea EKG data.

[Star Strider]

The power-line frequency filter is written for 60 Hz North America power frequencies. If you have 50 Hz power, you have to define the filter notch frequency for 50 Hz instead.

Try this:

d = designfilt('bandstopiir','FilterOrder',2, ...
               'HalfPowerFrequency1',49,'HalfPowerFrequency2',51, ...
               'DesignMethod','butter','SampleRate',Fs);

If you have some other power frequency, make the appropriate changes to the code to notch it out.

[Neal]

designfilt() function gives the error message, actually I am using MATLAB 2013. If u give me the full code what can I do after baseline removal, above I have attached EKG data.

[Star Strider]

The Code —

Fs = 1E+3;
Fn = Fs/2;
Wp = [48 52]/Fn;
Ws = [49 51]/Fn;
Rp =  5;
Rs = 50;
[n,Ws] = cheb2ord(Wp,Ws,Rp,Rs);
[z,p,k] = cheby2(n,Rs,Ws,'stop');
[sos,g] = zp2sos(z,p,k);

figure(1)
freqz(sos, 2^17, Fs)
set(subplot(2,1,1), 'XLim',[0 100])
set(subplot(2,1,2), 'XLim',[0 100])

This designs a stable Chebyshev Type II 50 Hz notch filter. A 60 Hz version has:

Wp = [58 62]/Fn;
Ws = [59 61]/Fn;

The rest of the code are unchanged.

[Neal]

I have a Sleep Apnea EKG collected from physionet whose sampling frequency is 100 Hz. Does it work for your above baseline removal code and 60 Hz notch filter that you have designed?

[Star Strider]

A signal with a sampling frequency of 100 Hz will not work with any of the filters I designed. The Nyquist frequency is 50 Hz, so you can only filter frequencies below that. You will have to design your own filters to work with that record.

[Neal]

Then how should I remove 60 hz noise from 100Hz sampling frequency EKG record?

[Star Strider]

You can’t.

You cannot design a filter for any frequency above about 49.8 Hz with a 100 Hz sampling frequency, since 50 Hz is the Nyquist limit.

This raises the possibility that you may not actually have any 60 Hz noise in your signal, or at least any you can actually identify. The reason is that 60 Hz is beyond the Nyquist frequency, so would either be eliminated by the hardware anti-aliasing filter in the ADC or be represented by a low-frequency aliased signal that you probably could not uniquely identify and eliminate without eliminating some of the information in your actual signal.

[Neal]

So, in this case how can I remove baseline?

[Star Strider]

You can remove low frequency baseline drift and d-c offset wtih an apppripriately-designed highpass filter. I would set the stopband at 0.5 Hz, passband at 1 Hz, and use a Chebyshev Type II filter design. (Depending on the frequency of the baseline variation, you might have to change these slightly.) You can use the code in my original Answer as a prototype for your filter.

[Neal]

Then, what should be my guess sampling frequency, and how can I choose the parameter of stopband-passband ripple, high passband cut-off frequency?

[Star Strider]

Your sampling frequency is the sampling frequency of your signal.

Use the parameters I used, except for the upper passband and stopband frequencies. For a highpass filter, those don’t apply. Specify 'high' for the filter type argument in the cheby2 call.

[Neal]

I used your above given parameter but there is no difference in filtered signal.

[Star Strider]

There may not be anything to filter out at those frequencies.

[Neal]

The below code has designed to remove low frequency, the original sampling frequency is 100Hz. I don't know this is right or wrong please verify it sir.

 Fs = 100;                                                           % Sampling Frequency (Guess)
Fn = Fs/2;                                                          % Nyquist Frequency                                                        % Sampling Interval
Wp = 1/Fn;                                                          % Passband
Ws = 0.5/Fn;                                                        % Stopband
Rp = 5;                                                            % Passband Ripple
Rs = 20;
[n,Wn] = cheb2ord(Wp,Ws,Rp,Rs);                                     % Chebyshev Type II Filter Order
[z,p,k] = cheby2(n,Rs,Ws,'high');                                          % Chebyshev Filter Transfer Function Coefficients
[SOS,G] = zp2sos(z,p,k);                                            % Convert to Second-Order-Section For Stability
figure(2)
freqz(SOS, 2^16, Fs) 
EKGf = filtfilt(SOS,G,EKG');

[Star Strider]

It appears correct.

You will get better results with this one change:

Ws = 0.8/Fn;

It is quite possible that there is no significant signal energy in the filter stopband, the reason it does not appear to have any effect. The filter should eliminate any d-c (constant) offset, regardless.

[Neal]

In MIT-BIH Arrhythmia database there are 48 records, these records are divided into normal and different types of arrhythmia EKG records. should I use the filter that you have designed for normal EKG in MIT-BIH Arrhythmia database.

[Star Strider]

With appropriate changes for the sampling frequency, it will work for most of them.

Whether you use it or another filter that you design is entirely up to you. It depends on what you want to do.

[Neal]

How to select the parameter of passband ripple(wp) & stopband attenuates (ws) in above filter design for same EKG?

[Star Strider]

They are arbitrary. The usual passband ripple is 1 dB, and the stopband ripple (or stopband attenuation) is 50 dB. You can always decrease the stopband attenuation further, however that increases the filter length, and so decreases the computational efficiency of the filter.

That Answer is a couple years old. A much more efficient filter design is:

Fs = 1000;                                                          % Sampling Frequency (Guess)
Fn = Fs/2;                                                          % Nyquist Frequency
Ts = 1/Fs;                                                          % Sampling Interval
Wp = [2  97]/Fn;                                                    % Passband
Ws = [1  99]/Fn;                                                    % Stopband
Rp =  1;                                                            % Passband Ripple
Rs = 50;                                                            % Stopband Ripple
[n,Wp] = ellipord(Wp,Ws,Rp,Rs);                                     % Elliptical Filter Order
[z,p,k] = ellip(n,Rp,Rs,Wp);                                        % Elliptical Filter Transfer Function Coefficients
[SOS,G] = zp2sos(z,p,k);                                            % Convert to Second-Order-Section For Stability
figure(1)
freqz(SOS, 2^16, Fs)

So I recommend the elliptical filter design. The Chebyshev filter that meets those specifications is 34-order, and the elliptical filter is 10-order.

[Neal]

How should I check the number of order returns?

[Star Strider]

That’s the ‘n’ output in ellipord (or cheb2ord).