Assessing the ability of substrate mapping techniques to guide ventricular tachycardia ablation using computational modelling

Background: Identification of targets for ablation of post-infarction ventricular tachycardias (VTs) remains challenging, often requiring arrhythmia induction to delineate the reentrant circuit. This carries a risk for the patient and may not be feasible. Substrate mapping has emerged as a safer strategy to uncover arrhythmogenic regions. However, VT recurrence remains common. Goal: To use computer simulations to assess the ability of different substrate mapping approaches to identify VT


Introduction
Recurrent ventricular tachycardias (VTs) are associated with an increased risk of sudden death in patients with ischemic heart disease (IHD) [1,2]. Catheter ablation has become an increasingly utilized treatment for post-infarction VTs [3]. However, success rate remains sub-optimal, particularly in structural heart disease (50-70% versus 80-90% in the structurally normal heart [4]) as ablation therapy is highly depended on the ability of current mapping systems to accurately locate areas of the heart muscle causing the VT.
VTs are often associated with circulating electrical wavefronts that are known to be formed at the so-called infarct border zone (BZ). The BZ is comprised of a heterogeneous admixture of remodelled myocardium and patchy fibrosis that may form a substrate for functional conduction block and reentry via channel isthmuses within the scar [5][6][7]. Efficacy of VT ablation in IHD relies on the ability of current mapping systems to identify such isthmuses and deliver adequate lesion sets to the critical arrhythmogenic sites [8].
Identification of the reentrant circuit often requires arrhythmiainduction. However, VT in patients with IHD is often hemodynamically non-tolerated or non-inducible [9]. In these cases mapping can still be performed to uncover the arrhythmogenic substrate during sinus or paced rhythm. Substrate-based ablation strategies consist of locating areas with abnormal electrograms or slow conduction as these are believed to correspond to critical VT sites [8,10]. However, electrogram characteristics have been shown to have limited specificity with the VT isthmus [11], which could result in over-or underestimation of ablation lesions [12]. Moreover, the delineation of the isthmus amid the scar can be challenging particularly in the presence of intramural substrates, where the re-entry circuit is deep into the ventricular wall [13]. Given these many limitations, the most effective substrate mapping strategy to guide VT ablation in the infarcted heart is still undetermined.
Computational modelling has been widely used to investigate mechanisms of post-infarction VTs [7,14,15], perform patient risk assessment [16], optimize novel arrhythmia mapping metrics [17,18] and more recently to guide personalized therapy planning [19]. However, there is often a disconnection between how simulated signals relate to clinically obtained signals as most in-silico studies use transmembrane potentials, instead of typical unipolar (UE) and bipolar (BE) electrograms measured by available electroanatomical mapping systems. Computational models constructed based on detailed imaging data have the capability to understand the origin of signals in relation to the detailed infarct structure, and thus guide mapping techniques in a safe, controlled and robust manner.
The goal of this study is to use a structurally detailed model of the infarcted left ventricle (LV) reconstructed from in-vivo late gadolinium enhancement magnetic resonance imaging (MRI) data to: 1) investigate the electrical signature (electrogram morphology and local conduction parameters) of intramural BZ structures within the scar; 2) relate the VT exit site to vulnerable regions detected by electroanatomical maps obtained with typical substrate mapping techniques as well as with the Reentry Vulnerability Index (RVI), an activation-repolarization metric, also aimed to identify vulnerable regions without inducing VT [18,20,21].

Geometrical model of the porcine LV
The biophysical simulations in this study were carried out in a geometrical finite element (FE) model built based on contrast-enhanced MRI scans of the in-vivo porcine post-infarction heart. Model generation and simulation protocols are described in detail below.

In-vivo MRI data
The geometrical model of the porcine heart was built based on lategadolinium enhanced cardiovascular MRI with an isotropic voxel resolution of 1 mm obtained seven weeks following myocardial infarction [22]. The infarct region was detected using a signal intensity thresholding approach as follows: the left ventricular (LV) myocardium, blood pool, aortic root and coronary ostia were manually segmented by a single observer within Seg3D (www.seg3d.org). Papillary muscles were included in the LV myocardial segmentation only if there was no visible blood pool between the muscle body and the endocardial surface. The signal intensity range within the segmented myocardium was used to calculate the threshold for scar at 60% of the maximum intensity within the image, and the threshold for heterogeneous tissue at 40% as previously reported and validated [23].  1. Pipeline used to construct the LV model from MRI stacks. A) A representative late gadolinium enhanced MRI slice of the infarcted porcine heart. B) Segmentation of the MRI data into blood pool (green), healthy myocardium (red), infarct BZ (gray) and scar (black). C) 3D computational model of the LV geometry. Endocardial cutaway (left) and epicardial (right) surfaces of the LV. S1-S2 stimuli were applied to cells located at the LV base. Fig. 1 illustrates the pipeline used in our group to convert MRI scans into geometrical models of the heart. The tetrahedral FE mesh was created based on the segmented MRI data using Tarantula (CAE Software Solutions, Eggenburg, Austria). The resulting LV model (Fig. 1C) consisted of 4 274 087 myocardial nodes connecting 23 663 301 FEs with an average edge length of 338 ± 127 μm [24]. A rule-based approach was employed to compute and assign fiber orientations into the LV model [25].

