R-DECO: an open-source Matlab based graphical user interface for the detection and correction of R-peaks

Pre-processing

The pre-processing consists of an enveloping procedure, which enhances the QRS-complexes and flattens the rest of the ECG (Varon et al., 2015). A visual explanation of the method is shown in Fig. 2.

First, the upper (U) and lower (L) envelopes are computed from the ECG signal by the secant method. This method selects the segment with the steepest positive and negative slope in a user-defined window with length t. Once U and L are obtained, they are used to derive a flattened version of the ECG signal (F): F = U − L. Since L is subtracted from U, the baseline is eliminated and only a positive signal, F, remains.

Decision

The locations of the QRS-complexes are found by detecting the peaks in the flattened ECG. These peaks are detected in three stages. First, all samples with an amplitude lower than the amplitude of the sample 80 ms further are selected. The 80 ms step size was experimentally defined. This results in the selection of the upward slopes. As a second step, only the upward slopes that are longer than the step size are selected in order to exclude small peaks. Finally, the maximum is selected in a window, with a length equal to the step size, that starts from the last selected sample of the upward slope. A graphical representation of this process is shown in Fig. 3.

On this selection of peaks, the adaptive thresholding procedure of the Pan–Tompkins algorithm is applied to define the peaks that correspond to the QRS-complexes.

Post-processing

The thresholding procedure generally produces satisfactory results for the detection of the QRS-complexes. However, some of the automatically generated RR-intervals might be physiologically unreasonable and need to be removed for further analysis. A slightly modified version of the search-back procedure as proposed in De Chazal et al. (2003) was used for this purpose.

Once the positions of the QRS-complexes are identified, the original ECG is used to find the exact location of the R-peaks. The search for an R-peak is performed up to 50 ms from the peak in the flattened signal. This extra search is necessary because the presence of large S-waves might shift the peak in the flattened signal toward the valley of the S-wave.

Evaluation on the physionet MIT/BIH arrhythmia database

We used the MIT/BIH arrhythmia database to evaluate the proposed algorithm (Moody & Mark, 2001). This dataset consists of 48 half-hour ECG signals, which were recorded in the Boston’s Beth Israel Hospital between 1975 and 1979. All recordings were annotated by two independent cardiologists who also made a distinction between normal and abnormal beats. In total, 110,122 heartbeats were annotated, of which 89,133 were labeled as normal. Each recording contains two channel ECG signals with a sampling frequency of 360 Hz. In most records, one channel is lead II and the other channel is V1. However, we only used the first channel for the evaluation.

As mentioned previously, the pre-processing consists of a flattening step of the ECG with a user-defined window width. To evaluate the sensitivity of the performance to the choice of the width we have tested multiple window widths. As can be observed from Fig. 4, comparable results were obtained for window widths between 250 and 350 ms.

In Table 1 we list the R-peak detection results of the proposed algorithm with an envelope width of 300 ms and without post-processing. We obtained an overall sensitivity of 99.60% and PPV of 99.69%. When including the post-processing, we obtained an overall sensitivity of 99.09% and PPV of 99.80%. These results are comparable with those in literature, especially with the Pan–Tompkins algorithm, which reaches a sensitivity of 99.76% and a PPV of 99.56% (Elgendi et al., 2014).

While these results are very promising, we can also observe that for some recordings only moderate detection results are achieved. This decrease in performance is generally due to loss of signal, unusual morphology or stretches of extremely irregular rhythms. For instance, recording 116 and 208 contain stretches where the signal is lost in the first channel. However, the recordings with the highest amount of false detections are 108, 203 and 210. Correctly detecting the R-peaks in recording 108 has been proven difficult for many algorithms (Pan & Tompkins, 1985). It contains a lot of baseline wander and additionally, very tall and sharp P-waves. These characteristics make it difficult to distinguish P-waves from R-peaks and thus result in a high false positive count. However, the highest amount of false positives is observed in recording 203 (92.25% sensitivity). This might be explained by the extremely high percentage of premature ventricular contractions (PVC) present in the recording, almost 15%. Since the envelope width was fixed during the detection process, one may assume that the performance could be improved if manual adjustments were permitted. This holds as well for other records with PVC’s, such as record 210.

