Granger Causality–Based Analysis for Classification of Fibrillation Mechanisms and Localization of Rotational Drivers

Supplemental Digital Content is available in the text.


Non-standard Abbreviations and Acronyms
The underlying mechanisms sustaining myocardial fibrillation remain unclear and there is a lack of consensus on a unifying mechanism. The anarchical hypothesis of multiple self-propagating wavelets and the opposing hierarchal hypothesis of organised spiral wave re-entry, referred to as rotors or rotational drivers (RDs), organised around a non-anatomical unexcited core termed a phase singularity (PS) point, continue to be debated as the mechanism sustaining fibrillation.
Although evidence exists from high spatiotemporal resolution optical mapping studies to support these conflicting mechanisms [1][2][3][4] , clinical translation to guide treatment remains challenging due to the poor spatiotemporal resolution of clinical mapping systems and a lack of suitable analysis tools. Thus, outcomes from catheter ablation remain poor, especially in persistent atrial fibrillation (AF), where recurrence remains as high as 40-50% 5 .
Multiple prominent investigators have proposed that ablation of sites localising RDs can terminate AF 6 and prevent ventricular fibrillation (VF) reoccurrence in survivors of sudden cardiac death 7 . Others have found no evidence for their existence in fibrillation mapping studies during cardiac surgery with high resolution electrode arrays 8 and proposed complex asynchronous endocardial-epicardial disassociation of fibrillatory conduction as a further mechanism.
A challenge in mapping fibrillatory mechanism is its dynamic nature. There is beat-tobeat variability in periodicity and amplitude of signals, and global wavefront propagation is nonuniform and temporally variable. RDs often demonstrate meandering and transient trajectories 9 .
Conventional mapping techniques, such as activation mapping which require annotation of a reference signal 10 and stable linear propagation are poorly suited. To overcome this issue, phase analysis is frequently used in fibrillation mapping instead. Phase analysis assigns a phase value between π and -π to the activation-recovery cycles of a given myocardial area 11 . Tracking these changes in phase allow for annotation of propagating wavefronts and localisation of PSs as areas devoid of a definitive phase. Accurate phase analysis however requires both global panoramic mapping and adequate spatial resolution 12 , which is currently not possible with clinicallyavailable tools.
We previously demonstrated that RD localisation from phase analysis of intracardiac electrograms (EGMs) acquired by multipolar catheters, including the 64-electrode basket catheters is inaccurate as it lacks sufficient resolution, requires significant interpolation, stitching of sequentially acquired data and is prone to generation of a high number of false positive RDs 13 . Other research groups have also highlighted the limitations of phase analysis with low spatial resolution 14 and demonstrated an increase in false positive RD detection rate with increasing inter-electrode distance and noise 15 . In addition, basket catheters can provide incomplete surface coverage and are susceptible to poor contact 16 . Alternatively, global cardiac mapping noninvasively is possible from the body surface using a multielectrode electrocardiographic imaging (ECGI) vest. ECGI utilises inverse solution algorithms for interpolating intracardiac EGMs.
Whilst there are some limitations with ECGI, such as correlation between surface and contact EGMs 17 and noise artefact, it does provide higher resolution mapping and has shown potential in a number of early studies where it has been utilised to map AF 18 , VF 7 and VT 19 mechanisms.
Here, we propose that granger causality (GC) analysis, originally an econometric tool designed to determine causal relationships between complex time series data 20 , can be repurposed as a novel tool to analyse fibrillation. Given GC analysis depends only on neighbouring causal relationships, we postulated that it could overcome the limitation of spatial resolution and sequentially acquired limited coverage data in fibrillation. We generated a range of fibrillation mechanisms in a rat VF model by modulating gap junction coupling and fibrosis, two factors implicated in cardiac remodelling. We hypothesised that GC-based analysis can be used to: (a) analyse temporal dependence of fibrillatory signals in neighbouring regions and determine the dominant propagating pattern, (b) quantify the global organisation and general mechanism of fibrillation, and (c) map stable RDs at low spatial resolution with limited coverage. Initially developed and validated against high-resolution phase analysis in a rat VF model, these novel GC-based analysis tools were further tested in previous VF optical mapping recordings of coronary perfused donor heart LV wedge preparations and finally adapted to analyse multi-electrode catheter recordings of persistent AF patients.

