Electrocardiogram arrhythmia detection with novel signal processing and persistent homology-derived predictors

1.Introduction

Cardiovascular diseases are among the leading causes of death per the World Health Organization and the Centers for Disease Control and Prevention [8,64]. Arrhythmias are heart rhythms other than normal sinus rhythm with a heart rate between 60 beats/minute and 100 beats/minute; that is, arrhythmias are heart rhythms that are either too fast, too slow, abnormal, and/or irregular. Most arrhythmias must be treated since they can either lead to 1) more chaotic electrical activity of cardiac muscle resulting in loss of cardiac output and/or 2) the formation of thromboemboli (e.g. as in atrial fibrillation) possibly resulting in stroke [40]. The overall prevalence of arrhythmias among adults is estimated to be around 2% with atrial fibrillation being among the most common arrhythmias [13,30]. The global prevalence of atrial fibrillation has been estimated to be about 0.51% [35].

The contraction and relaxation of cardiac muscle cells is driven by ion movement across cell membranes and must be coordinated in order for the heart to pump blood effectively. This ion movement is governed by an electrochemical potential comprised of 1) ion concentration gradients and 2) electric potentials. The depolarization and subsequent repolarization of cardiac muscle cells causes changes in electric potential on the body surface which can be measured non-invasively using an electrocardiogram (ECG). ECG analysis is important for accurate diagnosis, treatment, and prevention of cardiovascular diseases.

Topological data analysis (TDA) refers to a collection of methods concerned with quantifying ‘shapes’ of data which are invariant under continuous deformations such as stretching and twisting. The main tool of TDA is persistent homology which quantifies the homology of structures within the data which persist over a range of scales. Persistent homology has been applied to many tasks across various fields such as electroencephalogram analysis [3], genomics [4,7,14,37,44,50,57,63], classifying skin lesions based on images [10], and tumor segmentation on histology slides [49]. Cycle representatives – which will be described in Section 1.1 – of topological features have shown utility in various fields outside of ECG analysis such as analyzing structures on the atomic scale [47] and in structural engineering [24].

Several approaches to computer-aided ECG rhythm classification have been performed, including neural networks [5,15,17,21,25,39,46,48,51,56,61,62,66–68], wavelet transformation and independent component analysis [31,65], using higher-order statistics of wavelet-packet decomposition coefficients as features [32], and support vector machines using projected and dynamic ECG features [9]. An overview of TDA applied to cardiovascular signals has recently been performed [23]. In the field of computer-aided ECG analysis, TDA has been used to construct metrics of heart rate variability [11,20]. Additionally, the Mapper algorithm has been applied to predict the presence and severity of heart disease [2]. Computer-aided ECG rhythm classification methods which utilize TDA include neural networks with topological-based features [16,53], fractal dimension in tandem with neural networks [55], mapping ECG signals to a higher dimensional space prior to computing topological features [26,27,34,36,41], and utilizing a sliding window and Fast Fourier Transform to process the ECG signal prior to computing topological features [43]. These approaches construct topological predictor variables utilizing information directly derived from the birth and death radii statistics along with extra information such as heart rate, fractal dimension statistics, and persistent entropy.

To the author’s knowledge, constructing predictor variables for use in machine learning models to classify ECG rhythms based off information derived from cycle representatives has not yet been performed. Additionally, to our knowledge, there has been no computer-aided ECG analysis which utilizes only the N-most persistent topological features for use in rhythm classification, nor has there been an approach which introduces an isoelectric baseline into ECG signals to yield non-trivial topological features corresponding to P, Q, S, and T-waves (if they are present to begin with). Introducing an isoelectric baseline prior to computing persistent homology and utilizing the N-most persistent topological features and properties of their area-minimal cycle representatives for use in constructing predictor variables makes the approach taken here distinct from other combinations of TDA and machine learning described in the literature.