The noise tolerance of the algorithm was evaluated with the MIT-BIH Noise Stress Test Database (Moody, Muldrow & Mark, 1984). We observed that both median PPV and sensitivity remained around 100% above a signal-to-noise-ratio of six dB. From this threshold the performance of the algorithm decreased significantly.

From the analysis of the results, we could deduce two main factors that influence the results of the algorithm: (1) the envelope width and (2) the RR-post processing. The number of samples in the envelope is important, since it can be regarded as a filter of the RR-intervals. Smaller envelope widths might result in the enhancement of more peaks than only the R-peaks. This might be beneficial in the case of small R-peaks, but might also enhance artefact peaks. Larger envelope widths might cause adjacent R-peaks to be merged in the flattened signal. In practice, this might result in the failure of detecting premature heartbeats, which appear shortly after the previous heartbeat. In summary, a larger envelope width results in less false positives and more false negatives and the opposite is true for a small envelope width. A similar effect can be observed when the RR-intervals are post-processed. This increases the certainty of detection of the algorithm and thus results in less false positives. The downside is that, in the presence of abnormal rhythms, it also results in more false negatives.

Data

In the data panel, the user has the option to load data in two different ways: from a file or from the MATLAB workspace. Both can be accessed via the respective pushbuttons.

When a file is selected, the software visualizes a small segment of the ECG signal. If the signal is inverted, the user can indicate this and the file will be inverted before further analysis. The software also scans the selected file for the sampling frequency. If this is not found, the software prompts the user to manually indicate the sampling frequency.

Finally, the data panel also contains a reset button. This button allows the user to reset the entire GUI. It empties all plots, restores all default variables and deletes all results. If the user has not yet saved the current analysis, the user is prompted to confirm the reset action to prevent unwanted loss of information.

Filter

Since ECG signals can be contaminated with noise, filtering is often essential for further analysis. R-DECO provides three basic filters: high pass, low pass and a notch filter.

The high pass filter consists of a zero phase, second order Butterworth filter. The low pass filter consists of a zero phase, fourth order Butterworth filter. Finally, a zero phase notch filter is also included. The latter could be used to remove the power-line interference. Important to note is that filtering actions are always executed on the original signal to ensure repeatability.

In order to aid the user in the selection of appropriate frequency threshold(s), R-DECO is able to compute and display the power spectrum. An estimate of the power spectrum is computed using the Welch method (Welch, 1967). As a default, we used a window of 500 samples with 60% overlap.

To visualize the effect of the filtering in the frequency domain, R-DECO displays both the filtered and the original power spectrum. Furthermore, the effect of the filtering in the time domain can also be investigated by checking the “Show Original” checkbox. This will overlay the original signal on top of the filtered signal (Fig. 6).

Analysis period

In the Analysis Period panel the user has the possibility to define an analysis window. After pushing the “Define analysis period” button, the user has to select a window by clicking, dragging and releasing the mouse. The window is shown as a transparent patch over the data and can be enlarged, shrinked or moved with the mouse. After the initial window is drawn, the user can further modify the analysis window by changing the start time and/or the duration of the window.

When a window is selected, the user can accept the analysis window by pressing the “Apply changes” button. This prompts the x-limits of the graph to match the analysis window and disables the window modifications. From here on, all further analysis will be performed solely on the selected window.

All the above is very useful when the to-be-selected time period is known in advance, but this is not always the case. Sometimes the selection of the analysis window is depending on certain patterns in the tachogram, hence the tachogram has to be constructed first. Therefore, if R-peaks have been detected already, the user is prompted to indicate whether he/she would like to keep the detected R-peaks.