Methods
The data, analytic methods, and study materials are available from the corresponding author to other researchers for the purposes of reproducing the results or replicating the procedure upon reasonable request. Methods are described briefly here. For full details, please see the supplemental materials.

Ethical Approval
The animal work was performed in accordance with standards set out in the United Kingdom

Experimental Protocols
Eighteen Sprague-Dawley rats (250-300g) were humanely killed and the hearts were explanted, heparinised, and rapidly perfused ex-vivo on a Langendorff apparatus with Krebs-Henseleit solution, and stabilised for a 15 minute period before ex vivo optical mapping studies of transmembrane potential. In order to create a range of VF activity, in group one 8 of the hearts were acutely perfused with a gap junction uncoupler, carbenoxolone (0-50µM), which in our previous experiments produced increasing disorganised VF at increasing doses. In group two, the other 10 hearts had chronic patchy ventricular fibrosis which had been induced with ischaemiareperfusion cardiac surgery four weeks prior to the experiment (Supplementary Figure 1). No drugs were added to the perfusate, and in our previous experiments patchy fibrosis was found to sustain a more organised form of VF. Programmed electrical stimulation (PES) using a burst pacing protocol with the aid of Pinacidil 30µM was used to induce and sustain VF.

Optical mapping
Explanted hearts underwent optical mapping of the left ventricular anterior wall after VF induction. The transmembrane voltage was recorded from optical mapping fluorescence data using our custom made complementary metal oxide semiconductor (CMOS) camera (Cairns, Feversham UK) utilising the potentiometric dye RH237 (25µl of 1mg/ml dimethyl sulfoxide [DMSO]; Thermo-Fisher, Massachusetts, USA) and excitation-contraction uncoupler blebbistatin (10µmol/L, Tocris Bio-Sciences, Cambridge UK) in 160 x 128 pixel resolution for a 10 second duration. All our methods for filtering and analysing optical mapping fluorescence data have been previously described in detail 22,23 .

Organisational analysis
The processed optical mapping data were firstly analysed to quantify the degree of global organisation with two novel independent methods; causality pairing index (CPI) derived from global GC analysis and frequency dominance index (FDI) derived from the dominant frequency (DF) analysis. These indices are described below in detail.

Granger causality analysis
GC analysis is an econometric methodology for quantifying the causal dependence between two or more complex time series using a linear autoregressive model 20,24,25 . GC is a concept based on statistical prediction, whereby GC analysis statistically tests if a given time series signal A causes time series signal B, by analysing information contained in past values of time series A and determining if the information contained in signal A can predict signal B, beyond predictions from past values of signal B alone. In this work, we developed novel GC-based tools for fibrillation analysis adapted for low spatial resolution and limited coverage sequential mapping from optical fluorescence data of transmembrane potentials in a perfused rat VF model. GC analysis was used to measure the strength of the causal relationship between signals in neighbouring regions, and to quantify whether the fibrillation signal in one region over time could predict signals in another. GC vector mapping (described below) was performed based on quantified strengths of these relationships to determine the dominant propagating patterns.

Causality pairing index
CPI calculates global organisation of fibrillation from GC analysis, as described above. In this study, the temporal dependence structure between signals from different pixels was calculated by fitting a vector auto-regression model to a multi-variate signal. Thereafter, the CPI was measured by quantifying the percentage of possible pixel pairings between which there are propagational effects on a normalised scale of 0 to 1, where 0 is defined as no possible pairing having causal dependency and 1 where all possible pairings have causal dependency ( Figure 1A, 1B). The more pixel pairings that have a propagational effect between them above a specified threshold, the greater the level of global organisation in fibrillation.

Frequency dominance index
The FDI calculates the total level of global organisation by analysing all the DFs from all the signals within a fibrillating ventricle. The FDI is defined as the proportion of area in fibrillation occupied by the highest amplitude DF in the global spectrum relative to the area of all frequencies in the mapped area ( Figure 1C). The methodology for calculating DF has been previously described in detail 26 .

Down-sampling spatial resolution
After quantifying the level of global organisation in the VF data with CPI and FDI, we tested whether these indices were reliable and adaptable to lower spatial resolutions. CPI and FDI were benchmarked at decreasing spatial resolution (with data down-sampling) against two measures of global fibrillatory organisation, lps and nr/lr, generated from full spatial resolution phase processed data.