In Section 1.1, we give a brief overview of the aspects of persistent homology utilized in this study. Appendix A formalizes the intuition underlying persistent homology described in Section 1.1. The Methods portion is split into three parts: Section 2.1 describes the novel ECG processing pipeline, Section 2.2 describes the construction of predictor variables primarily based off the topological features of the processed ECG signal, and Section 2.3 describes the specific classification tasks along with the statistical models and evaluation metrics used. The Results/Discussion section presents the evaluation metrics and ROC curves for each statistical model used. The Conclusion section contains a brief comparison between the method proposed here and other methods which use TDA and machine learning for rhythm classification in addition to describing some future directions.

1.1.Intuition behind persistent homology

The background on persistent homology presented both here and in Appendix A is restricted to two-dimensional data and one-dimensional homology features. The methods discussed generalize to higher dimensions, but we restrict our focus to the relevant dimensions used in the ECG analysis presented here. A toy example dataset X and its persistent homology are used to build some intuition for persistent homology. The informal treatment of persistent homology described in this section is made rigorous in Appendix A.

Fig. 1.

Example dataset with persistence diagram. A: example dataset; B–E: radius 0.57, 0.81, 1.55, and 3.01 Geometric Čech Complex depicted in black, respectively; F: persistence diagram of equivalence classes of non-contractible loops.

Consider the set of points in the plane R2 shown in Fig. 1. Consider drawing a circle around each point, each with the same radius r. We will refer to the union of these circles as the Geometric Čech Complex of radius r, denoted C˘r(X), not to be confused with the Čech complex of radius r, which commonly refers to an abstract simplicial complex. Observe that for r<0.57, none of the circles comprising C˘r(X) overlap around a “void” of non-overlapping space. Furthermore, observe that for r∈[0.57,0.81), the circles comprising the smaller loop of points nearby the point (1,1) overlap such that there is a “void” of non-overlapping space enclosed by their region of overlap. Hence for r∈[0.57,0.81), there exists a non-contractible loop within C˘r(X). “Non-contractible” here means that the loop drawn around the void of non-overlapping space cannot be continuously deformed down to a single point without leaving C˘r(X); that is, the loop gets “stuck” on the void encircled by C˘r(X). This non-contractible loop can be continuously deformed to construct another non-contractible loop “stuck” around the same void. These two non-contractible loops are homotopic to one another. For example, the green and red loops in Fig. 1 are homotopic. The set of all possible non-contractible loops “stuck” around some void encircled by C˘r(X) forms an equivalence class of non-contractible loops, i.e. a set of non-contractible loops where any two non-contractible loops in the set are homotopic. In practice, rather than homotopy – of which is relatively straightforward to garner intuition in the context of TDA – we use a weaker but more technically-involved equivalence relation on loops called homology to utilize efficient algorithms such as Ripser [6] and GUDHI [38] in computing topological features. For a rigorous treatment of homotopy and homology, see [22].

For a given two-dimensional dataset X such that there exists a non-contractible loop ℓ within C˘r(X), we define the birth radius of the equivalence class of non-contractible loops containing ℓ as the smallest real number b such that some loop in C˘r(X) which is equivalent to ℓ and which is contained in the subset C˘b(X) of C˘r(X) exists. Similarly, we define the death radius of the equivalence class of non-contractible loops containing ℓ as the smallest real number d such that r⩽d and such that ℓ becomes contractible when regarded as a loop in C˘d(X). That is, the birth radius of an equivalence class of non-contractible loops is the smallest radius at which the equivalence class of that non-contractible loop forms, and the death radius is the smallest radius at which it vanishes (i.e., becomes contractible). For r∈[b,d], the equivalence class of non-contractible loops ‘persists,’ and this motivates the definition of the persistence of an equivalence class of non-contractible loops as the difference between the death radius and the birth radius. The two non-trivial equivalence classes of non-contractible loops in Fig. 1 have coordinates (0.57,0.81) and (1.55,3.01) in the persistence diagram and correspond to the subset of data clustered near (1,1) and the subset of data clustered near (8,8), respectively. Note that the larger loop-like structure of data in the upper-right corner of each subplot has a larger persistence than the smaller loop-like structure of data in the lower-left corner of each subplot (i.e. 3.01−1.55=1.46>0.81−0.57=0.24).

