Using machine learning to identify local cellular properties that support re-entrant activation in patient-specific models of atrial fibrillation

Abstract Aims Atrial fibrillation (AF) is sustained by re-entrant activation patterns. Ablation strategies have been proposed that target regions of tissue that may support re-entrant activation patterns. We aimed to characterize the tissue properties associated with regions that tether re-entrant activation patterns in a validated virtual patient cohort. Methods and results Atrial fibrillation patient-specific models (seven paroxysmal and three persistent) were generated and validated against local activation time (LAT) measurements during an S1–S2 pacing protocol from the coronary sinus and high right atrium, respectively. Atrial models were stimulated with burst pacing from three locations in the proximity of each pulmonary vein to initiate re-entrant activation patterns. Five atria exhibited sustained activation patterns for at least 80 s. Models with short maximum action potential durations (APDs) were associated with sustained activation. Phase singularities were mapped across the atria sustained activation patterns. Regions with a low maximum conduction velocity (CV) were associated with tethering of phase singularities. A support vector machine (SVM) was trained on maximum local conduction velocity and action potential duration to identify regions that tether phase singularities. The SVM identified regions of tissue that could support tethering with 91% accuracy. This accuracy increased to 95% when the SVM was also trained on surface area. Conclusion In a virtual patient cohort, local tissue properties, that can be measured (CV) or estimated (APD; using effective refractory period as a surrogate) clinically, identified regions of tissue that tether phase singularities. Combing CV and APD with atrial surface area further improved the accuracy in identifying regions that tether phase singularities.


Introduction
Atrial fibrillation (AF) is a supraventricular tachyarrhythmia, characterized by unco-ordinated activation of the atria, 1 deterioration of mechanical function, 2 and increased incidence of cardiovascular disease, stroke, and pre-mature death. 3 In drug-refractory patients, AF is commonly treated by radiofrequency catheter ablation (RFCA). 4 Central to RFCA for AF is pulmonary vein isolation (PVI), where the tissue surrounding the pulmonary veins that is prone to triggered spontaneous electrical activation is electrically isolated. Pooled single-procedure 12-month arrhythmia-free survival following PVI only in paroxysmal and in persistent AF is 66.6% and 51.9%, 5 respectively. Ablation targets beyond PVI are likely needed to improve RFCA success rates for AF patients. To improve success rates, additional lesion sets that aim to remove specific regions of tissue or substrate that support AF have been proposed. An early, and initially promising, approach for identifying AF sustaining tissue was targeting complex fractionated atrial electrograms (CFAE). 6 However, further larger studies found that the addition of CFAE ablation did not improve patient outcomes. 7 Subsequent re-entrant activation substrate ablation strategies have targeted focal or re-entrant activation patterns identified with basket catheters 8 or high frequency activation regions, 9 however, these approaches have yet to produce consistent results.
Many of the proposed ablation targets require high spatial and temporal resolution activation maps. While high-density catheter mapping systems are available, they only map from one region of the atria at a time. Catheters that map the whole atria are lower density or rely on non-contact mapping, both approaches can lead to artefacts. 10 This can make it hard to reliably identify ablation targets. However, implicit in strategies that target tissue that support reentrant activation, is that a combination of material properties and/or location leads to a specific region of tissue playing a critical role in sustaining the arrhythmia. We hypothesise that if re-entrant supporting tissue has specific characteristics, then these can be identified by local mapping of tissue properties.
Biophysical computer simulations of the atria provide a computational framework for linking cellular properties to emergent tissue scale properties and provide high spatial and temporal resolution signals across the entire atria for identifying re-entrant activation patterns and identifying the regions of tissue that support these activation patterns. Workflows for creating and validating models of the atria have been proposed based on CMR images or electro-anatomical maps. 11,12 These personalized models of the atria have developed to a point where they are being used to guide therapies 13 and virtual patient cohorts have been used to investigate patient-specific lesions, response rates, and experimental ablation procedures. 14,15 Virtual cohorts provide a testing environment for studying the local tissue properties that support re-entrant activation patterns. We aim to create a validated virtual cohort of patientspecific atrial models that capture patient-specific electrical heterogeneity to simulate AF and then use machine learning to identify the tissue properties associated with sustained reentrant arrhythmias in the left atria.