Phase Mapping
Rat VF phase processed analysis from full spatial resolution optical mapping was used for fully characterising the underlying mechanism of fibrillation, and for benchmarking our novel GCbased tools. Our methods for phase analysis and tracking of RDs have previously been described in detail 26,27 . A phase map of VF at each sampled time point was constructed and PS tagged using our algorithm. The edge of each wavefront was tracked in a 9x9 pixel window and maximum number of rotations [max (nr)] calculated. A minimum 2-rotation filter was used to threshold and define a significant RD and to construct RD heat maps from full spatial resolution data for validating our methodology for fibrillation analysis.

Phase characterisation of organisation and stability
From phase processed fibrillatory data, rotational activity was quantified by our metrics of organisation and stability (nr/lr and lps) and compared with CPI, FDI and a more widely used analysis feature in fibrillation literature, Shannon entropy (Shen) 28 . PS with < 2 were labelled 'non-significant' PSs and the number of locations (pixels) they occupied (lps) acted as a measure of global disorganisation, spatial meander and instability, whereby a large number of short lived meandering PSs would generate the highest value by this metric. PSs with ≥2 rotations were labelled 'significant' rotations drivers (nr) and the number of locations (pixels) they occupied (lr) over a fibrillatory recording was tracked, thus nr divided by lr acted as a measure of stability and global organisation, whereby RDs with high number of rotations localising to a small area would generate the highest values of this metric. These objective measures of fibrillation organisation (nr/lr and lps) calculated from high resolution phase analysis were correlated with our novel low resolution measures of fibrillation organisation (CPI and FDI). In addition, we correlated a conventional fibrillation analysis tool, Shannon entropy (Shen) 26,28 , with CPI and FDI.

GC Vector Maps
GC vector mapping was tested as a low spatial resolution tool adapted specifically for limited spatial coverage to determine dominant propagation and localise RDs independent of phase analysis. In this work, GC vector maps were generated from an 8 x 8 data grid at 25% of the full spatial resolution optical fluorescence data, whereby three pixels were discarded between each data point. Within this 8 x 8 data grid, an algorithm determined firstly if there were any signals with temporal causal dependence, thereafter it quantified the strength of the temporal causal dependence between these signals. A vector was then plotted only between the source signal and the signal where it exerted the greatest causal or propagational effect above a specified threshold ( Figure 1A). GC vector maps were benchmarked against full resolution phase processed data to determine ability to characterise fibrillation mechanism. Thus by applying the equation below to these 3D vectors outputs the circular interdependence value can be calculated.
CIV is quantified on a normalised scale of 0 to 1, whereby a stable RD present throughout the recording would generate a value of 1, and disorganised or propagating wavefronts in one direction will generate a value closer to 0.

Ex vivo human VF mapping
Our GC-based analysis tools were tested on human VF optical mapping data that was previously acquired as part of a separate study by Aras et al. (2019) and the methodology was reported in detail 21 . Briefly, we tested our GC-based analysis tools on 33 VF recordings from 12 representative de-identified human donor hearts. These recordings were 4-seconds in duration and taken from coronary perfused LV wedge preparations that had VF induced with 25µM pinacidil pre-treatment. The mean LV wedge dimensions were 7 cm x 3.5 cm x 1.8 cm (height x width x thickness).

Clinical in vivo AF mapping
In 16 patients presenting with symptomatic persistent AF for a first ablation procedure, EGMs  (Figure 2A, 2B). Disorganised fibrillation with a high number of meandering non-significant PSs had low FDI and CPI values. Conversely, nr/lr, a measure of global fibrillatory organisation that tracks presence of spatiotemporally stable RD from full resolution phase processed data, correlated positively with CPI and FDI at decreasing resolution ( Supplementary Figure 2A, 2B). Shannon entropy values, more conventionally utilised in fibrillation analysis, showed no statistically significant correlation with lps or nr/lr ( Figure 2C, Supplementary Figure 2C).
After determining the applicability of our novel organisational indices to low resolution data, we selected representative hearts along the organisational spectrum, to delineate whether the fibrillatory mechanism characterised by full resolution phase analysis correlated with the level of global organisation as determined with CPI and FDI with low spatial resolution data.
Hearts with the highest FDI had the most spatiotemporally stable RDs (Figure 3A, 3B), with high numbers of rotations and much fewer short-lived PSs in comparison to hearts with a lower FDI ( Figure 3C). The RD heat map showed well localised discrete areas harbouring the organised RD in hearts with a high FDI, whereas hearts with a low to intermediate FDI values did not show such areas ( Figure 3A).