Fig. 2.

Relationship between Geometric Čech Complex of radius r and geometric realization of radius r Vietoris Rips Complex. A–C: Geometric Čech Complex of radius 0.2, 0.5, 0.71 depicted in black, respectively; D–F: geometric realization of radius 0.2, 0.5, 0.71 Vietoris Rips Complex, respectively.

The cycle representatives of a given equivalence class of non-contractible loops {ℓα}α∈I (note that I is an uncountable indexing set) with birth radius b and death radius d are the subsets of the data which give rise to non-contractible loops in C˘r(X) with birth radius b and death radius d. For example, the cycle representatives of the equivalence class of non-contractible loops with birth radius 0.5 and death radius 22≈0.71 in Fig. 2 are given by {{a,b,c,d},{a,b,c,d,e}}. The Python package Homcloud can be used to identify cycle representatives which are optimal in some sense such as having the minimum number of points or spanning the minimum area among all cycle representatives [45]. Associating a single optimal cycle representative to each equivalence class of non-contractible loops is important 1) for reproducibility and 2) to select cycle representatives which more closely resemble the P, Q, S, and T-waves for the relevant equivalence classes of non-contractible loops.

2.Methods

The free and publicly available Shaoxing Hospital Zhejiang University School of Medicine electrocardiogram (ECG) database was used in this study [69]. This database consists of 10646 12-lead ECG signals, each spanning 10 seconds with a sampling frequency (i.e. the number of electric potential differences recorded per second) of 500 Hz, of which 10605 have non-empty Lead 2 signals. This study strictly utilizes Lead 2, i.e. the ‘rhythm lead’, so the term ‘ECG signal’ is henceforth used to refer to Lead 2 ECG signals. Each ECG signal is labeled with one of 11 rhythms by professional experts. The distribution of these 11 rhythms across the 10605 ECG signals is shown in Table 1.

ECG signals are typically characterized as 1-dimensional lists of real numbers of length F·tmax where F is the sampling frequency of the ECG machine (i.e. the number of electric potential differences recorded per second), tmax is the total amount of time (in seconds) over which the signal was gathered, and each real number in the list represents the signal amplitude at the given time index. In order to compute 1-dimensional topological features of an ECG signal, the ECG signal must be considered as a subset of R2. Therefore, rather than treat a given ECG signal S as a one-dimensional list with a sampling frequency F over a length of time tmax, we use the equivalent formulation of S given by S={(t,f(t))|t∈D}⊂R2 where D={iF|i∈{1,…,F·tmax}} represents the set of time indices and f:D→R defines the signal amplitude at each time index.

In the remainder of this section, we describe 1) ECG signal processing prior to extraction of topological features, 2) the construction of predictor variables derived from persistent homology, and 3) the statistical modeling approaches and evaluation metrics used. A flowchart providing an overview of our approach to arrhythmia detection is shown in Fig. 3.

Fig. 3.

Flowchart of ECG signal processing and arrhythmia classification.

2.1.Electrocardiogram signal processing

Given a raw ECG signal S={(t,f(t))|t∈D}⊂R2 with time domain D={hF|h∈{1,…,F∗tmax}} and signal amplitude given by f:D→R, the signal is first normalized by applying the transformation g:f(D)→[0,1] given by:

(1)g(f(t))=f(t)−min{f(D)}max{f(D)}−min{f(D)}.

The resulting signal Snormalized={(t,g(f(t)))|t∈D}⊂R2 has maximum amplitude