Anatomical analysis
An eikonal model was used to obtain quantitative measurements of wall thickness (WT) and percentage transmurality of scar and fibrosis. Briefly, the algorithm solves the eikonal equation on the ventricular domain with initial activations prescribed on the epicardial surface. Also, fixed propagation velocities v myo = 1.0 m/s in the longitudinal, transversal and normal to the fiber directions were enforced. The resulting scalar values represent the arrival times in the myocardium t myo . From the combination of t myo and v myo , the distance between the epicardial and endocardial surfaces can be computed, providing thus WT values throughout the LV: The amount of fibrosis in the post-infarction LV was assessed by solving a second eikonal problem with a slower velocity (v scar = 0.1 m/s) prescribed within the scar (shown in black in Fig. 1C). The percentage of fibrotic tissue FIB across the ventricular wall can then be estimated as follows: where t myo and t scar are the arrival times obtained after solving the eikonal models with v myo = v scar = 1.0 m/s, and v myo = 1.0 m/s and v scar = 0.1 m/s, respectively.

Prediction of electrogram amplitudes
Unipolar (UV) and bipolar (BV) voltage amplitudes within the porcine post-infarction model were estimated with the equations proposed by Glashan et al. [26]: where WT and FIB(%) were obtained by solving Eqs. (1) and (2).

Biophysical simulations
Cardiac electrical activity within the LV geometry was simulated according to the monodomain equations: ∂η where σ m = diag(σ ml , σ mt , σ mt ) is the effective bulk conductivity tensor [27], which combined with the fiber orientations (obtained with a widely used rule-based approach [25]) accounts for the cardiac anatomy and anisotropy; V m is the transmembrane voltage; β = 0.14 μm − 1 is the surface to volume ratio; I m is the transmembrane current density; C m is the membrane capacitance per unit area; I ion is the density of the total ionic current flowing through the cell membrane; and I stim is the stimulus current density.
Anisotropic conductivities values of σ ml = 0.143 S/m and σ mt = 0.064 S/m were assigned to the healthy myocardium to produce longitudinal and transverse velocities within the range recorded in the ventricles [28]. The ten Tusscher (TT) model [29] was used to simulate cellular membrane dynamics, such as channel gating and ionic concentrations (represented by η in Eq. (7)) given by I ion . The maximum conductances of the rapid (I Kr ) and slow (I Ks ) delayed rectifier currents were increased by 50% to shorten the AP duration (APD) from 299 ms (epicardial phenotype) to 255 ms (basic cycle length of 1000 ms). While the resulting APD is still above the average 220 ms reported porcine ventricular myocytes [30], it was enough to shorten the cardiac wavelength allowing for sustained monomorphic VT induction in the post-infarction model. Experimental data obtained at the chronic phase of myocardial infarction (the time point in which the MRI data used in this study was obtained) suggest that remodeling processes promoting conduction slowing as opposed to APD prolongation dominate [6]. Moreover, in a previous in-silico investigation we showed that in simulations with prolonged APD, reentries were often not sustained because of the subsequent block caused by long-lasting repolarization (down-regulated I Kr and I Ks ) [7]. Thus, macroscopic electrophysiological properties of the BZ of healed infarcts were represented here by: 1) reducing the maximum conductance of fast sodium current (I Na ) to 25% of its control value [7]; and 2) prescribing isotropic tissue conductivity with values of σ ml = σ mt = 0.064 S/m [31]. The choice to only decrease I Na and combine it with isotropic tissue conductivity in the BZ to account for the presence of interstitial and patchy fibrosis as well as fiber disarray, all reported alterations in the BZ of healed infarcts [6]. The scar was represented as being nonconducting (i.e., by imposing no-flux boundary condition at its interface).

Calculation of unipolar electrograms
The set of monodomain Eqs. (5)-(7) describes only the potential across the cell membrane V m . However, potentials resulting from a current flow from the cardiac cells to the extracellular space can still be computed by assuming that the heart is immersed in a uniform volume conductor of infinite extent and conductivity σ b = 1.0 S/m (taken as the measure reported for blood [33]). Unipolar electrograms (UEs) can be estimated according to the following equation [34]: where φ e are the UEs at a given field point; r is the distance vector between field points and the source term βI m given by Eq. (5) within the intracellular domain Ω defined by the FE mesh; and d = 7.14 μm (average cell radius) is the electrode-tissue distance added to r in order to ensure convergence of the integral.

