Late-Gadolinium Enhancement Interface Area and Electrophysiological Simulations Predict Arrhythmic Events in Patients With Nonischemic Dilated Cardiomyopathy

Objectives This study sought to investigate whether shape-based late gadolinium enhancement (LGE) metrics and simulations of re-entrant electrical activity are associated with arrhythmic events in patients with nonischemic dilated cardiomyopathy (NIDCM). Background The presence of LGE predicts life-threatening ventricular arrhythmias in NIDCM; however, risk stratification remains imprecise. LGE shape and simulations of electrical activity may be able to provide additional prognostic information. Methods Cardiac magnetic resonance (CMR)-LGE shape metrics were computed for a cohort of 156 patients with NIDCM and visible LGE and tested retrospectively for an association with an arrhythmic composite endpoint of sudden cardiac death and ventricular tachycardia. Computational models were created from images and used in conjunction with simulated stimulation protocols to assess the potential for re-entry induction in each patient’s scar morphology. A mechanistic analysis of the simulations was carried out to explain the associations. Results During a median follow-up of 1,611 (interquartile range: 881 to 2,341) days, 16 patients (10.3%) met the primary endpoint. In an inverse probability weighted Cox regression, the LGE–myocardial interface area (hazard ratio [HR]: 1.75; 95% confidence interval [CI]: 1.24 to 2.47; p = 0.001), number of simulated re-entries (HR: 1.40; 95% CI: 1.23 to 1.59; p < 0.01) and LGE volume (HR: 1.44; 95% CI: 1.07 to 1.94; p = 0.02) were associated with arrhythmic events. Computational modeling revealed repolarization heterogeneity and rate-dependent block of electrical wavefronts at the LGE–myocardial interface as putative arrhythmogenic mechanisms directly related to the LGE interface area. Conclusions The area of interface between scar and surviving myocardium, as well as simulated re-entrant activity, are associated with an elevated risk of major arrhythmic events in patients with NIDCM and LGE and represent novel risk predictors.

the CMR results. Mortality status was verified from the UK Health and Social Care Information Service. Cause of death was established from death certificates, postmortem and medical records.

Image Analysis
Entropy was measured by the standard Shannon Entropy [3], being applied to the native grey-scale image within the segmented fibrosis region [4], providing a measure of fibrotic disorder present. Units: dimensionless.
Volume was computed by summing together all pixels that contained LGE multiplied by slice thickness. Units: cm 3 .
Interface Area was measured by extracting the border between myocardium and LGE and adding together every borders' arclength multiplied by slice thickness. Note that the interface between the scar and the valvular annuli was not considered, neither was the interface between scar and epicardial or endocardial surfaces.

Units: cm 2
Transmurality was measured by a ray tracing method in which 580 rays emanated from the central pixel of the blood pool. The slice-specific transmurality score was then defined to be the mean fraction of LGE along each ray that intersected scar, while the patient's transmurality score was the mean of the slice-specific scores among slices with LGE. Units: dimensionless.
Number of Components was computed by a connected component algorithm, defined to be the total number of 4-connected regions of LGE in all slices.
Radiality was measured in each slice as the angular variance of all LGE pixels with respect to the central blood pool pixel, and for each patient as the mean score of the slices containing LGE. Units: normalised angle (radians/2π).
All metrics were computed with an in-house Python script which is freely available at www.github.com/GabrielBalabanResearch/lgemri_scar_metrics. Figure S1 shows histograms of all computed LGE metrics across the full patient cohort.

Calculation of the Inverse Probability of Weights
Inverse probability weights (IPW) were calculated for each LGE metric using a non-parametric method [5]. In brief the method solves an optimization problem that generates a set of weights which maximally de-correlate the target variable with the potential confounders, while maximizing the empirical likelihood of observing the data. We implemented the IPW method in Python and used the L-BFGS-B solver from Scipy to solve the optimization problem. Our potential confounder set consisted of all baseline variables listed in Table 1 (main article) that were not LGE metrics, as well as ICD or CRT receipt during follow-up. Figure S2 shows the distributions of the resulting IPW; the mode of each distribution is very close to 1, and most weights have value less then 4. In Figure S3 we visualise the Pearson correlations between the LGE metrics and other baseline variables. All correlations were reduced to 0 after reweighing.
For the simulated reentries variable, it was not possible to calculate IPW using all baseline variables due to the large number of patients (124) with simulated reentries = 0. Instead, we focused on a limited set of the most important confounders, moderate alcohol excess, and NYHA class (which were significantly different in the event and non-event patient groups) and ICD receipt. The distribution of the simulated reentry variable, as well as it's IPW and Pearson correlations are shown in Figure S4.