max{g(f(D))}=1 and minimum amplitude min{g(f(D))}=0. Since equivalence classes of non-contractible loops are not scale-invariant, this normalization is necessary for the magnitude of persistent homology-derived statistics to be comparable across ECG signals.

Next, an isoelectric baseline is included in Snormalized in order to form ‘loop-like’ structures with non-trivial topological properties in the ECG signal corresponding to the P, Q, S, and T-waves (if they are present). The inclusion of this baseline emphasizes the shape of the P, Q, S, and T-waves (if they exist to begin with), as illustrated in Fig. 4. This is done by inserting the baseline value computed as the median of g(f(D)) at the beginning of the signal and between every pair of consecutive time indices, doubling the number of points of the signal while still spanning the same amount of time. More explicitly, after the inclusion of the isoelectric baseline to Snormalized, we obtain the signal:

(2)Sprocessed={(t,h(g(f(t))))|t∈E},E={i2F|i∈{1,…,2·F·tmax}},h:[0,1]→[0,1]:g(f(iF))↦median{g(f(D))}if i is oddg(f(iF))if i is even.

Note the appearance of highly-persistent equivalence classes of non-contractible loops around birth radius 0.005 once the isoelectric baseline is included in the normal sinus rhythm ECG signal in Fig. 5. Also observe in Fig. 5 that for a rhythm such as atrial fibrillation with the property of absent/attenuated P-waves, the isoelectric baseline does not produce additional highly-persistent equivalence classes of non-contractible loops to the same extent that it does for rhythms with normal wave-shape such as normal sinus rhythm.

Fig. 4.

Illustration depicting the effect of the isoelectric baseline on the persistence diagrams of ECG signals. PD: persistence diagram. A–B) normal sinus rhythm ECG signal without baseline and corresponding PD. C–D) same as A–B but with the isoelectric baseline included. Note the cluster of topological features that appeared and the P, S, and T-waves their area-minimal cycle representatives correspond to. E–F) atrial fibrillation without baseline included and corresponding PD. G–H) same as E–F but with the isoelectric baseline included. H1 features: equivalence classes of non-contractible loops.

The onset of each QRS-complex in the processed ECG signal Sprocessed is identified using Zong, Moody, and Jiangs’ approach of “passing Sprocessed through a low-pass filter, applying a transformation with a non-linear scaling factor to enhance the QRS-complexes and suppress unwanted noise, and applying adaptive thresholds to the signal to determine the onset of each QRS-complex” [70]. An illustration of the preprocessing transformations applied to a raw signal is shown in Fig. 5.

Fig. 5.

Depiction of preprocessing transformations applied to a normal sinus rhythm ECG signal. A) raw ECG signal with normal sinus rhythm. B) normalized ECG signal with maximum amplitude 1 and minimum amplitude 0. C) normalized ECG signal with isoelectric baseline included and R-waves identified.

2.2.Construction of predictor variables

Each equivalence class of non-contractible loops with birth radius b and death radius d corresponds to a set of subsets of Sprocessed given by Yprocessed={S⋆⊂Sprocessed|(birth radius of S⋆=b) and (death radius of S⋆=d)}. That is, there may be multiple subsets of data which generate a given equivalence class of non-contractible loops. Equivalently, for a single point in a persistence diagram, there may be multiple subsets of data such that the Geometric Čech complex births and vanishes the given equivalence class of non-contractible loops with the same birth and death radii. The Python package Homcloud is used to compute a single unique area-minimal cycle representative S⋆ from Yprocessed for each equivalence class of non-contractible loops in the signal Sprocessed [45]. Given an equivalence class of non-contractible loops with centroid coordinates of the area-minimal cycle representative (T,A), the effective centroid coordinates (x,y) are computed as

• x=tR−T where tR is the time-coordinate of the onset of the subsequent QRS-complex.

• y=A−baseline1−baseline where baseline represents the amplitude value of the isoelectric baseline median{g(f(D))}.