Pacing protocol
The LV model was initialized with stabilized single-cell model states obtained after pacing the TT model for 100 cycles with a basic cycle length of 1000 ms. If not stated otherwise, electrical activity was initiated at the base of the LV model chosen to roughly represent typical pacing through the coronary artery. Three S1 beats followed by a premature S2 with a shorter coupling interval (CI) were simulated.
VT induction: VT was induced in the LV model by shortening the S2 CI until a sustained electrical reentry was detected.
Fast pacing: fast pacing is considered here as the shortest S2 CI which did not induce reentry.

LV model
Simulated action potentials (APs) of all nodes in the LV model were stored and post-processed. VT dynamics were analyzed with maps constructed from ATs extracted from APs. ATs were taking as the time at which V m crossed − 20 mV with a positive gradient. The VT exit site was detected visually from the isochrone maps and is considered here as the "reference solution".

Mapping grid
While APs are readily available for all points in computational models, this is not the case in the clinical setting. Typical electroanatomic mapping systems only allow operators to record electrograms (rather than APs) within the cardiac chamber, limiting the size of the data collected with multipolar catheters [21]. In order to replicate clinically-relevant conditions, a mapping grid resembling the inner LV chamber was created by coarsening the endocardial surface shown in Fig. 1C. The resulting grid (1 073 recording electrodes with an average spacing of 3 829 ± 751 μm) is shown in Fig. 2A. UEs (φ e in Eq. (8)) were calculated for each electrode in the mapping grid. BEs were computed by taking the difference between UEs from electrodes sharing the same triangular surface as illustrated in Fig. 2. The BE with the largest amplitude amongst the three was assigned to the center of the triangular surface. The surface data were then interpolated back onto the nodes (electrodes) of the mapping grid.

Feature extraction
The following features were extracted from electrograms in the mapping grid ( Fig. 2A): • Unipolar voltage (UV): the amplitude of φ e , calculated as the difference between its minimum and maximum values; • Bipolar voltage (BV): the amplitude of the bipolar signal calculated as the difference between its positive and negative peaks; • Fractionation index (FI): the number of separated activation events in the time derivative dφ e /dt was calculated based on a modified template-based algorithm [35]. Briefly, the cross-correlation of a template (derived from φ e and dφ e /dt) and the signal under investigation was performed. The number of peaks in the cross-correlation function (CCF) corresponds to the number of activation events. Peaks with less than 20% of the maximum CCF or separated by less than 10 ms were ignored;

RVI calculation
Along with traditional substrate mapping techniques, the capability of the Reentry Vulnerability Index (RVI) [18,20,39] to detect the VT exit site was also assessed in this study. Reentry initiation around a circuit requires that a wavefront traveling along a line of block finds tissue that has already regained excitability, enabling its reactivation. The RVI is a activation-repolarization metric that provides a point-by-point quantification of this principle. RVI values are computed by taking the time interval between the regaining of excitability (repolarization) of tissue and the arrival of a subsequent wavefront at the same site [39]. The RVI value between a pair of recording electrodes i and j is given by: The metric in Eq. (9) has been extended and optimized to construct a global map that can be used to identify spatial regions with high susceptibility to reentry [18,20]. The algorithm consists of four steps: 1) ATs and RTs of the S2 beat are calculated for all recording electrodes; 2) for a given recording site i, all other sites j that are activated later than site i and lie within a search radius (R = 16 mm) are found; 3) the RVI for each recording site pair RVI(i, j) as in Eq. (9) is calculated; 4) the minimum RVI among all RVI(i, j) is found and associated with node i to create an electroanatomical map. During normal propagation, RVI maps resemble the distribution of AP duration of the tissue. In the presence of block, on the other hand, absolute RVIs are smaller than the APD, becoming negative if a reentry is induced [18].

Detection of the VT exit site using substrate mapping
Substrate maps constructed based on all features described above were used to assess their capability of localizing the VT exit site during fast (non-arrhythmogenic) pacing. Potential ablation targets with each map were detected following previous studies [11,21,40,41]: • Signal amplitude: UEs or BEs with voltages smaller than prescribed cut-offs of 0.5 mV and 1.0 mV; • Signal fractionation: UEs with the highest 5% FI; • AT gradient: the spatial gradient of the AT map was computed and regions with the highest 5% values were selected; • RVI: areas with lowest 5% RVI values.
Specificity and sensitivity of all maps were assessed by calculating the minimum (d min ) and average (d avg ± STD) distances between each point within the detected region and the actual exit site obtained from VT mapping [21,41]. Finally, Pearson correlation coefficients between simulated metrics and ventricular anatomical properties (WT and fibrosis) were also computed.