Computational Simulations
Computational simulations were performed as part our previous study on 2D scar shapes and microstructures [6]. In summary, our method for running computational simulations was as follows: segmented images containing LGE were processed as individual short-axis slices and meshed (computational geometry algorithms library CGAL: https://doc.cgal.org/) into triangular finite element models with maximum edge length 250um. A monodomain representation was used to simulate electrical activity with a human ventricular cell model [7] implemented in the Cardiac Arrhythmia Research Package [8]. Realistic myo-fibre architecture was incorporated using a rulebased method [9]. Conductivities were tuned to match experimentally observed conduction velocities [10], with conduction suitably modulated within scar regions, as performed in our previous study [6].
Patchy distributions of fibrosis were represented using the widely-used percolation method [6,11].
Briefly, areas defined as scar in image segmentations had triangular elements or mesh edges randomly removed, with a probability that depended on the normalized local image intensity on the LGE image.
In this manner, less compact regions of fibrosis (characteristic of NIDCM [12]), could be incorporated into the models, allowing conduction of activation waves through the scar substrate, albeit via slowed and convoluted pathways [13].
Simulated programmed electrical stimulation was performed from an endocardial pacing location determined to be in a consistent location with respect to scar in each model (obtained by minimizing a cost function [6]). A steady-state train of 4 S1 pulses at 500ms coupling intervals was followed by up to 5 shorter extra-stimuli, the timing of each determined by the local effective refractory period (ERP). Following each stimulus sequence, 600ms of activity was simulated and any reentrant transmural reentrant activations were recorded, that is activation waves which reversed direction and reactivated the original pacing location. A detailed analysis of activation wavefront dynamics was performed to identify the mechanistic role of different LGE metrics defining the scar in unidirectional conduction block and reentry initiation. A total of 10 simulations were run per image slice, testing 10 scar microstructure combinations, with 5 different scar density levels and 2 fibrosis types (interstitial and replacement fibrosis represented by mesh triangle and mesh edge removal respectively), for further details see [6]. The pacing location, and assignment of electrical and ionic properties to tissues, were kept constant in each image slice. Figure 4 in the main manuscript shows 2 simulations highlighting the specific role of the scarmyocardial interface in the genesis of reentry. Below, we provide a more in-depth description of these simulations.

Simulation 1:
Panel C (1950ms) shows the wavefront 100ms following the final S1 beat. Here, the wavefront passes easily through the healthy myocardium at the epicardial and endocardial tissue edges (as the scar is largely mid-myocardial in this example). However, activation of the scar itself does occur, albeit greatly delayed due to the tortuous conduction pathways the wave is forced to take through the fibrotic tissue. The panel C (2470ms) shows the wavefront 120ms after the final S4 of the stimulus train. Due to the reduced diastolic interval, conduction blocks at along the interface shown by the green parallel lines. However, it continues to propagate along the healthy myocardium in the sub-epicardial region, and conduction occurs (albeit slowly) across the interface in this area, as before.
Panel C (2620ms) shows that as activation propagates slowly through the scar, it reenters back through the initial site of block at the scar interface.

Simulation 2:
Panel D (1980ms) shows the wavefront 130ms following the final S1 beat. Here, the wavefront initially activates the more excitable healthy (sub-endocardial) myocardium, and slowly spreads and uniformly activates the fibrotic area extending through the midwall and towards the epicardium. Panel D (2530ms) shows the wavefront 80ms after the final S4 of the stimulus train. A strong heterogeneity in repolarisation is witnessed throughout the scarred region (following the S3 beat at 2100ms) as well as between the scar and the surrounding healthy myocardium, causing very slow, tortuous propagation through the scar and a zone of complete conduction block. Panel D (2900ms) then shows that as the healthy myocardium fully recovers from successful activation by the S4 beat it, may be re-activated by the wavefront that is still slowly 'hiding-out' in the scar.