The computation of the effective centroid coordinates of an area-minimal cycle representative is depicted in Fig. 6. The equivalence classes of non-contractible loops with centroid time coordinate T larger than the largest of all onsets of the QRS-complexes are not considered to ensure that the effective centroid time coordinates can always be computed. This effectively trims Sprocessed to end with a point representing the onset of a QRS-complex. Furthermore, all equivalence classes of non-contractible loops with area-minimal cycle representative centroid amplitude coordinate A larger than 1−baseline2 where baseline=median{g(f(D))} are discarded to obtain a larger proportion of highly-persistent equivalence classes of non-contractible loops corresponding to clinically-relevant subsets of ECG signals such as P, Q, S, and T-waves. For example, the computation of the effective centroid time coordinate for area-minimal cycle representatives that represent P-waves is a proxy of the clinically-relevant PR-interval. The computation of the effective centroid amplitude coordinate normalizes the amplitude coordinates of centroids of area-minimal cycle representatives across signals with differing isoelectric baselines.

Fig. 6.

Computation of the effective centroid coordinates of an area-minimal cycle representatives. A) processed ECG signal with normal sinus rhythm and R-waves, an area-minimal cycle representative corresponding to a P-wave, and an area-minimal cycle representative corresponding to a T-wave identified. B) zoomed-in region depicting the computations of the effective time-coordinates of the two area-minimal cycle representatives. C) zoomed-in region depicting the computation of the effective amplitude-coordinate of the area-minimal cycle representative corresponding to the T-wave.

The persistent homology of the processed signal Sprocessed is then computed, and the N most persistent equivalence classes of non-contractible loops are used to construct predictor variables for use in rhythm classification for N∈{5,6,…,29,30}. Specifically, for each of the N-th most persistent equivalence classes of non-contractible loops, the persistence, birth radius, effective time-coordinate of the centroid of the area-minimal cycle representative relative to the subsequent QRS-complex, effective amplitude-coordinate of the centroid of the area-minimal cycle representative relative to the isoelectric baseline, and Shannon entropy of the vector (a,b,c,d,e)sum((a,b,c,d,e)) where a=persistence, b=birth radius, c=death radius, d=centroid time-coordinate, and e=centroid amplitude-coordinate are used as predictor variables. Additional predictor variables include the mean and standard deviation of the persistences, birth radii, area-minimal cycle representative centroid time coordinates, and area-minimal cycle representative centroid amplitude coordinates of the N-most persistent equivalence classes of non-contractible loops along with the mean and standard deviation of the RR-intervals. Lastly, the total number of R-waves, the total number of equivalence classes of non-contractible loops, and the Shannon entropy of the normalized distribution of all persistences are also used as predictor variables. Note that death radii statistics are not included as predictor variables since their inclusion would introduce undesired collinearity due to the persistence of a given equivalence class of non-contractible loops being the difference between the death radius and the birth radius.

2.3.Statistical modeling and evaluation

Three different binary classifications are carried out:

• Atrial Fibrillation vs. Non-Atrial Fibrillation

• Arrhythmia vs. Normal Sinus Rhythm

• Arrhythmias with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia

For each of the three binary classifications, Logistic Regression, Linear Discriminant Analysis, Quadratic Discriminant Analysis, Naive Bayes, Random Forest, Gradient Boosted Decision Tree, K-Nearest Neighbors, and Support Vector Machine with Linear, Radial, and Polynomial Kernel Models are constructed. For background on the theory and/or implementation of these statistical models, see [28]. Stratified 5-fold cross-validation is performed, and in each of the 5 folds, the true positives (TP), false positives (FP), false negatives (FN), and true negatives (TN) are recorded in a confusion matrix like that shown in Fig. 7 for each statistical model used. The mean and standard deviation of the F1-Scores, Accuracies, Sensitivities, Specificities, Positive Predictive Values (PPVs), and Negative Predictive Values (NPVs) across the five folds are recorded. Definitions of these evaluation metrics can be found in [52].

