Functional feature extraction and validation from twelve-lead electrocardiograms to identify atrial fibrillation

Deep learning methods on standard, 12-lead electrocardiograms (ECG) have resulted in the ability to identify individuals at high-risk for the development of atrial fibrillation (AF)^1,2. These studies have the potential to aid with screening and risk stratification on a population level; however, the process remains a “black box” and does not help clinicians in understanding the electrocardiographic changes at an individual level. In particular, the findings on sinus rhythm ECGs from patients who eventually develop AF remain unknown. In this paper, we propose a data-driven approach to identify novel, subtle ECG features on sinus rhythm tracings that can provide insights into electrocardiographic abnormalities that precede the development of AF. In addition, subtle longitudinal changes of the sinus rhythm ECG features over time can provide insights into subclinical disease that precedes AF. These longitudinal changes can be surrogates for atrial arrhythmias and have the potential to be used for dynamic monitoring.

Deep learning methods have demonstrated great potential in disease detection and classification^{3,4,5,6,7,8,9,10,11,12}. An advantage of the deep learning methods is that they do not require feature extraction and work directly with the raw data. A known limitation of current machine-learning methods is that it is challenging to understand the rationale behind their results. Deep learning also requires tens of thousands of samples in the training dataset, which often stems from electronic health records. It is not feasible for these methods to leverage the strengths of traditional cohort studies that are designed to collect data serially across time and adjudicate rigorously clinical endpoints. In longitudinal settings, standard deep learning methods use data points as the basic unit of analysis to define predictors at different time points. Recurrent neural networks (RNN)^13,14 carry forward the predictors defined in the past and expand the number of predictors at each time. Consequently, the number of predictors can expand exponentially. Long Short Time Memory (LSTM)¹⁴ and its variants¹⁵ are proposed to reduce the number of predictors by ignoring the predictors from the distant history. Since the predictors are defined separately at different time points and are not directly comparable, it is not possible to derive the longitudinal changes of these predictors that are critical in understanding the subtle disease progression leading to the clinical onset.

In comparison to deep learning, feature extraction aims to create a reduced and meaningful representation of the raw data as a set of features, which can then be used or applied in subsequent analysis. Functional principal component analysis (fPCA)¹⁶ is a non-parametric data-driven approach for feature extraction in which a continuous curve is the basic unit of data analysis. The derived features are weighted averages of the entire curve. Standard principal component analysis (PCA) has been widely used in ECG signal processing^17,18. A few principal components can represent a majority of the variation in the original data, which can be used for the purpose of data compression^19,20,21,22, dimensionality reduction^23,24 and feature extraction^25,26. Clinical applications include classification and detection of arrhythmia²⁷, ischemia^28,29,30,31, long-QT syndrome³² and other T-wave morphology³³. PCA treats the measurements of ECG data at discrete time points as standard multivariate data, ignoring the temporal order of the data. In contrast, fPCA treats an ECG tracing as a continuous curve measured as discrete timepoints, taking into account its temporal order, which is particularly important when trying to have a clinical interpretation of the derived functional principal components. In addition, the features extracted using fPCA have intuitive interpretations, which is a significant advantage compared to deep learning methods. The extracted features can also be evaluated longitudinally at an individual level across serial ECG tracings. This process allows for explicit modeling of within-subject feature changes longitudinally and may represent a surrogate of subclinical disease progression. Through a longitudinal modeling, this approach is particularly powerful for detecting subtle changes over time, which can potentially lead to ahead-of-time prediction using the within-subject longitudinal trends.

In this manuscript, we apply the fPCA to raw ECG tracings collected in the Chronic Renal Insufficiency Cohort (CRIC)³⁴, an ongoing prospective study of individuals with chronic kidney disease (CKD). We evaluate the raw ECG tracings that were collected at the study baseline (2003–2008) to extract features, which are weighted averages of the ECG tracing of the entire heartbeat cycle, and screen those that are associated with AF. We also assess the change of selected ECG features over time and their prospective association with incident AF. Finally, we validate the selected ECG features in a separate group of participants who were enrolled in CRIC more than a decade after the initial phase (Phase III, 2013–2018).

In this study, we identify four ECG features related to the P-wave amplitude, QRS complex and ST segment. We show that both their baseline values and three-year changes are associated with the development of AF, which are validated using independent data from participants enrolled in Phase III of the CRIC study.