Methods
In this article, we generate personalized computer models of the left atrium of a cohort of 10 patients. The model was fitted to activation measurements recorded across the atria at up to 16 sites, as described by Corrado et al. 12 To trigger AF, we applied a burst pacing protocol on 12 different locations in the proximity of the pulmonary veins. On simulations that presented a self-sustaining (SS) AF, we evaluate the phase (PS) singularity map to classify the tissue capable of tether PS. Finally, we trained a SVM classifier on maximum local conduction velocity (CV max ), and action potential duration (APD max ) to identify portions of tissue that tether PS.
We provided a detailed description of the methods in Supplementary material online.

Virtual patient cohort
We recruited 40 patients; 10 patients with AF had sufficient evenly spaced recording points in the left atria to allow us to create a model. These included seven paroxysmal AF (PAF) and three persistent AF (PsAF) cases.

Creating and validating the model
In this article, we generate personalized computer models of the left atrium of a cohort of 10 patients. The model was fitted to activation measurements recorded across the atria at 42-100 recording sites per case, as described by Corrado et al. 12 Patient specific, spatially varying model parameters were fitted using an S1-S2 pacing protocol applied through the coronary sinus. Recordings were made at 14 ± 3 sites across the atria. Two model parameters (the conductivity and s close ) were analytically mapped to CV max and APD max allowing these parameters to be fitted directly to functional measurements. The remaining cellular depolarization and repolarization parameters (h min , s open , s in ) were then fit to measured local activation times and the effective refractory period. 12 The predicted activation times across the whole atria, when using locally fitted material properties, are compared against measured activation times used to fit the model (Figure 1, top). The model was validated by comparing predicted with clinically measured activation times during an S1-S2 pacing protocol applied at the high right atrium (Figure 1, bottom). The validation data were not used in training the models. The distribution of the fitted parameters for each patient for each parameter is presented in Supplementary material online, Figures S1-S5.

Induction of atrial fibrillation and atrial tachycardia
To systematically evaluate the capacity of each atria to sustain AF, we applied a burst pacing protocol from one of 12 sites, 3 around each pulmonary vein (Figure 2A), and measured the type (AF, meandering rotor, or macro-re-entrant), frequency, and duration of sustained activation. Sixty-six percent of simulations terminated within 1 s and were classified as non-triggering, 8% terminated between 1 and 3 s, 9% terminated between 3 and 40 s, and 17% were classified as sustained activation and remained in activation after 80 s ( Figure 2B). About 7.5%, 4.2%, and 5% of simulations that led to sustained AF with a frequency of 6.05 ± 0.98 Hz, meandering rotor with a frequency of 3.89 ± 0.23 Hz, or stable-re-entrant activation with a frequency of 4.77 ± 0.4 Hz, respectively.

What's new?
• Left atrial size and local measurements of conduction velocity (CV) and effective refractory period (ERP) identify tissue that can sustain phase singularities.
• An support vector machine classifier can be used to identify tissue that can sustain phase singularities.
• Larger atrial size and local short action potential duration and slow CV correlates with re-entrant activation patter tethering.
Machine learning to identify tissue supporting AF

Transient changes in activation dynamics
Self-sustained activation in the model is initiated from a model of sinus rhythm. There is likely to be a transitory period, where initial conditions are still important, at the start of the simulation. Similarly, there may be a change in frequency prior to termination of SS activation patterns. We measured the frequency of activation in 500 ms windows. In non-terminating simulations, the activation pattern reached the activation frequency achieved in the final 5 s of simulation after a median of 303-341 ms in AF, meandering rotor and macro re-entry cases. We saw no change in the activation frequency in the window when activation terminated and the prior window. The temporal evolution of the activation frequency is plotted in Supplementary material online, Figures S11-S14, for each patient and for each pacing site. For our analysis, we consider activation times from 500 ms through to, whichever comes first of, termination or 80 s.