Fig. 7.

Confusion matrix.

The optimal hyperparameters for the Random Forest, Gradient Boosted Decision Tree, K-Nearest Neighbors, and Support Vector Machines with Radial and Polynomial Kernel Models were chosen as the hyperparameters which yielded the largest mean F1-Score across all folds in 5-fold stratified cross validation. The grid search spaces of hyperparameters for the relevant models are:

• Random Forest:

  • ∗ number of trees∈{500,1250,2000,3000}

  • ∗ number of variables randomly sampled∈{int(0.25·T),int(0.5·T),int(0.75·T),T} where T is the total number of predictor variables.

• Gradient Boosted Decision Tree:

  • ∗ number of trees∈{500,1250,2000,3000}

  • ∗ interaction depth∈{5,10,15,20}

• K-Nearest Neighbors:

  • ∗ K∈{1,2,3,4,5,6,7,8,9,10}.

• Support Vector Machine with Radial Kernel:

  • ∗ cost=1.

  • ∗ γ∈{0.5,1,2,3,4,5}.

• Support Vector Machine with Polynomial Kernel:

  • ∗ cost=1.

  • ∗ degree∈{2,3,4,5}.

For each of the three binary classifications, the relative influence of the predictor variables in the statistical model yielding the largest mean F1-score across the five folds is quantified using the methods described in Section 8.1 of “Greedy Function Approximation: A Gradient Boosting Machine” by Friedman [18].

3.Results/discussion

The mean and standard deviation across the five folds for the binary classifications of (i) Atrial Fibrillation vs. Non-Atrial Fibrillation, (ii) Arrhythmia vs. Normal Sinus Rhythm, and (iii) Arrhythmia with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia with the hyperparameters yielding the largest F1-Score are shown in Tables 2, 3, and 4, respectively. The results corresponding to the top-performing model with respect to each evaluation metric are displayed in bold. Observe that the Gradient Boosted Decision Tree Model outperforms all other models with respect to F1-Score and Accuracy across each of the three binary classification tasks, closely followed by the Random Forest Model. The maximum mean F1-Score attained by the Gradient Boosted Decision Tree Model across the five folds was 0.967, 0.839, and 0.943 for binary classification of Atrial Fibrillation vs. Non-Atrial Fibrillation, Arrhythmia vs. Normal Sinus Rhythm, and Arrhythmia with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia, respectively. The corresponding mean Accuracy attained by the Gradient Boosted Decision Tree Model across the five folds was 0.946, 0.946, and 0.921 for binary classification of Atrial Fibrillation vs. Non-Atrial Fibrillation, Arrhythmia vs. Normal Sinus Rhythm, and Arrhythmia with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia, respectively. The Gradient Boosted Decision Tree and Random Forest models outperformed all other models with respect to the area under the Receiver-Operator Characteristic Curves (AUC) for all three classification tasks as seen in Fig. 8, Fig. 9, and Fig. 10. This may be due to heterogeneity of the data; regardless, in computer-aided ECG analysis, interpretability of statistical models may be less important than the performance of said models, rendering more support in favor of ensemble and tree-based modeling approaches given their favorable performance.

Recall that TDA quantifies the ‘shape’ of data. Thus, the motivation behind presenting the classifications of both (i) Arrhythmia vs. Normal Sinus Rhythm and (ii) Arrhythmias with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia is to illustrate how the results are improved when TDA is used to classify two groups that primarily have different shapes, not frequencies. With this in mind, it may not be surprising that the presented TDA approach performs much better when classifying arrhythmias when the only two arrhythmias characterized solely by abnormal periodicity (assuming the individual has at most one rhythm as is the case in the data used in this study) – i.e. tachycardia and bradycardia – are not considered to be part of the arrhythmia group.

Fig. 8.

Receiver operator characteristic curve for classification of atrial fibrillation vs. Non-atrial fibrillation.