Simulation of ablation lesions
Experiments with canine hearts showed that ablation lesions can vary from 8 mm to 14.8 mm in diameter and 5.9 mm-11.2 mm in depth [42]. In the post-infarction LV model, the ablation was created with a radius of 10 mm to ensure similar lesion penetration through the ventricular wall. Tissue residing within the lesion was set to be nonconducting as the necrotic scar.

VT mapping
VT was induced in the LV model by S1-S2 stimulation as described above. Fig. 3A shows V m maps of the (last) S1 and S2 beats as well as the resulting sustained reentry in the LV model. VT was induced after a S2 beat with a CI of 350 ms had initially blocked at the isthmus's mouth proximal to the stimulus site (t = 2570 ms). As can be seen in Fig. 3A, the wavefront travels around the scar entering the infarct region from a distal site at t = 2720 ms. The wavefront then propagates back towards the base exiting through the isthmus (t = 2940 ms) from where it reenters the myocardium. See Supplemental Video 1 for further details of the block and reentry within the LV model.
Supplementary video related to this article can be found at https://d oi.org/10.1016/j.compbiomed.2021.104214 Fig. 3B presents AT maps of the first cycle of the VT seen from both epicardial and endocardial surfaces of the LV. The VT exit site is indicated by a white double circle. Note that the scar is not fully transmural and that the isthmus is located in the epicardium. The bottom/right panel of Fig. 3B shows ATs computed in the recording grid (based on UEs rather than APs). The scar region (in black) was determined by mapping the fibrotic elements of the computational mesh onto the recording grid. The exit site indicated in the mapping grid is considered as the "reference solution" in this study.

WT and percentage transmurality of scar and fibrosis
Quantitative measurements of WT and percentage transmurality of scar and fibrosis were performed on the LV. These results are presented in Fig. 4. Note that WT varied between 2 and 16 mm with thinnest regions located at the LV apex. As can be seen in Fig. 4B, the infarct scar is very heterogeneous across the septal wall. The isthmus sustaining the VT in Fig. 3 can bee seen in light blue between the more scarred areas in yellow.

Substrate mapping
All electroanatomical maps in this section were computed during fast non-arrhythmogenic pacing, i.e., in absence of VT. The in-silico data presented below was obtained after applying a premature S2 beat with a longer CI (400 ms) than the one required to induce the VT shown in Fig. 3 (CI = 350 ms). Unlike in the scenario depicted in Fig. 3A, the S2 beat is not able to reenter the myocardium despite of the block at the isthmus' proximal mouth (see Supplemental Video 2). Figs. 5-6 show maps constructed based on features extracted from both UEs and BEs (see Fig. 2). The correlation between features extracted from computed electrograms and LV anatomy (WT and scar transmurality) is presented in Fig. 7.
Although simulated signals with small amplitudes were found in and around the scar, they lacked specificity for the identification of exit sites. Areas with tissue thinning, such as the apex, or with large amounts of fibrosis also resulted in small voltages. Simulated UVs and BVs had similar low correlation coefficients ( Fig. 7) with WT and fibrosis.
The scar size estimated by removing signals with a prescribed cut-off is shown in black. Removing UEs with amplitudes smaller than 0.5 mV resulted in a smaller scar area when compared to the removal of BEs (mid panels). In both UE and BE maps, increasing the cut-off value to 1.0 mV further increased the scar region. However, the electrical signature of the scar resulting from removal of BEs was overall more similar to the reference solution in Fig. 3B (mapping grid) than the one obtained by applying cut-offs to UEs.
Additional maps computed based on fractionation and time of activation of UEs are shown in Fig. 6. Signals with higher FI were found near the border of the infarct scar. Note that non-fractionated signals were present within the scar core. These were removed after applying the 0.5 mV cut-off to UEs (mid panel). The highest correlation (0.73) was observed between scar transmurality and UE fractionation (see Fig. 7). Fig. 6B AT maps resulting from the premature S2 beat. When all signals are considered (bottom/left panel), isochrone crowding can be seen next to the area with latest AT. As UEs are removed (from left to right), the scar region increases and zones of slow conduction are no longer seen. As can be seen in Fig. 7 AT gradient as well as the RVI had the lowest correlation with all LV features.