GC mapping to localise driver regions
After establishing that high global organisation in fibrillation, as measured by CPI and FDI, correlated positively with the existence of localised stable RDs, we tested whether GC vector mapping could identify causally dependent neighbouring regions to localise areas harbouring RDs at low spatial resolutions independent of phase analysis. In a representative heart with high global fibrillatory organisation, optical fluorescence data was down-sampled to 25% of full spatial resolution for GC vector mapping. GC vector mapping identified a number of regions with causal dependence and localised an area harbouring a spatiotemporally stable RD with a circular interdependence of GC-vectors as shown in Figure 4A, 4B. Optical fluorescence of the transmembrane potentials along this region showed repetitive sequential activation over time ( Figure 4C). At 12.5% of the full spatial resolution RD regions could no longer be localised accurately with GC mapping (results not enclosed).
We further validated RD localisation with GC-vector mapping performed at 25% spatial resolution against full spatial resolution phase analysis in a further three hearts classified as organised from high FDI and CPI measures.  Figure 5A). In regions harbouring stable RDs, the CIV was 2.6 fold higher than regions without RD (0.91±0.05 versus 0.35±0.06, p=0.0002, n=3, Figure 5B).

GC-based analysis of human VF
We further tested and validated our low spatial resolution adapted GC-based tools in optical mapping of human VF in LV wedge preparations to test applicability to a larger spatial scale and benchmarked these tools against full resolution phase-analysis. As with rat VF, a spectrum of VF mechanisms were found, ranging from fibrillation driven by organised RDs to completely chaotic activity. At decreasing spatial resolutions, lps correlated negatively with CPI as before; 50% spatial resolution: F(1,10)=24.4, p<0.0001, R 2 =0.42 and 25% spatial resolution: F(1,10) =21.3, p=0.0001, R 2 =0.38 ( Figure 6A). As before, the optical fluorescence data was downsampled to 25% of full spatial resolution for GC vector mapping. Figure 6B shows two representative hearts with low and high global organisation respectively. In a representative heart with a high global fibrillatory organisation, as quantified by CPI, GC vector mapping identified a regions with high CIV, harbouring a stable RD and this correlated with the same region identified by phase analysis on the RD heatmap. Similarly, in a representative heart with low global fibrillatory organisation, as quantified by CPI, GC vector mapping showed random vector distribution with low CIV values and this correlated with multiple wavelet driven fibrillation on phase analysis with no stable RDs.

GC vector mapping of intracardiac EGMs from AF mapping
After developing and validating these novel GC-based tools for use with low resolution data from optical mapping of transmembrane voltage, we adapted our methodology to process intracardiac EGMs in human persistent AF acquired sequentially using a 20-pole AFocusII TM mapping catheter. We measured global fibrillatory organisation and localised RDs in persistent AF with our novel indices. Processing EGMs for GC analysis requires different considerations for signal processing than optical fluorescence data. Firstly, a three-dimensional spatial map was constructed from spatial correlates of the electrodes and corresponding bipoles. Each signal underwent high and low band pass filtering before been rectified and down-sampled. Causal dependence between bipoles was established from rectified down-sampled data ( Figure 7A). representative paced data as expected is low.

GC vector mapping quantifies global AF organisation and identifies areas harbouring RDs
The EGM recordings during AF mapping were processed and causality maps plotted from 20seconds of recording. The threshold value for a RD positive site was determined from a receiver operating characteristic (ROC) curve (Supplementary Figure 3). Figure 8A Figure 8E]. In patients with low CPI fibrillation the underlying mechanism was chaotic activity with no RDs, whereas in patients with higher CPI values there was a higher number of RDs.