Fig. 9.

Receiver operator characteristic curve for classification of arrhythmia vs. Sinus rhythm.

Fig. 10.

Receiver operator characteristic curve for classification of arrhythmia with morphological changes vs. Sinus rhythm with bradycardia and tachycardia treated as non-arrhythmia.

The relative influence [18] of each predictor variable in the top-performing model with respect to mean F1-score (i.e. Gradient Boosted Decision Tree model) across the five folds in the classifications of Atrial Fibrillation vs. Non-Atrial Fibrillation, Arrhythmia vs. Normal Sinus Rhythm, and Arrhythmia with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia are shown in Fig. 11, Fig. 12, Fig. 13 in Appendix B. Atrial fibrillation is characterized by (1) absent/attenuated P-waves and (2) irregularly irregular frequency, so it is not surprising that the standard deviation of the RR-interval holds most influence for the classification of Atrial Fibrillation vs. Non-Atrial Fibrillation. Note that 9 of the 15 most influential predictor variables in the classification of Atrial Fibrillation vs. Non-Atrial Fibrillation stem from area-minimal cycle representatives and that 44115=38.3% of all predictor variables stem from area-minimal cycle representatives. In the case of Arrhythmia vs. Normal Sinus Rhythm, 8 of the 15 most influential predictor variables stem from area-minimal cycle representatives and 2876=36.8% of all predictor variables stem from area-minimal cycle representatives. Lastly, for the classification of Arrhythmia with Morphological Changes vs. Sinus Rhythm with Bradycardia and Tachycardia Treated as Non-Arrhythmia, 11 of the 15 most influential predictor variables are derived from area-minimal cycle representatives while 2466=36.4% of all predictor variables are derived from area-minimal cycle representatives. These results demonstrate that the predictor variables constructed from properties of the area-minimal cycle representatives are highly influential in the top-performing model regardless of the arrhythmia classification task.

The methods used in other studies that approach computer-aided ECG rhythm classification through a combination of TDA and machine learning are summarized in Table 5. Due to the wide range of classification tasks performed and evaluation metrics used in these studies, the classification tasks and evaluation metrics are not shown in Table 5 to avoid (i) presenting misleading comparisons and (ii) subjectivity in choosing the results from other studies to present. These other studies use a variety of databases [19,42] and sometimes a sample size on the scale of tens or hundreds, in addition to having longer – and consequently more informative – signals compared to the database used in this study [69]. Another factor to consider when comparing analyses of TDA and machine learning in ECG rhythm classification is the fact that different ECG databases often have signals labeled with different rhythms that may not be found in other ECG databases. The approach presented here attains similar results as these previous studies with respect to classification outcomes while utilizing a novel ECG signal processing pipeline and topological predictor variable construction, particularly with respect to using information derived from area-minimal cycle representatives.

4.Conclusion

The method presented here differs from other methods utilizing TDA and machine learning in three main ways:

There have been people working on computer-aided ECG analysis since the invention of the ECG machine. Over the past 20 years, there have been many machine learning approaches taken, yielding encouraging results. Some of these methods have involved TDA. Regardless of the type of method taken in computer-aided ECG analysis and the goodness of the evaluation metrics, we must take care to not rush to replace ECG interpretation by skilled health care professionals, however tempting the potential time and cost savings may be. In addition to the obvious danger of automated arrhythmia classification algorithms missing a harmful arrhythmia that a skilled healthcare professional would not have missed, bells and whistles from automated arrhythmia detection algorithms can lead to unnecessary medical staff fatigue and an increase in stress and adverse outcomes in hospitalized patients [12,29,33,54,58,59].

The data used in this study are free and publicly available at https://figshare.com/collections/ChapmanECG/4560497/2 [69]. The code used in this study is free and publicly available and can be found on GitHub: https://github.com/hdlugas/ekg_tda_arrhythmia_detection.