Detection of the VT exit site
Next, the electroanatomical maps in Figs. 5-6 had their ability to detect the VT exit site assessed. The effect of combining voltage cut-offs with the other metrics under investigation here was also assessed by recomputing d avg and d min after removal of UEs with amplitudes smaller than the cut-off. In Fig. 8 vulnerable regions estimated solely by voltage cut-offs are shown in black, while vulnerable regions detected by the other metrics are highlighted in red. The exit site is indicated in all maps by a white double circle as base for comparison.
As can be seen in Fig. 5 (left to right), increasing the voltage cut-off from 0.5 mV to 1.0 mV also increases the size of the estimated vulnerable region (in black). Increasing the cut-off to 1.0 mV resulted in a decrease in d min for both UVs (from 18.7 mm to 4.5 mm) and BVs (from 9.0 mm to 6.1 mm). d avg of UVs, on the other hand, increased while remaining similar for the BVs. Fig. 8A and B highlight regions with highest 5% UE fractionation and AT gradient values, respectively, computed from all maps in Fig. 6. The minimum distance between the exit site and the closest fractionated signal was d min = 14.7 mm d min decreased while d avg became larger when voltage cut-offs were applied (from top to bottom).
Areas with highest AT gradients were found around the scar (Fig. 8B  mid panels). d min between the recording site with largest AT gradient was 13.1 mm when all signals were taken into account for the analysis. Removing UEs reduced d min to 3.2 mm and to 0.0 mm, for voltage cutoffs of 0.5 mV and 1.0 mV, respectively. d avg also decreased when a cut-off of 1.0 mV was applied.
The ability of RVI maps to detect the VT exit site was investigated in        When all signal were considered, the RVI performed better (smaller d min and d avg ) than all other metrics. The removal of UEs increased d min . Note that when a 0.5 mV cut-off was considered, small RVI values were found at more distant regions (apex) which in turn increased d mean from 28.2 mm to 36.1 mm. Next, a simple ablation was performed to the post-infarction heart model followed by a test of VT non-inducibility. The RVI map in Fig. 8C (top panel) was used to guide the ablation since this metric, when used alone (no voltage cut-offs), identified regions closest to VT exit site. To assess the success of the simulated ablation, the exact S1-S2 pacing protocol which successfully induced sustained reentry in Fig. 3 was repeated in the ablated model. As can be seen Supplemental Video 3, VT was no longer inducible in the newly ablated model.

Influence of the pacing location
In order to quantify the impact of wavefront incidence on the electroanatomical maps in Fig. 8, the LV was paced at the apex (see Supplemental Video 4 for more details) and all electrogram features used to construct the maps were recomputed. Fig. 9 presents maps highlighting the vulnerable regions detected after the S2 beat (CI = 400 ms).
Supplementary video related to this article can be found at https://d oi.org/10.1016/j.compbiomed.2021.104214 UVs, on average, remained similar when the pacing location was changed from the base (2.8 ± 1.3 mV) to the apex (2.6 ± 1.4 mV). The estimated ablation regions resulting after applying voltage cut-offs of 0.5 and 1.0 mV can be appreciated in Fig. 9 (black areas). It can be seen that the size of the estimated vulnerable regions increased when changing the stimulus from the base (Fig. 8) to the apex (Fig. 9). This was also reflected in d avg of UVs which increased from 29.7 ± 8.1 mm (base, Fig. 5A) to 34.2 ± 12.6 mm (apex, map not shown) with a cut-off 0.5 mV d min , on the other hand, decreased from 18.7 mm to 8.0 mm when the stimulus location was moved from the base to the apex.
The combination of voltage cut-offs and the other metrics was also assessed by recomputing d avg and d min after removal of UEs. Overall, the spatial distribution of the vulnerable regions estimated by all metrics (red areas in Figs. [8][9] were affected by change in pacing location. When all signals were considered in the analysis, a decrease in d min was observed for all metrics, but the RVI (Fig. 8C and 9C).

Discussion
In this study we employed computer simulations to investigate: 1) the relation between clinically significant signals and post-infarction scars; and 2) the ability of metrics derived from these signals to guide ablation without inducing VT. Simulations within a structurally detailed model of the infarcted LV showed that low UVs and BVs were found in and around the scar. However, BEs provided a better electroanatomic picture of the scar. Moreover, as expected, our results suggest that while signal amplitude cut-offs can provide information about the infarct anatomy, they alone are not capable of identifying the VT exit site. Electrogram fractionation had the highest correlation with scar transmurality. Among all metrics assessed, when all signals (including dense scar regions) were considered, the RVI identified the vulnerable region closest to the VT exit site. The VT exit site could only be accurately detected (d min = 0) when areas with large AT gradients were computed after applying a voltage cut-off on UEs (1.0 mV). Finally, all metrics were shown to be affected by the direction of the wavefront, suggesting that multi-site pacing may be necessary to detect all vulnerable regions.