Discussion
In this study, we developed novel GC-based tools for fibrillation analysis adapted for use with low spatial resolution data acquired with limited coverage sequential mapping, independent of The clinically available fibrillation analysis tools, most widely available of which is phase analysis, considered to be the "gold-standard", have significant limitations. In AF simulation studies, we previously demonstrated that the spatial resolution of multiple commercially available clinical catheters including the AFocusII TM used in this study is prone to generating a large number of false positive RD detections from phase processed data and is ineffective at locating RDs with meandering trajectories. In perfused heart AF mapping studies, phantom RDs often appear as spatial resolution is reduced 38 . Phase mapping additionally requires careful consideration of a number parametric inputs; such as rotational thresholds for labelling RDs and average fibrillation cycle lengths for accurate analysis 39 .
GC-based analysis negates the issues generated by low spatial resolution phase analysis and dynamic nature of fibrillation by considering only causal interdependence of neighbouring areas over time rather than attempting to construct panoramic videos of wavefront dynamics within a chamber from interpolated, temporally stitched and sequentially acquired data to determine mechanism. Furthermore, by determining only causal relationships in neighbouring regions over long periods of mapping in fibrillation, where activation patterns are dynamic and vary beat to beat, GC mapping also addresses a key limitation of activation mapping in fibrillation, which is annotation of relative activation times from complex electrograms. GC analysis was initially developed as an econometric tool, it has been used to characterise the causal connectivity between different brain regions based on neuroimaging data 40,41 and is widely used in climate science to establish causal relationships between two parameters 42 .
By quantifying the circular interdependence of GC-vectors that characterises a RD with CIV through measurements of the directionality of cross-product vectors relative to the centre, we provide a non-subjective measure for localising RDs that does not rely on user interpretation.
Windowing CIV in time segments can quantify the temporal stability of RDs, accounting for their transient and meandering nature and determine their overall contribution to the fibrillation mechanism.
We found a direct correlation between global fibrillatory organisation and the underlying mechanism in both rat VF and human VF, and therefore classifying this organisation with FDI or CPI in low spatial resolution without detailed mapping may help classify the predominant mechanism and guide treatment strategy. For instance, the operator may choose to pursue detailed mapping and RD ablation in patients with high FDI and CPI values only. Some clinical studies support the existence of a spectrum of organisation and mechanisms in fibrillation. For instance in VF mapping studies of patient undergoing cardiac surgery stable RDs, meandering RDs and multiple wavelets were all found as predominant mechanisms in different patients 43 .
Similarly in AF, non-invasive mapping with ECGI has shown coexistence of a number of mechanisms and varying fibrillation complexity 18 . The global fibrillatory organisation and incidence of RDs were low in many of our persistent AF mapped patients. This suggests that mechanism guided ablation may only be suitable in a select number of persistent AF patients. This study has a few important limitations. GC-vector mapping was used to analyse twodimensional fibrillatory data and may not reflect the transmural propagation in fibrillation.
EGMs show sharp deflection in fibrillation, rather than a sinusoidal waveforms seen from optical fluorescence. GC-based analysis is more dependent on determining causal relationships over time than timing of local activation, further work is needed to evaluate whether GC vector maps from EGMs and optical fluorescence are comparable. GC based analysis is used to determine both global fibrillatory organisation with CPI and to localise RDs with GC vector maps, thus, these methods are not mutually exclusive and may influence each other. One of the limitations of human VF mapping in this study was that the analysis was not performed in a whole intact ex vivo perfused heart. Nevertheless, the volume of an LV wedge preparation was clearly sufficient to sustain VF, and the volume of myocardium in the wedge preparation greater that the wavelength volume of the fibrillation 21 and therefore a suitable preparation to study fibrillatory dynamics and to validate of the GC-based analysis tools in large hearts.

Conclusion
In summary, we present novel methodologies based on GC analysis for measuring global fibrillatory organisation and mapping RDs. The techniques presented here are optimised for sequential mapping with limited spatial resolution and coverage, and were developed and

determination -R 2 and p value is indicated, n =16)
What is Known?
• Pre-clinical studies have implicated multiple competing mechanisms for sustaining myocardial fibrillation.
• Clinical translation to guide treatment in patients with AF and VF survivors remains challenging due to the poor spatial resolution of clinical mapping systems and a lack of suitable analysis tools.

What this Study Adds?
• Granger causality (GC) analysis, originally an econometric tool for quantifying causal relationships between complex time-series, was developed in rat VF, and validated in human VF and AF as a novel fibrillation mapping tool.
• GC-based fibrillation analysis can measure global fibrillation organisation, characterise dominant propagating patterns and map rotational drivers using low spatial resolution sequentially-acquired data.