Importance of pacing site on duration of self-sustaining activation
Pacing from different veins caused different durations of activation patterns, with pacing from the left superior PV leading to twice as many sustained activations compared with the right infer pulmonary vein ( Figure 2C). The ability of an atria model to sustain AF cannot be determined by simulating burst pacing to induce AF from a single site. We consider all activation patterns initiated from all the pacing sites in our analysis.

Distributions of phase singularities
To test if the type of SS activation pattern supported in each atrium model was consistent, regardless of pacing location, we plotted the phase singularity maps for each SS activation pattern. This identifies if regions of tissue that tether re-entrant activation patterns were consistent, regardless of pacing site. We plot the unwrapped phase singularity density maps, starting from 75 s after pacing to remove potential transitory effects of moving rotors, for cases with sustained activation patterns ( Figure 3). In some cases (6 and 9), the PS map is similar for different pacing sites, while for other cases (7 and 8), the PS maps are different for different pacing sites. Supplementary material online, Figure 15 plots the combined phase singularity maps for each clinical case. There is not a unique PS map for a given atrial anatomy and set of material properties. We can only determine if a region has the potential to support a PS, not that it will support PS.
Are specific parameter sets associated with self-sustained activation?
We classified each case as self-sustaining if at least 1 of the 12 simulations produced a self-sustained activation pattern for the full 80 s and non-triggering or self-terminating (NTþST) otherwise. We compared the distribution of fitted parameters between patients (Supplementary material online, Figure S10) and when grouped across patients ( Figure 4A) for these two outcomes. CV max , s in , h min , and s open have considerable overlap in parameter distribution while APD max is shorter in atria with SS activation.
Are specific parameter sets associated with phase singularities?
We classified as 'tethering' the portions of the tissue presenting a PS density maps at 1SD above the mean of the PS 16 and as 'normal' otherwise. Figure 4B plots the value of the parameters that characterize the substrate that support PS (red) compared with other tissue. The CV max value shows a separation between the two classes (median values: 60 cm/s on the substrate with PS; 160 cm/s in the rest of the tissue). Slower values of CV max are associated with regions that anchor PS. Machine learning to identify tissue supporting AF i15 How do APD max and CV max vary in each patient? Figure 4C plots the distribution of CV max and APD max for each clinical case; we obtained each distribution by fitting a Gaussian Kernel and choosing the empirical standard deviation as the bandwidth. We see the high distribution of APD max in four of the five cases that exhibit SS activation patterns. The majority of cases has a comparable distribution of CV max regardless of their ability of the case to support a SS activation pattern. The distribution of cell model properties on all the 10 clinical cases are provided in Supplementary material online.

Do tissue properties identify PS location?
We have demonstrated that APD max and CV max are correlated with SS activation and determining the location of PS. To test how well these two material properties support PS tethering, we trained a support vector machine (SVM) on the CV max and APD max values at every point. The SVM was trained on material properties at 3 568 549 data points, with 21 391 classified as tethering and 3 547 158 as normal.
The large data set (3 568 549 data points) and the complex correlations between the local input features (CV max and APD max ) and global features (atrium size) that produce a PS are not readily analysed using classical statistics.  The best classifier achieves an accuracy of 91% (sensitivity = 0.96, specificity = 0.91) indicating that APD max and CV max play an important role in the location of PS tethering.
Does atrial size impact the importance of CV max and APD max on PS Location?
We have previously shown that effective size, which combines an APD surrogate (ERP), CV, and area values, predicts response to PVI in persistent AF cases. 17 To test if atrial size impacts the regions that can support PS, we retrained the SVM including each patients atrial size as an input feature. We added atrial size into the CV max and APD max SVMs. When atrial surface area CV max and APD max are taken into account, the classifier achieves an improved accuracy of 95% (sensitivity = 0.99, specificity = 0.95).

