How to filter ECG and detect R peaks

Question

I am trying to extract R peak from raw ECG data and some samples are seemed to be interfered by EMG. I used the lib provided by biosppy with python, biosppy.signal.ecg, it seems that the lib does,

1. Default filter is firwin with filtfilt, mode is bandpass with 3-45Hz, the filter order is 0.3 * sample_rate, in my case, sample_rate is 500Hz, and order is 150
2. The R peak extract algorithm is P.S. Hamilton, "Open Source ECG Analysis Software Documentation", E.P.Limited, 2002

The result seems not stable when handling some samples interfered by EMG,

Between 2 and 8, the R peaks seem to be not correct. I'm not familiar with DSP and what can I do to get the raw ECG from the composed signal?

What I tried(I read some papers for ECG filtering) and the issue encountered,

• Baseline Wander, some papers suggest to use firwin with highpass filter to remove this kind of noisy, and the frequency can be 0.5Hz, I tried this both with firwin and butter, the filtered signal was more smooth than before, but the EMG interfered still exists
• Powerline Interference, the bandstop frequency is 50Hz or 60Hz in most cases, since the firwin using frequency from 0.5Hz to 45Hz, this kind of noisy also filtered
• EMG Noisy, some paper mention to use moving average to remove this kind of noisy since this noisy is complex and may overlap with the raw ECG, some suggest to use a 8 window moving average. But in my case, the noisy last a long time with many samples, even the window set to 150, the EMG also exists

What I want is to get the stale R peaks without the interference? if the EMG is not easy to remove, can I remove the invalid R peaks and do interpolation after R peaks detected, such as,

If the R peaks detected is,

[439, 433, 433, 245, 193, 150, 414, 307, 135, 300, 432, 430, 429]

The first 3 and last 3 indices seem to be correct with the correct raw PQRST wave, but other indices are the invalid ones interfered by EMG, remove these indices, and insert 3 new indices with the value 433?

Thanks for your kindly help.

Answer 1

You should use py-ecg-detectors Siply install by doing

pip install py-ecg-detectors

Then you can use for instance the well known Pan Tompkins algorithm to find the R-peaks

Here I used an ECG recording from the Apnea-ECG Database

from ecgdetectors import Detectors
import matplotlib.pyplot as plt
import pandas as pd
fs=100 # sample freq

heartbeat = pd.read_csv("/myworkdir/mitbih-database/100.csv")

detectors = Detectors(fs)

r_peaks_pan = detectors.pan_tompkins_detector(heartbeat.iloc[:,2][0:1000])
r_peaks_pan= np.asarray(r_peaks_pan)

plt.plot(heartbeat.iloc[:,2][0:1000])
plt.plot(r_peaks_pan,heartbeat.iloc[:,2][0:1000][r_peaks_pan], 'ro')

If you want you can also add a R-peak correction algorithm:

def R_correction(signal, peaks):

    num_peak=peaks.shape[0]
    peaks_corrected_list=list()
    for index in range(num_peak):
        i=peaks[index]
        cnt=i
        if cnt-1<0:
            break
        if signal[cnt]<signal[cnt-1]:
            while signal[cnt]<signal[cnt-1]:
                cnt-=1
                if cnt<0:
                    break
        elif signal[cnt]<signal[cnt+1]:
            while signal[cnt]<signal[cnt+1]:
                cnt+=1
                if cnt<0:
                    break
        peaks_corrected_list.append(cnt)
    peaks_corrected=np.asarray(peaks_corrected_list)            
    return peaks_corrected 

corrected_R_peak=R_correction(heartbeat.iloc[:,2][0:1000],r_peaks_pan)
plt.plot(heartbeat.iloc[:,2][0:1000])
plt.plot(corrected_R_peak,heartbeat.iloc[:,2][0:1000][corrected_R_peak], 'ro')