Ablation of VTs versus substrate ablation in the post-infarction heart
The mechanisms leading to the initiation and maintenance of reentrant waves in the presence of an infarct scar have been one of the most important topics of research in computational electrocardiology [7,[43][44][45][46]. Previous studies have employed simpler geometrical models to investigate characteristics of post-infarction VTs [45,46]. Here, a  structurally detailed model of the infarcted LV was used to address the clinical question of whether typical substrate mapping approaches can identify the VT exit site. Mapping of cardiac tissue sustaining reentry requires VT induction as illustrated in Fig. 3. However, VT inducibility poses a risk for the patient as haemodynamic compromise can occur during the procedure. Unlike in the structurally normal hearts, VTs in patients with IHD are often hemodynamically unstable [9]. Another limitation of VT mapping is that VT inducibility may not be achieved. Indeed, VT in the LV model could not be induced by pacing at the apex. Absence of unidirectional conduction block (a pre-requisite for reentry formation) may be one reason for non-inducibility. This directional dependency of the substrate has been demonstrated in a simulation study to assess mortality risk in post-infarction patients [16], where the virtual hearts were paced from different locations. However, multi-site pacing may not be possible due to time constraints of the ablation procedure. If the VT can not be induced or is not hemodynamically tolerated, ablation can be guided by the identification of the arrhythmogenic substrate during sinus rhythm or by controlled pacing protocols.
Another factor contributing to the unsatisfactory ablation outcomes is infarct anatomy [13]. While more common in patients with non-ischemic cardiomyopathy, intramural or epicardial VT isthmuses, such as the one in Fig. 3, have also been reported in patients with IHD [47]. Critical VT sites have been shown to be located in areas with more than 25% scar transmurality [48]. In the model used here, the VT isthmus, when seen from the endocardium, corresponds to the areas of light blue (̃30% transmurality) between the more transmural scar (yellow and red) in Fig. 4B. The presence of such components poses a challenge for current mapping systems. Since most clinical VTs are ablated endocardially, all metrics investigated here were based on features obtained from endocardial signals.

Relation between low-voltage areas and the infarct scar
The relation between tissue structure and electrograms has been previously investigated in experimental and computational studies [26,49,50]. Electrogram amplitude (voltage maps) is often used as a metric to distinguish between healthy myocardium from scar tissue.
Glashan et al. [26] proposed equations to estimate UV and BV based on the amount of fibrosis and WT in non-ischemic cardiomyopathy. Although their approach may not hold for ischemic cardiomyopathy, it was nonetheless included here as an important basis for comparison since the fundamental biophysical arguments upon which their formula are based (reduced electrogram amplitude with reducing viable tissue under the electrode due to either fibrotic infiltration or thinned myocardial wall) should hold, to some degree, for ischemic scar. While estimated UVs (6.9 ± 0.6 mV) are below values measured in remote myocardium without scar (8.32 ± 3.39 mV) [51], the minimum estimated BVs in Fig. 5 (> 3 mV) are above standard values used to detect scar tissue in patients with IHD (0.5 mV) [8,10]. This discrepancy is likely to be due to the fact that the authors estimated the parameters in Eqs (3) and (4) using multiple linear regression performed on ex-vivo data of patients with non-ischemic heart disease [26].
Marchlinski et al. [40] have demonstrated that areas with voltage above 1.5 mV are consistent with normal tissue while sites with amplitudes smaller than 0.5 mV are associated with the presence of a dense scar. Moreover, voltages between 0.5 and 1.5 mV are believed to represent surviving myocardial bundles interspersed with fibrosis in the BZ [40], and have been thus been used to locate conducting channels supporting VT [8]. Here, voltage cut-offs 0.5 and 1.0 mV were applied to simulated UEs and BEs (black areas in Fig. 5) to characterize the infarct anatomy. A cut-off of 1.0 mV instead of the typical 1.5 mV was used because the latter would lead to an overestimation of scar/BZ regions. This is because the simulated signals were obtained with Eq. (8) assuming that the LV is immersed in a uniform volume conductor with conductivity σ b = 1.0 S/m. Increasing σ b would increase simulated UVs and BVs, but would not change their distributions. In accordance with experimental and clinical evidence, signals with small amplitudes were found in and around the scar (Fig. 5). Furthermore, the results suggest that BEs might offer a better electroanatomical picture of the scar than UEs. Ablation of all low voltage areas would result in scar homogenization, an approach that has been shown to be effective at reducing VT recurrence in IHD [52].

Fractionated electrograms
Another approach to substrate ablation is to identify and ablate abnormal electrograms as they are thought to represent the VT isthmus [11,53]. Abnormal electrograms have been shown to originate from the propagation across surviving viable myocardium separated by fibrosis following complex zig-zag patterns [50,54,55]. The presence of such tortuous conduction pathways produce abnormal signals characterized by multiple deflections representing multiple depolarization of the different surviving myocardial fibers under the recording catheter. These abnormal electrograms are typically classified as fractionated, split or late potentials [53]. Here fractionated signals are defined as UEs with multiple deflections in their time derivative [35].
As shown in Fig. 6A, the electroanatomical map constructed based on signal fractionation could reveal areas of uniform UEs (blue) surrounded by fractionated ones (green/red areas). The highest correlation (0.73) was also observed between UE fractionation and scar transmurality (see Fig. 7). Such strong correlation is likely because of the presence of surviving myocardial bundles in areas of more scar. Critical VT sites have been shown to be located in areas with more than 25% scar transmurality [48]. In the model used here, the VT isthmus, when seen from the endocardium, corresponds to the areas of light blue (̃30% transmurality) between the more transmural scar (yellow and red) in Fig. 4B. However, there was no fractionated signals near the exit site (the closest one was located 5.5 mm away, Fig. 9). Moreover, d min decreased when a voltage cut-off of 0.5 mV was applied suggesting that the fractionation metric might benefited from the removal of UEs. In fact, Desjardins et al. [51] showed that UVs are more useful to detected intramural scars, such as the one in Fig. 1, than BVs. However, fractionated signals are known to have smaller amplitudes than uniform ones [8] and voltage cut-offs may result in an undesired removal of signals that could otherwise provide important information about the substrate. As shown in Fig. 8A and 9A (top to bottom) the size of the regions with highest FI decreased to a few points once all UEs with voltages < 1.0 mV were removed. Furthermore, in experiments with swine with healed infarction, Anter and colleagues [11] showed that VT critical zones consistently occur in areas with abnormal electrograms with voltages lower than 0.55 mV during sinus rhythm. However, the authors also showed that these abnormal signals are not specific for isthmus sites. As can be seen in Fig. 6, fractionated UEs were also observed near the base as well as in other regions of the LV.

Activation mapping
Isochronal activation mapping is another common metric used for ablation during sinus rhythm. Activation maps are used for visualization of slow conduction within scar to identify deceleration zones, which have been associated with sites critical for reentry [11,56,57]. Indeed, the VT critical zone was shown to correspond to a region of steep activation gradient [11]. In our simulations, most of the zones with slow conduction were found to be within and around the scar, but some deceleration can be seen near the exit site (Fig. 6B). Here, areas marked for ablation were identified as having the highest 5% gradients in ATs (Fig. 9B). While alone this approach was not able to identify the exit site, when applied in combination with a signal cut-off of 1.0 mV it could accurately detected the exit of the VT circuit (d min = 0 in Fig. 8B, bottom panel).

Activation-repolarization metrics
The initiation and maintenance of electrical reentries depend on the interaction between the front (activation) and back (repolarization) of the cardiac wavefront [58]. Despite the evidence about the link between spatio-temporal repolarization dynamics and arrhythmogenesis [39], repolarization is often neglected during the mapping of the arrhythmogenic substrate. The RVI is a metric proposed to overcome this limitation by combining both ATs and RTs to identify vulnerable regions without the need to induce VT [20].
The RVI has been thoroughly assessed and optimized in in-silico studies for its use in the clinical setting [17,18]. Due to its ability to predict VT sites of origin in patients with right ventricular pathology [41] and more recently in patients with scar-related VT [21], the RVI was included in our analysis and compared to traditional substrate mapping techniques. Since the RVI is an absolute time difference, its interpretation depend on factors such as presence of conduction block and activation-repolarization dynamics that are patient-specific and can change from beat to beat. In order to cope with this limitation, we followed the approach of Orini et al. [21] and sought for areas representing the bottom 5% of all RVI values. Although the RVI had the lowest correlation with scar transmurality (see Fig. 7), it could identify regions closest to VT exit site (d min = 3.2 mm) when all signals, including those recorded on the dense scar, were used to construct the map.
A possible explanation for the better performance of the RVI metric at identifying the exit site is the use of repolarization to detect the presence of unidirectional conduction block. As can be appreciated in Supplemental Video 2, the S2 beat blocks at the proximal mouth and travels around the scar propagating back towards the base where it blocks again. Indeed, conduction block was not present when the LV was paced from the apex (Supplemental Video 4), which resulted in an increase in d min from 3.2 mm to 9.2 mm. According to the definition of the RVI, regions of conduction block are associated with RVI values close to zero [20,39]. Even though the block occurred at the epicardium, its presence could still be detected from ATs and RTs taken from endocardial UEs. This finding suggests that the RVI can detect regions susceptible to conduction block even in infarcts with non-transmural scars. The ability to identify components of the VT circuit, including those found deep within the scar such as the one in Fig. 3, is clinically relevant because infarct anatomy is another factor contributing to the unsatisfactory ablation outcomes [13].

Direction of the wavefront
Substrate-based ablation strategies rely on features extracted from electrograms to obtain information of the underlying tissue structure/ physiology. However, amplitude and morphology of both unipolar and bipolar electrograms is influenced by many factors such as electrode size and spacing, conduction velocity and direction of the activation wavefront [59]. In particular, the influence of wavefront incidence on electrograms has been investigated using multi-site pacing [50] and by rotating the position of the electrodes [60]. The effects of wavefront incidence on substrate mapping techniques was assessed here by changing the initiation of cardiac propagation from the base of the LV (Fig. 8) to its apex (Fig. 9). Moreover, a multicenter study conducted in patients with infarct-related VT showed that the location as well as the magnitude of the conduction slowing changed depending on LV activation [57]. Overall, all metrics investigated here were affected by the change in pacing location. This suggests that a combination of maps created with different pacing locations may be necessary for the accurate identification of vulnerable regions. However, this would lead to an increase in the time of the ablation procedure.

Limitations
The main limitation of this study is the use of only one post-infarction heart model. However, the use of single scar anatomy model allowed a detailed analysis of typical metrics used for substrate ablation strategies. Further, the model was constructed based on MRI data from a swine heart. The swine heart shares several anatomic and physiologic characteristics with the human one and they are often used in pre-clinical studies as their scars and VT dynamics are known to be similar to humans [61]. Furthermore, the ability of models to reliably capture specific cardiac macro-and microstructures is limited by the resolution of current MRI scanners. Although the data used here has a higher (1 mm voxel) resolution than typical clinical MRI resolution, it still limits the accuracy of the reconstructed BZ structural features. Nevertheless, clinical studies have reported critical VT channels with a mean diameter of 7.4 mm [53] in swine and 9.5 mm in patients [62] which are above the resolution used here.
In addition, the VT exit site, taken as the endocardial breakthrough during a stable reentry, was considered here as the "reference solution" for the analysis. However, the exit site in the LV model as well as in many clinical cases is likely 3D and not just on the surfaces of the heart. The aim of the present study was just to compare what is being done clinically. Clinical electroanatomical mapping has its limitations but addressing them is out of the scope of this study. Also, while the maps constructed in Figs. 8-9 could be used to guide ablation, there are a number of uncertainties in accurate simulation of the biophysics of radiofrequency thermal lesion creation within the model, along with accurate representation of the realities of clinical ablation delivery which make faithfully modelling this process challenging.
State-of-the-art models, such as the one used in our study, only accounts for patient-(animal)specific ventricular geometry, not electrophysiology. A complete parametrization of biophysically detailed models would require acquisition of invasive data. Thus, personalization/parametrization of models is a topic of its own merits. Here, our goal was not to produce a like-for-like "patient specific" model but a generic model to highlight key mechanisms. Specifically, to address the clinical question of whether typical substrate mapping approaches can identify the VT exit site.
Another limitation is the use of a Monodomain instead of a Bidomain model. A full Bidomain model would require the 3D heart model to sit within a perfusing bath surrounding the tissue and filling the cavity, which is computationally significant. The approach chosen here is computationally cheaper and has been shown to yield similar results of a full-blown Bidomain solve [63].
Moreover UEs were computed assuming an infinitely small electrode tip. While electrodes with very small diameter (̃18 μm) exist [64], electrode sizes routinely used in the clinical setting are larger (millimeter scale). Electrodes with larger tip sizes may hamper the identification of ATs in fractionated signals due to attenuation effects [55]. Further, the amplitude of the BEs depends on interelectrode distance, which was not fixed (3 829 ± 751 μm). Instead, the largest amplitude amongst three BEs sharing the same triangular surface was considered here (see Fig. 2). Although this choice might have affected BVs in some regions, the known directional sensitivity of BEs was lessened by this approach. Furthermore, cut-offs were applied here to UEs rather than BEs as the former have been shown to be more useful than the latter at detecting intramural scars [51]. Finally, ATs and RTs in this study were obtained from noise-free UEs. RTs might be particularly challenging to be extracted from clinically-measured signals, which could hinder the performance of the RVI. However, our previous study demonstrated that even in the presence of noisy AT and RT maps, vulnerable regions could still be localized by the RVI [18].

Conclusions
In this study we used a structurally detailed 3D model to investigate clinically relevant features of signals used to guide catheter ablation of post-infarction VTs. Simulation results within the infarcted LV showed that: 1) BEs may better delineate scar anatomy than UEs; 2) voltage cut-offs alone are not capable to detect the VTs isthmus; 3) electrogram fractionation correlates with scar transmurality; 4) among all metrics, the RVI, when used alone (i.e., without applying the cut-offs), could identify regions closest to VT exit site, but it was slightly outperformed by the AT gradients combined with voltage cut-offs; This finding suggests that activation-repolarization metrics may help improve ablation therapy without the nuisances and risks of inducing VT even in infarcts with non-transmural scars. However, all electroanatomical maps were shown to be affected by the direction of the wavefront, which should be taken into account by VT substrate mapping strategies.