R-peak detection

The execution of the R-peak detection algorithm, as described in the first section, is initiated when the “Detect Peaks” button is pushed. However, before the actual algorithm is executed, the user is able to adjust the default parameters of the algorithm.

In Fig. 7, an epoch of 10 s of the ECG signal and its respective enveloped signal is shown as an example of the flattening procedure. The user can adjust the envelope size to the desired value and enact the changes by pressing the “Apply” button. Additionally, an estimation of the average heart rate can be defined to compute the boundaries of normal-to-normal RR-intervals. Lastly, the user can select the automatic post-processing of the RR-intervals.

Since some devices have built-in QRS detection algorithms and some researchers have their own preferred QRS detection algorithm, the software allows to load R-peak locations. These will be displayed the same way the R-peaks of the algorithm are displayed.

R-peak correction

In case of HRV studies, only the normal-to-normal RR-intervals need to be taken into account. This can be achieved by selecting the ectopic removal option. This option corrects ectopic beats, without altering the normal RR-intervals.

After finishing the detection process, either by detecting or loading the R-peaks, it is still possible that not all R-peaks are accurately detected. In R-DECO, the user can make manual and (semi-)automatic adjustments to the R-peak locations.

The manual methods are: add, adjust and delete. These methods can be activated by selecting the specific radiobuttons or via a context menu, which is linked to each R-peak. These manual methods allow the user to correct wrong or missing annotations either in all or in individual leads.

• Add: When this radiobutton is active, new R-peaks can be added by clicking in the ECG graph. The program selects the signal sample that is closest to the mouse position in a symmetric window of 300 ms around the mouse position. Upon mouse release, a new R-peak is added and the tachogram is adapted.

• Adjust: While hovering over the ECG graph, the R-peak closest to the mouse location is selected. After clicking on the desired R-peak, it can be moved by dragging the mouse. The movement of the selected R-peak is restricted by the previous and next R-peak. While adjusting the R-peak location, the tachogram is automatically updated. Upon mouse release, the new R-peak location is saved.

• Delete: While hovering over the ECG graph, the R-peak closest to the mouse location is selected. After clicking on the desired R-peak, more to-be-deleted R-peaks can be selected by dragging the mouse. Upon mouse release, the selected R-peaks are removed and the tachogram is adapted.

The three (semi-)automatic R-peak correction methods are: cross-correlation, peak conversion and maximum search.

• Cross correlation: For this method, a symmetrical window of 300 ms around each R-peak is selected. Then, all heartbeats are normalized by subtracting the mean and dividing it by the standard deviation. Then, a trimmed average QRS-complex is computed of all the positive and “negative” R-peaks. In this work a “negative” R-peak is understood as the absence of an R-peak or the presence of a very prominent S-wave, also described as RS-complex.

• The user is prompted to select either the positive or the “negative” average heartbeat. If necessary, the user can also adjust the location of the R-peak on the selected template. This is all graphically displayed as shown in Fig. 8. Finally, the cross-correlation of every heartbeat is computed with the trimmed average and the R-peak is re-located, based on the highest correlation value.

• Peak conversion: The absolute amplitude of every R-peak annotation is compared with the previous and next sample’s absolute amplitude. If it is bigger than the previous and smaller than the next, the location of the R-peak will be shifted forward, until an extremum is obtained. If it is the other way around, the location of the R-peak will be shifted backwards.

• This functionality avoids the “jumping” of R-peaks from, for example, an actual R-peak to a pre-mature ventricular contraction, which might be the case when window search is applied. Furthermore, pressing the button more than once will not affect the relocation after one correction.

• Maximum search: Firstly, a symmetric window of 300 ms around each R-peak is selected. Secondly, the user is prompted to select either the maximum, minimum or absolute maximum. Based on this selection, the respective extremum within the window is selected as new R-peak location.