Discussion
We have created a validated virtual cohort of 10 (7 paroxysmal and 3 persistent) AF patients with heterogeneous material properties fitted to each individual patient's electrophysiology measurements.
We show that the capacity to initiate an arrhythmia is dependent on the pacing site, after a transitory phase the activation rates remain stable, termination can still occur up to 40 s after initiation, PS maps are not unique for a given atrium anatomy and material properties, shorter APD max is associated with self-sustained activations, the CV max is associated with locations that support PS, that a combination of APD max and CV max can identify PS tethering locations with 91% accuracy and that accounting for atria size improves the identification of PS tethering locations to 95%.

Initiating and analysing simulated atrial fibrillation
To reduce computational cost, atrial arrhythmia simulations have also adopted relatively short simulation durations, often as low as 5 s 13 for identifying sustained arrhythmia. Here we show that 9% of self-terminations occur after 3 s and self-terminations can still occur even 40 s after initiation, motivating the use of longer simulations for clinically meaningful durations. However, 80 s is a long time to perform a simulation. For the 120 simulations that we performed, only one case would have been misclassified as self-sustained if we had adopted a 30 s simulation duration, in line with the duration used for defining an episode of AF, 18 and this may represent a reasonable balance of clinical relevance and computational tractability.
We considered that the activation patterns may organise or slow down prior to termination. However, consistent with clinical studies, 19,20 we found that there was no change in frequency in the models prior to terminations. We did, however, observe a transitory period following initialization of AF of 330 ms. This may be an underestimate of the transition period in clinical studies and simulations studies using complex models, where slower changes in ionic homoeostasis, which are not accounted for in the mMS model, may cause a longer transitory period.

Re-entrant and ectopically sustained activation
Consistent with clinical observations that ectopic activations start from the left superior pulmonary veins, 21 we find that pacing from the left superior pulmonary veins leads to more sustained activation patterns. In our models, only 50% of the virtual cohort cases could support a self-sustained activation pattern. We did not account for sustained ectopic events following the initial burst pacing in the models, as we aimed to only simulate the effect of the electrophysiological substrate to support re-entrant activation patterns. In many cases, particularly in paroxysmal atria, AF is initiated and potentially driven by ectopic beats, often originating in the ostium of the pulmonary veins. 21 The inability of our atrium models to sustain AF may represent cases where their clinical AF was supported by ectopic beats and not an arrhythmogenic substrate. These patients may represent cases where PVI is more likely to be effective.
At the same time, the ability to sustain AF may be highly dependent on attributes not included in the models such as endocardial and epicardial dissociation, fibrosis, or the right atrium. We did not have access to AF induction studies in these cases, the location of ablation lesions delivered, or ablation outcome. The addition of both more sophisticated models and additional therapy and outcome data from patients would allow us to test what factors contribute to simulated sustained activation and if cases that cannot sustain activation have AF driven by pulmonary vein ectopic beats.

Local tissue properties that support phase singularities
Previously re-entrant activation patterns have been associated local fibrotic regions measured using MRI in simulation 22,23 and clinical 24 studies, however, others found no link between re-entrant activation and fibrosis. 25 Alternately, other groups have suggested regions with low voltage and CV, 26 or regions with CFAE 27 could support reentrant activation. In this paper, we assumed the effects of fibrosis would be captured in electrical measurements 28 and consequently be encoded in the fitted cell models. The classifier identified tissue associated with tethered PS with 91% accuracy and comparable F1 score. Experimental noise may impact the SVM classifier. As we trained and tested the SVM classifier using model parameters fitted to experimental data that were subject to noise, it is not expected that the accuracy of the SVM will degrade when applied to measurements with a similar level of noise. To estimate how increased noise may impact the classifier, we added a 10% white Gaussian noise into testing data set. This elevated level of noise, beyond what is seen clinically caused the accuracy of the classifier to decrease from 95% to 93%.
We also showed how patient attributes beyond tissue properties can be used to improve identification of PS tethering locations. The addition of patient atrial size further increased the ability of the SVM to identify PS tethering sites to 95% accuracy. This is consistent with many previous imaging studies that have identified atrial size or shape as potential contributors; these are used as indexes of risk, although not as a factor in guiding therapy.