The CRIC study was reviewed and approved by the Institutional Review Board (IRB) at all participating clinical sites: University of Pennsylvania IRB, Johns Hopkins IRB, The University of Maryland Baltimore IRB, University Hospitals Cleveland Medical Center IRB, MetroHealth IRB, Cleveland Clinic Foundation IRB, University of Michigan Medical School IRB, Wayne State University IRB, University of Illinois at Chicago IRB, Tulane Human Research Protection Office IRB and Kaiser Permanente Northern California IRB. Each clinical center prepared a site-specific informed consent form following the guidelines of their local IRB and applicable regulations for informed consent. All participants provided written informed consent. The current study was approved to use deidentified data by the IRB at the University of Pennsylvania (IRB protocol # 849552).

The CRIC study was established in 2001 by the National Institute of Diabetes and Digestive and Kidney Diseases and is an ongoing cohort study of participants with CKD. In Phase I of the study, 3939 individuals were recruited through seven clinical centers nationwide between 2003 and 2008³⁴. A second wave of recruitment occurred during phase III of the study between 2013 and 2015 when an additional 1686 individuals with CKD were recruited. Study participants are followed through annual clinic visits, and a series of demographic and health-related data are collected according to study protocols, e.g., clinical history, kidney function etc. Cardiovascular endpoints including heart failure (HF), myocardial infarction (MI), stroke, and AF underwent a rigorous adjudication process that has been detailed previously³⁴. In this manuscript, we focus on identifying ECG features that are important markers of AF, the most common arrhythmia affecting CKD patients.

Our study cohort was restricted to those individuals who had no history of any cardiovascular disease at the time of recruitment into the cohort. As such, 2623 participants who were recruited during Phase I, and 1116 participants, who were recruited during phase III, comprised our study population (Supplemental Fig. 1).

Electrocardiogram

All twelve-lead ECGs were recorded annually as part of routine clinical study visits using standard procedures³⁵ and identical equipment (MAC 1200; GE Medical Systems Information Technologies, Milwaukee, WI). Digitally recorded ECGs were transmitted to the ECG reading center located at Wake Forest University. After being visually checked for quality, the ECGs were automatically processed using the 2001 version of the GE Marquette 12-SL program. The ECG tracing was digitized at 500 Hz for ten seconds that includes approximately 10 heartbeat cycles, from which the median across cycles within-subject is calculated for each lead, and is referred to as the raw tracing throughout the manuscript.

Each ECG raw tracing included 600 data points, which represents a 1.2 s recording at a frequency of 500 Hz. The raw ECG tracing was pre-processed in a few steps. First, we subtracted each raw ECG tracing by its mean value between 0.1 and 0.2 s prior to the onset of P-wave, so that all tracings have an initial value of zero. We then took the square root of the raw values while keeping the sign of the original value and smoothed the raw tracing using B-splines. The final data included the recordings between 0.2 to 1 s of the raw tracing after excluding the flat region that corresponds to either the start or end of a cardiac cycle.

Feature extraction and selection

Our initial focus was to extract ECG features associated with the development of AF in participants recruited during Phase I of the CRIC study. To define the candidate features, we applied fPCA to the raw ECG tracing using the R package fda³⁶ for each lead separately. We selected the first eight principal components (referred to as ECG features) from each lead that explained most of the variation (i.e., >90%) in the raw tracing. An orthogonal transformation, i.e., varimax rotation, was then applied to the eight principal components from each lead to improve the interpretability. Consequently, a total of 96 candidate ECG features were extracted from the raw 12-lead ECG tracing. Their values are the corresponding fPCA scores, calculated as the weighted averages of the ECG tracing. The weight functions for calculating the fPCA scores are then fixed in all subsequent analyses.

We compared the candidate features of normal sinus rhythm ECG tracings between individuals who developed AF and those who remained free of any cardiovascular event (including incident HF, MI, stroke, or AF) or death during the study follow-up. To maximize differences of the cross-sectional comparisons in defining ECG patterns, we selected the sinus rhythm ECG that was most proximal within a 2-year timeframe prior to the development of AF (N = 294) (Supplementary Fig. 1A). We also evaluated baseline ECGs from those individuals who did not develop any cardiovascular event during follow-up (N = 1333) (Supplementary Fig. 1A). A total of 1627 ECGs were used for the cross-sectional feature extraction and selection. Since the 12 leads are correlated, the features extracted from different leads may represent similar underlying signals. An additional selection procedure was used to select the most predictive subset of features. We compared the derived ECG features between those who developed AF and those who did not using Student t test and selected features that had a p-value less than 0.0005. This allows for correction of an overall type one error rate less than 0.05 with 100 simultaneous comparisons. Forty-seven features were selected from this univariate analysis. Logistic regression with a stepwise selection strategy was then used to select the final list of ECG features with a p-value threshold of 0.001 for both entering and remaining in the model. To visualize the selected features, we plotted the average ECG tracings that are representative of each feature between the AF and non-AF individuals. The average ECG tracings were calculated as the weight function multiplied by the corresponding average principal component scores for each group, adding the population average ECG tracing. Additionally, we examined the association of selected ECG features with AF through a multivariable logistic regression model that included all selected features. The odds ratio (OR) (and 95% confidence interval) of AF occurrence associated with each feature was reported as per one standard deviation change in the corresponding feature value.

To examine if the selected ECG features appear similarly across all twelve leads, for each selected feature, we applied the same weight function estimated from the cross-sectional analysis to all twelve leads. The selected features were then visually compared between the AF and non-AF individuals in all twelve leads.

Changes of the selected features and their associations with AF

Our other objective was to evaluate if longitudinal changes in the final set of selected ECG features can be a precursor of AF. As an example, Fig. 1 provides the ECG tracings of two individuals, each with two repeated ECG tracings. One individual developed AF during the study follow-up and the other one did not. In our analysis, we evaluated prospectively the changes in the selected final set of ECG features between the baseline and three-year followup. Among those who were recruited in Phase I of the CRIC study, 2357 participants had two or more ECGs collected within a three-year window from study entry (Supplementary Fig. 1B). We used the first and last ECGs to calculate changes of selected features, which were then evaluated as risk factors for incident AF. We applied the same weight function estimated from the cross-sectional comparison to both the baseline ECG and those collected 3 years later to compute the corresponding fPCA scores. The difference of the two fPCA scores within each individual was then computed.

a An individual who developed atrial fibrillation (AF) during the study follow-up. b An individual who remained healthy throughout the study follow-up. The black and red lines represent the ECG tracing collected at the study entry and follow-up, respectively.

Full size image

Similar to the previous feature visualization, we reviewed graphical depictions of the selected features between the two, longitudinal ECG recordings stratified by AF incidence status. To study the risk of the ECG feature change and AF, we fit multivariable Cox proportional hazards model. The model included both baseline ECG features and their 3-year changes as covariates. In the model, participants were censored at the end of study follow-up, date of withdrawal or death, whichever occurred first. The hazards ratio (and 95% confidence interval) was reported as per one standard deviation change in the respective covariate.

We were also interested in examining if the selected ECG features and their changes over time can be potential surrogates for subclinical disease progression. We evaluated a group of representative baseline demographic and clinical risk factors (age, race, sex, smoking status, body mass index (BMI), systolic and diastolic blood pressure, diabetes, estimated glomerular filtration rate (eGFR) and urine protein to creatinine ratio (uPCR)) that have previously demonstrated independent associations with AF. We assessed whether the association between these demographics and clinical risk factors and AF can be explained by the selected ECG features and their longitudinal changes. The analysis was done using R package mediation³⁷ with modifications to allow for both multiple exposure (i.e., baseline clinical risk factors) and intermediate variables (i.e., selected ECG features and their change over time). We fit two separate parametric survival models for time to incident AF with Weibull distribution that included baseline demographics/clinical risk factors and ECG features respectively. We then estimated the linear combination of baseline demographics/clinical risk factors and ECG features from the two fitted models. The mediation package was used to estimate the percent of variation between the linear combination of baseline demographics/clinical risk factors and AF that can be explained by the linear combination of the selected ECG features.

Validation of selected ECG features using participants recruited one decade later (Phase III of the CRIC Study)

Since the raw ECG tracings from phase 3 participants were not used for model development, we validated our phase 1 findings in those recruited more than a decade later. In Phase III of the CRIC study, ECGs were collected in 1077 participants, and 909 participants had at least two ECGs collected within three years from study entry (Supplementary Fig. 1C). We assessed fPCA scores at baseline and 3-years follow-up. Similar to the phase 1 analysis, a multivariable Cox proportional hazards model for time to AF was fit that included the selected ECG features at baseline and their 3-year changes. Finally, we examined how much the association of baseline clinical risk factors with AF is explained by selected ECG features evaluated through both recordings in Phase III participants.

Description of study participants

There were 3939 individuals recruited during Phase I of the CRIC study and had a median follow-up period of 7.6 years [Interquartile range (IQR): 4.3–13.0 years]. 1627 participants did not have any history of cardiovascular disease and were used for ECG feature extraction and selection. Of these, 294 individuals developed AF and the remaining 1333 individuals were free of AF and other cardiovascular events (Supplementary Fig. 1A). The average age at the time of enrollment was 55 (standard deviation (SD): 12 years and 49% were female (Supplementary Table 1). In addition, 50% were white, 34% were Black and 12% were Hispanic. Diabetes was present in 33% of individuals. The mean eGFR was 50 mL/min/1.73 m² (SD: 18) and median urine protein to creatinine ratio (uPCR) was 0.097 (interquartile range (IQR): 0.048–0.45)).

As part of the longitudinal assessment that evaluated the changes in selected ECG features and risk of AF, we evaluated 2357 participants from CRIC Phase I (Supplementary Fig. 1B). For this analysis, we expanded our study population to include those who may have developed other cardiovascular events as we were prospectively evaluating the selected features and incident AF risk. The characteristics of these participants were similar to the more restricted phase 1 population that was used for the feature extraction and selection analysis (Supplementary Table 1). In addition, there were 909 participants recruited during CRIC Phase III and included in the external validation analysis for the selected ECG features (Supplementary Fig. 1C). Compared to those recruited in Phase I, these participants were older (mean (SD): 64.3 (8.3). There were also more Black participants (49.0%). The mean eGFR at baseline was higher (mean (SD): 56.7 mL/min/1.73 m² (SD: 13.7)) and uPCR was similar (median: 0.12; IQR: 0.055–0.41).

ECG feature extraction and selection using CRIC phase I data

After applying functional PCA to all 12 leads, there were 8 ECG features extracted from each lead. These accounted for more than 90% variation in the raw tracings for each lead. After univariate screening and stepwise regression, four features were selected representing components of the P-wave, QRS complex and ST segment from leads II, aVF, V2 and aVL, respectively (Fig. 2). Those who developed AF had lower P-wave amplitude, lower QRS amplitude and lower slope in the ST segment, compared to the non-AF individuals.

The features are a the P wave (II.5), b QRS complex (aVF.2), c ST segment (V2.6), and d ST segment (aVL.6). The red and blue lines correspond to those who developed AF and healthy individuals respectively. The lines are calculated by adding to the mean ECG curve the product of the principal component weight functions (displayed in the top panel) with the corresponding average principal component scores for each group.

Full size image

For the four selected features, we evaluated the corresponding features in the other eleven leads (Supplementary Fig. 2 through 5). Differences in the P-wave related feature were more pronounced in leads I, II, III, aVR and aVF compared with the other leads (Supplementary Fig. 2). The other three selected ECG features also demonstrated greater differences between AF and non-AF cases in specific leads (Supplementary Fig. 3–5).

The association between each ECG feature on the initial ECG tracing was associated with AF after adjusted for the other features [OR (per one standard deviation change) for features related to P wave amplitude 1.72 (95% CI: [1.51, 1.97]), QRS complex 1.59 (95% CI: [1.38, 1.85]), and two separate regions of the ST segment 1.39 (95% CI: [1.21, 1.60] and 1.87 (95% CI: [1.61, 2.17])] (Table 1).

We evaluated the correlation between the four selected features and other existing P wave measures from GE Muse program that include P wave amplitude in Lead II, P wave duration in Lead II, P prime wave amplitude in Lead V1, PR interval and the left and right atrial enlargement using the NOVA code. The largest correlations were observed between our selected P wave amplitude (II.5) feature and P amplitude and duration in Lead II (0.47 and 0.38 respectively). The correlations between the other selected ECG features and the MUSE variables are relatively small (Supplementary Table 2). In addition, the standardized mean difference of the four selected features between those who had atrial enlargement (either left or right) and those who did not range between 0.16 (ST segment (V2.6) feature) and 0.30 (P wave amplitude (II.5) feature), all of which were relatively small. When including these five measures in the multivariable logistic regression, the OR for the P wave amplitude (II.5) feature was attenuated (OR: 1.34; 95% CI: [1.09, 1.63]). The ORs for the other selected ECG features remain largely unchanged (Supplementary Table 3).

Prospective evaluation of selected ECG features and their changes on AF risk in CRIC phase I data

Individuals who developed AF had both a lower P-wave amplitude at baseline and a greater reduction in the P-wave amplitude over time compared to those that did not develop AF (Fig. 3). Similarly, there was both a lower QRS amplitude and ST-segment slope in baseline ECGs and a greater reduction in these parameters among individuals who developed AF compared to those who did not. When evaluating the model that combined both baseline and longitudinal changes in ECG parameters, the 4 ECG features assessed at baseline and changes in the P-wave and ST-segment were associated with incident AF risk after adjustment for the other features. Changes over time in the QRS amplitude and the other ST segment feature were not independently associated with the development of AF (Table 2, left column). Further, 12.8% (95% CI: [10.1%, 16.0%]) of the association between baseline clinical risk factors and AF was explained by the selected ECG features and their changes.

The features are a the P wave (II.5), b QRS complex (aVF.2), c ST segment (V2.6), and d ST segment (aVL.6). The ECGs measured at baseline and follow-up are represented by dotted and solid lines respectively. The red and blue lines correspond to those who developed AF and healthy individuals respectively.

Full size image

When including the five Muse variables and their changes from baseline to follow up visit in the multivariable Cox regression, the HR for the baseline P wave amplitude (II.5) feature became nonsignificant (Supplementary Table 4). The association for the other three features at baseline remained statistically significant (Supplementary Table 4). Additionally, the HR for the change of P wave amplitude (II.5) feature remains largely unchanged (HR: 1.22 [1.08, 1.38]).

External validation of ECG features in CRIC phase III data

Similar findings were observed when validating the findings in the phase 3 data. In particular, 3 of the 4 ECG features that were evaluated in baseline tracings among phase 3 participants and the change in the P-wave feature were independent risk markers for incident AF risk. (Table 2, right column). Both the baseline P-wave feature and its change over time had the greatest estimates of risk: (HR for baseline P-wave feature: 1.74 95% CI: [1.42, 2.13]; HR for change of P-wave feature: 1.47, 95% CI: [1.18, 1.83]). Similar to the participants enrolled in phase I, ECG features and their changes explained 14.5% (95% CI: [7.9%, 22.0%]) of the association between clinical risk factors and incident AF.

Using functional principal components analysis, we extracted baseline ECG features in sinus rhythm that were associated with incident AF, independent of well-established P-wave measures that are known markers of atrial health. One of these features, P-wave amplitude, was a dynamic marker of risk as reduction in its magnitude over a 3-year interval of serial ECG recordings was associated with incident AF. In particular, one standard deviation of the 3-year decline in P-wave amplitude was associated with a 29% increased risk of AF.

Our agnostic findings highlight the significance of the P-wave amplitude and complement prior studies that demonstrated associations between select P-wave abnormalities and AF^{38,39,40,41,42}. By pursuing a longitudinal analysis on annual ECGs over a 3-year duration, our analyses also suggest that progressive atrial disease, which is evidenced by a larger reduction in atrial voltage than what is seen with time alone, is present among individuals who develop AF. When evaluated sequentially over time, this feature could serve as an AF surrogate and could alert at-risk patients for potential stroke prevention therapies such as anitcoagulation. In our analysis, we found that the four selected ECG features and their changes within three-years accounted for more than 12% of the association between a set of clinical risk factors (including demographic, comorbidity and kidney function measures) and AF incidence in both Phase I and III of the CRIC data. Expanding the longitudinal assessment beyond 2 ECGs may help to increase the signal to noise ratio of selected features and improve their ability to serve as risk markers. This work will be examined in a planned follow-up study.

Deep learning methods have demonstrated remarkable efficacy in detecting arrhythmia and AF^12,43 by analyzing ECG tracings. These studies have evaluated thousands of raw ECG tracings at a single time point and developed classification schemes to identify individuals at the greatest risk of developing AF. However, these prior studies provide limited insights into the specific measures of the tracing that are responsible for identifying risk. In contrast, our proposed method extracts features through functional principal component analysis, which enables straightforward interpretation by highlighting specific waveforms or components of the ECG tracing that are associated with the outcome of interest. These features are easily visualized by examining the corresponding weight functions. In addition, the fundamental unit of analysis in our proposed method is the functional principal component estimated using cross-sectional data, which remains constant when computing the longitudinal change of corresponding fPCA scores. This allows us to model subtle longitudinal changes of features over time. We identified the P-wave amplitude as the most significant feature associated with incident AF. We were also able to evaluate systematically its change in serial ECGs over years as a marker of AF risk. Although the average reduction of p-wave amplitude over three-years is only 2.3% in those recruited in Phase I (1.2% and 8.0% in the non-AF and AF individuals respectively), we were still able to demonstrate a significant association between the P-wave amplitude change and prospective risk of AF. We also identified subtle changes in ventricular repolarization that involved the slope of the ST-segment and delay of T-wave onset. These repolarization changes may represent abnormalities in myocardial health in the CKD population and be markers of AF risk. Further, the likelihood that these repolarization measures represent subclinical disease is likely limited, especially when compared to P-wave amplitude, since they did not prove to be independent, dynamic markers of AF risk in our external validation subset. Overall, our findings suggest that subclinical changes in cardiac electrophysiology can potentially be used for ahead of time prediction of AF. Although deep learning methods like LSTM can be used to model longitudinal data by tracking the up-to-date information on selected features, its ability to model the feature change over time is limited. Consequently, subtle changes like what we observed in the P-wave amplitude may not be identified using deep learning methods.

Our proposed fPCA method of ECG feature extraction from the raw tracing belongs to the large category of unsupervised feature extraction. It is a nonparametric data adaptive approach for feature extraction that assumes no structure on the raw ECG tracing and does not require the labeling of different phases of a cardiac cycle, an advantage compared to other methods, e.g., mathematical modeling, that require strong and rigid parametric assumptions⁴⁴. In addition, fPCA is very effective in reducing the dimensionality of the raw ECG tracing. For example, in the analysis of the CRIC ECG data, it requires only eight features from each lead to account for more than 90% of the total variation observed in the raw ECG tracing.

There are a few limitations in this study. First, in the analysis of ECG raw tracing, we did not account for phase variability caused by different heart beat cycles among individuals. As a result, the extracted features may not accurately represent the raw tracing. To address this limitation, future research could explore the use of curve registration techniques to explicitly remove the phase variation in ECG raw tracings before feature extraction. Second, the evaluation of longitudinal changes in ECG features was based on only two ECGs from the same individual. Including a larger number of serial ECGs can provide a more accurate characterization of feature changes over time and strengthen the association with AF. Third, the feature is extracted one lead at a time, without the consideration of correlation across different leads. Future work may consider multivariate fPCA^45,46,47 to simultaneously analyze all 12-lead tracings. It will be interesting to see if additional insight can be obtained by explicitly exploring the between lead correlations. Last, although we validated the ECG features and their changes using an independent set of individuals with CKD recruited during a later phase of the CRIC study, it may be beneficial to conduct external validation in real clinical settings, e.g., using electronic health records.

In the external validation, we will apply the same weight function for each identified feature from the current study to ECG raw tracings collected from an external data source to compute features and compare them between those who developed AF and those who did not. Furthermore, we will assess the longitudinal change of each feature and examine if it is associated with the incident risk of AF in a prospective setting. If the selected features are significant in the external cohort, we would consider it as a success validation, although the magnitudes can be slightly different. Regarding the selection of an external data source, we will first restrict the study population comparable to that of the CRIC study, e.g., including only patients with CKD and no history of cardiovascular disease. We will consider the adjustment of baseline covariates to make the external cohort as comparable to the CRIC cohort as possible. If these features are validated in the CKD population, we will then expand it to a more general population.

In conclusion, we have proposed a framework for data driven feature extraction using raw ECG tracings, and we discovered P-wave amplitude and its longitudinal change as risk markers for the development of AF. The extracted features have intuitive interpretations. The proposed method is straightforward and applicable to other cardiovascular disease endpoints.