Model critique
We have created a model that reflects the complexity of the available clinical data, we have aimed to capture the data that are available in the catheter lab. This represents a trade-off between model personalization and sophistication. In our simulations, we have not attempted to model fibre distribution, biophysical cellular physiology, endo-epicardial dissociation, fibrosis, the right atria, or wall thickness. We have assumed that the effects of tissue change due to fibrosis, remodelling, or physiological heterogeneity will be captured by the measured activation times and will be encoded within the spatially varying cell model parameter values. We have not included the right atria as the source of AF is predominantly found in the left atria and the right atria is not routinely mapped during ablation procedures. While we recognise that many of the attributes we have not included in the model may play a role in AF, we were able to train and validate our model on an independent data set and did not find any evidence that we needed to include these extra attributes to replicate the clinical observations. This approach represents an alternate strategy for patient-specific atrial models from purely image-based approaches, where cellular properties are assumed common between patients and are based on MRI image intensity values.
We created and validated our model using a systematic S1-S2 pacing protocol. We did not validate the model against patient-specific electrical recordings of AF data. We do not have AF recordings for these patients. Validating patient-specific simulations of AF is challenging. Boyle et al. 13 have indirectly validated their model by showing that in 10 cases the model can identify critical ablation sites. Comparison of simulated PS maps with clinical measurements is challenging due to disagreement on the number of PS locations that should be identified, 9,29 recording devices have limited coverage and the use of open and proprietary PS identification tools. We do not see one standard PS distribution in our simulations. We see high densities of PS, consistent with 29 (Case 2), meandering rotors, consistent with 30 (Cases 7 and 9), and stable rotors consistent with 8 (Cases 6 and 8). In this paper, we presented a limited number of cases as each personalized model requires a high-quality data set. Increasing the number of cases modelled, in particular for persistent AF or cases with a sustained AF, may provide a promising approach for identifying the mechanisms underpinning these different re-entrant activation patterns.

Limitations
We modelled the atria as a shell with homogenous thickness. Some of the impacts of heterogeneous wall thickness on local conduction velocity will be captured by regional variations in conductivity in the i18 model. However, the model we have does not capture three potentially important impacts of wall thickness on atrial electrophysiology modelling. First, gradients in wall thickness will lead thicker tissue activating thinner tissue, or vice versa, depending on the wave direction. This will lead to direction-dependent electronic loading differences that will have activation wave direction-dependent effects on the conduction velocity that cannot be captured by our proposed modelling framework. Experiments in sheep indicated that gradients in wall thickness may impact the tethering of re-entrant activation patterns. 31 However, simulations in human left atrium, including wall thickness, did not find an impact of wall thickness gradients on re-entrant activation tethering in the left atrium. 23 Second, wall thickness raises the potential for distinct activation patterns to occur on the endocardium and epicardium, this would explain increased complex fractionated electrograms in thicker tissue. 32 Different activation waves on endocardial and epicardial tissue may lead to measurement artefacts that compromise the fitted model parameters. Third, our model cannot capture transmural differences in fibrosis, which simulation studies have shown could play a role in AF activation patterns. 33

Conclusions
We have created a virtual patient cohort encoding patient-specific electrophysiology heterogeneity. We have used the models combined with machine learning to show that combining local CV and APD with left atrial size can be used to identify tissue that can sustain phase singularities. Using a mathematical model that captures patient-specific dynamic atrial electrophysiology response we found that larger left atrium atrial size and local short APD and slow CV were correlated with re-entrant activation patter tethering. This is consistent with previous studies, supporting the critical mass hypothesis for sustain AF, that have correlated left atrium size and measures of left atrium wave length and ERP with AF. 17,34,35 Importantly our study indicates that it is not only the average but also the local CV and APD within the context of atrial size that determines the ability of an atria to tether re-entrant activation patterns. This may provide an alternate method for identifying ablation targets or guiding anti-arrhythmic drug therapies that adapt CV and APD, to remove these regions.

Supplementary material
Supplementary material is available at Europace online.