Cardiac re-entry dynamics & self-termination in DT-MRI based model of Human Foetal Heart

TL;DR

This study uses anatomically realistic DT-MRI based models of the human fetal heart to investigate how heart geometry and fiber anisotropy influence cardiac re-entry behavior and self-termination, revealing significant effects of anisotropy.

Contribution

It introduces a novel simulation approach using DT-MRI data to analyze re-entry dynamics in fetal hearts, highlighting the impact of fiber anisotropy on re-entry stability and termination.

Findings

Fiber anisotropy changes re-entry from pinned to anatomical in 2D slices.

In 3D models, anisotropy leads to re-entry self-termination.

Heart geometry influences re-entry dynamics significantly.

Abstract

The effect of heart geometry and anisotropy on cardiac re-entry dynamics and self-termination is studied here in anatomically realistic computer simulations of human foetal heart. 20 weeks of gestational age human foetal heart isotropic and anisotropic anatomy models from diffusion tensor MRI data sets are used in the computer simulations. The fibre orientation angles of the heart were obtained from the DT-MRI primary eigenvalues. In a spatially homogeneous electrophysiological mono domain model with the DT-MRI based heart geometries, we initiate simplified Fitz-Hugh-Nagumo kinetics cardiac re-entry at a prescribed location in a 2D slice, and in the full 3D anatomy model. In a slice of the heart, the MRI based fibre anisotropy changes the re-entry dynamics from pinned to anatomical re-entry. In the full 3D MRI based model, the foetal heart fibre anisotropy changes the re-entry dynamics from a persistent re-entry to the re-entry self-termination.

Equations

$$\frac{\partial \mathbf{u}}{\partial t}=\mathbf{f}(\mathbf{u})+\nabla \cdot \hat{\mathbf{D}}\nabla \mathbf{u},$$

$$\displaystyle\qquad\qquad{{\vec{n}}\cdot\hat{\mathbf{D}}\nabla\mathbf{u}}\bigg{|}_{G}=0,$$

$$f(u,v)$$

$$g(u,v)$$

$$P_{jk}=P_{\perp }{\delta }_{jk}+\left(P_{\parallel }-P_{\perp }\right)f_j{f}_k,$$

Full text

††preprint: AIP/123-QED

The effect of the heart anisotropy and anatomy on cardiac re-entry dynamics, although difficult to demonstrate in experiment, is well appreciated Bishop et al. (2010); Bishop, Vigmond, and Plank (2011); Bishop and Plank (2012); Fukumoto et al. (2016), and has been studied in simplified mathematical and computer models Fenton and Karma (1998); Pertsov et al. (2000); Wellner et al. (2002); Rodriguez, Eason, and Trayanova (2006); Dierckx et al. (2013). The BeatBox Antonioletti et al. (2017) High Performance Computing (HPC) cardiac electrophysiology computer simulation environment allows direct incorporation of the high resolution DT-MRI heart anatomy data sets into the biophysically and anatomically realistic computer simulations. In the BeatBox in-silico model, the anisotropy of the tissue is switched “on” and “off” to allow for comparison between the anatomically realistic isotropic and anisotropic conduction, in order to see the specific pure anatomy effects, as well as the interplay between the anisotropy and anatomy of an individual heart. In this paper, we present the DT-MRI based anatomy and myofiber structure realistic computer simulation study of cardiac re-entry dynamics in the in-silico model of the human foetal heart Pervolaraki et al. (2013). We demonstrate that, in a 2D slice of the heart, the realistic fiber anisotropy of the tissue changes cardiac re-entry dynamics from pinned into fast anatomical re-entry. In the full 3D DT-MRI based model, depending on the initial location of the re-entry, the isotropic geometry of the heart might sustain a perpetual re-entry even with a positive filament tension; while the same positive filament tension re-entry initiated at the same location of the foetal heart with the realistic fiber anisotropy self-terminates within seconds. Generally, time of re-entry self-termination depends on the re-entry initial position, while the role of the heart anisotropy is to speed up the re-entry self-termination.

I Introduction

Since the hypothesis over a century ago that cardiac re-entry underlies cardiac arrhythmias Mines (1913); Garey (1914) , and the much later confirmation of the hypothesis in cardiac tissue experiment Allessie, Bonke, and Schopman (1973); Pertsov et al. (1993), the re-entry (aka spiral wave in 2D, cardiac excitation vortex in 3D), its origin and its role in sustained arrhythmias and fibrillation, as well as a possibility of its effective control and defibrillation, have been an object of extensive theoretical study and modelling Wiener and Rosenblueth (1946); Balakhovsky (1965); Krinsky (1968); Panfilov, Rudenko, and Pertsov (1984); Davydov et al. (1988); Keener (1988); Ermakova, Pertsov, and Shnol (1989); Biktashev and Holden (1994); Biktashev (1998); Fenton and Karma (1998); Pertsov et al. (2000); Wellner et al. (2002); Biktasheva and Biktashev (2003); Biktashev, Barkley, and Biktasheva (2010); Biktashev, Biktasheva, and Sarvazyan (2011); Biktasheva, Dierckx, and Biktashev (2015). From experiment, it is an established point of view that cardiac arrhythmias are due to a complex combination of electrophysiological Bosch and Nattel (2002); Workman, Kane, and Rankin (2008); Kushiyama et al. (2016), structural Pellman, Lyon, and Sheikh (2010); Eckstein et al. (2011); Takemoto et al. (2012); Eckstein et al. (2013), and anatomical MacEdo et al. (2010); Anselmino et al. (2011) factors which sustain cardiac re-entry Gray, Pertsov, and Jalife (1996); Wu et al. (1998); Nattel (2002); Yamazaki et al. (2012).

The specific effect of the heart anisotropy and anatomy on cardiac re-entry dynamics is well appreciated Bishop et al. (2010); Bishop, Vigmond, and Plank (2011); Bishop and Plank (2012); Fukumoto et al. (2016), and has been studied in simplified mathematical and computer models Fenton and Karma (1998); Pertsov et al. (2000); Wellner et al. (2002); Rodriguez, Eason, and Trayanova (2006); Dierckx et al. (2013).

The anisotropic discontinuities in the heart muscle have been commonly seen as a substrate for rise of cardiac re-entry due to the abrupt change in conduction velocity and wavefront curvatureFenton and Karma (1998); Spach (2001); Smaill et al. (2004). On the other hand, extensive mapping of cardiac myocyte orientation in mammalian hearts has shown that the transmural fiber arrangement, including the range of transmural change in fiber angle in ventricular wall, was consistent within a species, and varied between species(Hunter et al., 1997, p. 173). So that changes in anisotropy seen in healthy hearts can facilitate initiation of arryhthmias.

The recent advance in DT-MRI technology and High Performance Computing (HPC) allows the obtained DT-MRI data sets with the detailed heart anatomy and myofiber structure to be directly incorporated into the anatomically realistic computer simulations Antonioletti et al. (2017), so that the anisotropy of the tissue in the in-silico model can be switched on and off to allow for comparison between the anatomically realistic isotropic and anisotropic conduction in order to see specific anatomy effects as well as the interplay between the anisotropy and anatomy of an individual heart.

In this paper, we present the raw DT-MRI based anatomically and myofiber structure realistic computer simulation study of cardiac re-entry dynamics in the in-silico model of human foetal heart.

The raw DT-MRI image data Pervolaraki et al. (2013) was segmented into the tissue/non-tissue pixels based on the MRI luminosity threshold, followed by the calculation of the fiber angles at each voxel from the diffusion-weighted DT-MRI images. This very basic segmentation might be seen as a limitation of the study from the cardiac physiology point of view. However, the purpose of our study is not to provide results of immediate physiological or clinical relevance: for these we currently simply have not enough data. Rather, from the non-linear science point of view our rationale is to use the raw DT-MRI data “as is” as an example of an unaltered nature provided medium to study a re-entry dynamics. Although the DT-MRI yields three eigenvalues, the second and the third are often harder to distinguish, so we used only the primary eigenvector to define the local fibre orientations in the simulation study. The focus of the paper is to demonstrate the effect of a real mammalian heart anatomy and anisotropy on a re-entry dynamics. The available MRI data of a foetal heart provide an excellent oportunity for such study. The objectives for the use of foetal heart MRI data are: whether the anatomical settings of the although foetal but a real heart might support a positive filament tension re-entry, and what would it be the role of a real heart anisotropy in that case. So, here we demonstrate that the real heart anisotropy enhances re-entry self-termination.

We demonstrate that, in a 2D slice of the heart, the realistic fiber anisotropy might change the re-entry dynamics from pinned to anatomical re-entry.

In the full 3D DT-MRI based model, depending on the location of the re-entry initiation, the isotropic geometry of the heart might sustain perpetual re-entry even with a positive filament tension kinetics. While the same positive filament tension re-entry initiated at the same location of the foetal heart with the realistic fiber anisotropy self-terminates within seconds. Time of re-entry self-termination depends on the re-entry initial position. Anisotropy of the real heart speeds up re-entry self-termination, and in this sense has a rather anti-arrhythmogenic effect. The geometry and anisotropy of the heart together ensure the fastest self-termination of cardiac re-entry.

The novel significance of our findings is that we demonstarte that the real life heart anisotropy might have a rather anti-arrhythmic function as it facilitates fast self-termination of cardiac re-entry.

II Methods

II.1 DT-MRI based anatomy model

The DT-MRI data sets of the $128\times 128\times 128$ voxels size, with voxel resolution of $\sim 100\mu m$, of ethically obtained 143 days of gestational age (DGA) human foetal heart Pervolaraki et al. (2013), were converted into the BeatBox Antonioletti et al. (2017) regular Cartesian mesh .bbg geometry format, containing the DT-MRI cartesian coordinates of the heart tissue points together with the corresponding components of the diffusion tensor primary eigenvectors Antonioletti et al. (2017). The .bbg file is an ASCII text file, each line in which describes a point in a regular mesh in the following format:

[TABLE]

Here x, y, z are integer Cartesian coordinates of a DT-MRI voxel, status is a flag with a nonzero-value for a tissue point, and fibre_x, fibre_y, fibre_z are $x$-, $y$- and $z$-components of the fibre orientation vector at that point. To reduce the size of the .bbg files, only the tissue points, that is points with nonzero status need to be specified, because the BeatBox solver will ignore the void points with zero status in any case. Although the original DT-MRI images data sets had $128\times 128\times 128$ voxels size, the actual dimensions of the foetal heart minimum bounding box were $67\times 91\times 128$, with $181070$ tissue points.

The raw DT-MRI anatomy data Pervolaraki et al. (2013) were segmented into the “tissue”/“not tissue” pixels discretion based on the MRI luminosity threshold, with the cartesian fiber angles at each voxel obtained from the diffusion-weighted DT-MRI images. Only this basic segmentaion of the raw DT-MRI anatomy data Pervolaraki et al. (2013) was taken into account in the computer simulation of cardiac re-entry dynamics, so we shall refer to it as the raw DT-MRI based anatomy model.

In the 2D model, the fibres vectors were projected into the plane, in order to construct the 2D diffusivity tensor.

Fig. 1 shows the cross section of the 143 days of gestational age (DGA) foetal heart with already formed intramural laminar structure and more irregular epicardial, endocardial, and septal fibers, see Figure 4 of Pervolaraki et al (Pervolaraki et al., 2013, p. 5) for the color-encoded fractional anisotropy (FA) and all the three components of the fiber angles in human foetal hearts. The DT-MRI based foetal heart model offered a unique opportunity to see if the 20 weeks of gestation age intramural heart structure was capable to support cardiac re-entry, as it would not be possible for the re-entry to pin to the endocardial fine features which were yet to be developed later, such as e.g. the pinning to pectinate muscles junction with crystae terminalis reported in adult human atria Wu et al. (1998); Yamazaki et al. (2012); Kharche et al. (2015a).

II.2 Cardiac Tissue Model

To investigate the effects of anatomy on cardiac re-entry dynamics we used monodomain tissue model with non-flux boundary conditions

[TABLE]

where $\mathbf{u}({{\vec{r}}},t)={(u,v)}^{T}$, ${{\vec{r}}}$ is the position vector, ${\mathbf{f}}({{\vec{r}}},t)={(f,g)}^{T}$ is the Fitz-Hugh-Nagumo Winfree (1991) kinetics column-vector

[TABLE]

with the parameter values $\alpha=0.3$, $\beta=0.71$, $\gamma=0.5$, which in an infinite excitable medium support a rigidly rotating vortex with positive filament tension Biktashev, Holden, and Zhang (1994). The simplified FHN model was intentionally chosen for this study in order to fully eliminate the possible effects of a realistic cell excitation kinetics, such as e.g. meander Winfree (1991), alternansKarma (1994), negative filament tensionBiktashev, Holden, and Zhang (1994), etc., and in order to enhance and highlight the pure effects of the heart anatomy and anisotropy on the cardiac re-entry outcome. $\hat{\mathbf{D}}=\mathbf{Q}\hat{P}$, where $\mathbf{Q}=\mathrm{diag}(1,0)=\begin{bmatrix}1&0\\ 0&0\end{bmatrix}$ is the matrix of the relative diffusion coefficients for $u$ and $v$ components, and $\hat{P}=[P_{{j}{k}}]\in\mathbb{R}^{3\times 3}$ is the $u$ component diffusion tensor, which has only two different eigenvalues: the bigger, simple eigenvalue ${P_{\parallel}}$ corresponding to the direction along the tissue fibers, and the smaller, double eigenvalue ${P_{\perp}}$, corresponding to the directions across the fibres, so that

[TABLE]

where ${\vec{{f}}}=\left({f}_{k}\right)$ is the unit vector of the fiber direction; ${\vec{n}}$ is the vector normal to the tissue boundary $G$. In the isotropic simulation, ${P_{\parallel}}$ and ${P_{\perp}}$ values were fixed at ${P_{\parallel}}={P_{\perp}}=1$ (corresponding 1D conduction velocity 1.89). In the anisotropic simulations, ${P_{\parallel}}$ and ${P_{\perp}}$ values were fixed at ${P_{\parallel}}=2$, ${P_{\perp}}=0.5$ (corresponding conduction velocities 2.68 and 1.34 respectively). All the conduction velocities have been computed for the period waves with the frequency of the free spiral wave in the model, i.e. 11.36. With the isotropic diffusivity (${P_{\parallel}}={P_{\perp}}=1$) equal to the geometric mean between the faster and the slower anisotropic diffusivities (${P_{\parallel}}=2,{P_{\perp}}=1/2$), the isotropic conduction velocity 1.89 was almost exactly the same as the geometric mean 1.89 of the faster and slower (2.68 and 1.34 respectively) anisotropic conduction velocities, chosen in order to minimize the maximal relative difference between the isotropic and anisotropic propagation speeds.

All the computer simulations presented here were done using the BeatBox Antonioletti et al. (2017) software package with the explicit time-step Euler scheme, on the Cartesian regular grid with space step discretization $\Delta x=0.1$, time step discretisation $\Delta t=0.001$; 5-point stencil for isotropic, and 9-point stencil for anisotropic Laplacian approximation in 2D simulations; 7-point stencil for isotropic, and 27-point stencil for anisotropic Laplacian approximation in 3D simulations. The re-entry was initiated by the phase distribution method Biktashev and Holden (1998): in the 2D simulations, at a prescribed location of the cross section of the DT-MRI based anatomical model; in the 3D simulations, at a prescribed location of the full DT-MRI based whole heart anatomical model.

III Results

III.1 2D MRI-based “slice” simulations

In the 2D simulations, Fig. 2, a counter-clockwise re-entry was initiated by the phase distribution method Biktashev and Holden (1998), with the initial center of rotation placed at the prescribed location $x_{0}=40,y_{0}=60$ in the 2D cross section of the DT-MRI based anatomical model shown in Fig. 1.

In the Fig. 2(a-b), it can be seen that in both isotropic and anisotropic 2D simulations, at $t=0$, there was identical location of the initial re-entry rotation center: roughly in the middle of the slice, in the vicinity of the septum cuneiform opening.

Fig. 2(a) shows the isotropic dynamics of the re-entry, that is with the fiber orientation data “turned OFF”, so that only the geometry of the isotropic homogeneous slice affects the dynamics of the re-entry. While it is known that in an infinite medium the chosen FHN kinetics parameter values $\alpha=0.3$, $\beta=0.71$, $\gamma=0.5$ produce rigidly rotating spiral Winfree (1991), the anatomically realistic boundaries of the foetal heart cause the drift of the re-entry. The re-entry does not terminate because of the resonant reflection from the inexcitable boundaries Biktashev and Holden (1994), but after the transient first rotation around the septum cuneiform opening, the tip of the re-entry firmly pins to the sharp lower end of the cuneiform opening, see Fig. 2(a).

Fig. 2(b) shows the anisotropic dynamics of the re-entry, that is with the fiber orientation data “turned ON”, so that both the anatomically realistic geometry and the anisotropy of the otherwise homogeneous slice of the heart affect the dynamics of the re-entry, causing its drift. In the anisotropic slice, the re-entry also does not terminate at the inexcitable boundaries, but after the faster than in the isotropic case, see the $a.u.$ time labels in the Fig. 2(a-b), transient first rotation around the septum cuneiform opening, the anatomically realistic anisotropy of the medium turns the initial spiral wave into the fast anatomical re-entry around the septum cuneiform opening, see Fig. 2(b).

III.2 3D Whole heart MRI-based simulations

In the 3D whole heart MRI-based simulations shown in the Fig. 3and Fig. 4, a counter-clockwise excitation vortex was initiated by the phase distribution method Biktashev and Holden (1998), with the initial position of the transmural vortex filament (yellow line) at the prescribed location along the $x$ axis at $y_{0}=40,z_{0}=60$. It can be seen in Fig. 3 isotropic, and Fig. 4 anisotropic 3D simulations that, at $t=0$, there was identical initial location of the filament of the excitation vortex: that is transmurally, roughly in the middle through the ventricles of the heart.

Fig. 3 shows the isotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned OFF”, so that only the geometry of the otherwise isotropic homogeneous foetal heart affects the dynamics of the vortex. It is known that the chosen FHN kinetics parameter values $\alpha=0.3$, $\beta=0.71$, $\gamma=0.5$ produce rigidly rotating vortex with the positive filament tension Biktashev, Holden, and Zhang (1994), which, depending on the topology, either collapses or straightens up between the opposite boundaries of the excitable medium. In the 3D anatomically realistic isotropic simulations of the foetal heart, the anatomically realistic boundaries of the heart cause drift of the excitation vortex, and, depending on the initial position of the vortex filament, vortices with the positive filament tension tend to collapse. However, there exist initial locations of the excitation vortex, which although result in the drift of the vortex, still do not lead to the expected collapse of the vortex with positive filament tension. One of such outcomes is shown in the Fig. 3. Here, following the geometry of the heart, after a very short transient, the initial vortex filament breaks into the two short pieces, each of which finds its own synchronous perpetual pathway in the “isotropic” foetal heart, resulting in the seemingly perpetual cardiac re-entry, which failed to self-terminate within the extended simulation time, see Fig. 3.

Fig. 4 shows the anisotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned ON”, so that both the anatomically realistic geometry and the anisotropy of the otherwise homogeneous foetal heart affect the dynamics of the initial vortex.

In the 3D whole heart MRI-based simulations shown in the Fig. 5

and Fig. 6,

a counter-clockwise excitation vortex was initiated by the phase distribution method Biktashev and Holden (1998), with the initial position of the transmural vortex filament (yellow line) at the prescribed location along the $y$ axis at $x_{0}=40,z_{0}=60$, that is perpendicular to the initial orientation of the vortex filament shown in Fig. 3 and Fig. 4. It can be seen in Fig. 5 isotropic, and in Fig. 6 anisotropic 3D simulations, that at $t=0$, there was identical intial location of the filament of the excitation vortex: that is transmurally, roughly in the middle through the ventricles of the foetal heart, and perpendicular to the initial orientation of the vortex filament shown in Fig. 3 and Fig. 4.

Fig. 5 shows the isotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned OFF”, so that only the geometry of the otherwise isotropic homogeneous foetal heart affects the dynamics of the vortex. Here, contrary to the expectation for the positive filament tension vortex to always contract, the organising filament first transiently extends intramurally along the tissue walls, before finally breaking up to the two ring-like pieces, each of which quickly contracts and terminates at the opposite base and apex regions of the heart.

Fig. 6 shows the anisotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned ON”, so that both the anatomically realistic geometry and the anisotropy of the otherwise homogeneous foetal heart affect the dynamics of the initial vortex leading to its really fast termination at the apex of the heart.

In the raw DT-MRI model simulations shown in Fig. 3, Fig. 4, Fig. 5, and Fig. 6, it can be seen that although the organising filament of the vortex could not get through into the accidental “leftover” piece of tissue adjacent to the apical region, the piece got activated and might have served as an artificial “capacitor” affecting dynamics of the re-entry. In order to check whether this might be the case, we edited the original raw DT-MRI model by removing in the MRI the foreign piece, and repeated the whole heart isotropic and anisotropic simulations from the same two orthogonal initial locations of the re-entry, similar to the shown in Fig. 3, Fig. 4, Fig. 5, and Fig. 6.

In the 3D whole heart “edited” MRI model simulations shown in the Fig. 7

and Fig. 8,

a counter-clockwise excitation vortex was initiated by the phase distribution method Biktashev and Holden (1998), with the initial position of the transmural vortex filament (yellow line) at the prescribed location along the $x$ axis at $y_{0}=40,z_{0}=60$. It can be seen in Fig. 7 isotropic, and in Fig. 8 anisotropic 3D simulations, that, at $t=0$, there was identical initial location of the filament of the excitation vortex: that is transmurally, roughly in the middle through the ventricles of the foetal heart, similar to the initial location of the vortex filament in the raw DT-MRI simulations shown in Fig. 3 and Fig. 4 .

Fig. 7 shows the isotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned OFF”, so that only the geometry of the otherwise isotropic homogeneous foetal heart affects the dynamics of the vortex. Here, following the geometry of the heart, the organising filament of the initial vortex also breaks into the two short pieces, each of which also finds its own synchronous pathway similar to the beginning of the raw DT-MRI simulation shown in Fig. 3. However, this time, after a few rotations, the two re-entries find their end in their also almost synchronous termination of the filaments in the base region of the foetal heart, see Fig. 7.

Fig. 8 shows the anisotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned ON”, so that both the anatomically realistic geometry and the anisotropy of the otherwise homogeneous foetal heart affect the dynamics of the initial vortex. Here, the anisotropy of the heart also causes the significant transient distortion of the organising filament of the initial vortex, followed by its fast drift towards the apex, and the ultimate termination at the AV border before a completion of a single rotation, very similar to the raw DT-MRI simulation shown in Fig. 4. However, this time, without the “leftover” piece “incidental capcitor” effect, there is just a bit faster repolarisation of the whole heart than it was in the presence of the “incidental capcitor” in the raw DT-MRI simulation shown in Fig. 4. For the comparison of the re-entry termination times, and the whole heart recovery times, see FIG. 11.

In the 3D whole heart “edited” MRI model simulations shown in the Fig. 9

and Fig. 10,

a counter-clockwise excitation vortex was initiated by the phase distribution method Biktashev and Holden (1998), with the initial position of the transmural vortex filament (yellow line) at the prescribed location along the $y$ axis at $x_{0}=40,z_{0}=60$. It can be seen in Fig. 9 isotropic, and in Fig. 10 anisotropic 3D simulations, that at $t=0$, there was the identical initial location of the filament of the excitation vortex: that is transmurally, roughly in the middle through the ventricles of the foetal heart, perpendicular to the initial location of the vortex filament in the “edited” MRI simulations shown in Fig. 7 and Fig. 8, and similar to the initial location of the vortex filament in the raw DT-MRI simulations shown in Fig. 5 and Fig. 6 .

Fig. 9 shows the isotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned OFF”, so that only the geometry of the otherwise isotropic homogeneous foetal heart affects the dynamics of the vortex. Here, again contrary to the expectation for a positive filament tension vortex to always contract, the organising filament first transiently extends intramurally before breaking up into the two ring-like pieces, each of which quickly contracts and terminates at the opposite base and apex regions of the heart, identical to what can be seen in the raw DT-MRI simulation shown in Fig. 5. So that, this time, for this particular orientation of the initial re-entry, the “leftover” tissue “incidental capacitor” effect does not seem to play any role in the outcomes of the isotropic “heart geometry only” raw DT-MRI simulations shown in Fig. 5 as opoosed to the outcome of the identical initial re-entry location in the “edited” MRI simulations shown in Fig. 9.

Fig. 10 shows the anisotropic dynamics of the excitation vortex, that is with the fiber orientation data “turned ON”, so that both the anatomically realistic geometry and the anisotropy of the otherwise homogeneous foetal heart affect the dynamics of the initial vortex, which, in the absence of the “incidental capacitor” effect, results in the fastest possible termination of the re-entry at the apex of the heart, before the vortex first rotation ever started. The re-entry termination time here is more than twice shorter than in the raw and “edited” MRI isotropic simulations shown in Fig. 5 and Fig. 9, shorter than in the analogous simulation with the “incidental capacitor” effect shown in the Fig. 6, and times shorter than in any of the simulations of the re-entry with the perpendicular initial location shown in the Fig. 3, Fig. 4, Fig. 7, and Fig. 8 .

In Fig 11, we have summarized the results of the simulations shown in Fig. 3, Fig. 4, Fig. 5, Fig. 6, Fig. 7, Fig. 8, Fig. 9, and Fig. 10, with the re-entry termination time shown in arbitrary time units under each respected whole heart model and initiation cite panel. It can be seen that the realistic anisotropy of the heart causes at least twice faster termination of re-entry. It also can be seen that indeed the present in the raw DT-MRI model leftover piece of tissue connected to the apical region of the heart has served as an artificial “capacitor” affecting the dynamics of the re-entry, and significantly prolongated life time of the re-entry initiated at particular locations/orientation respective to the “capacitor”.

Finally, the 3D anatomically realistic simulations of the foetal heart show that the realistic anisotropy of the heart causes the fast transient distortion of the vortex filament, and the typical fast drift towards the apex area of the inexcitable boundary of the heart, which ultimately results in the fast self-termination of the excitation vortex, see Figs. 3-10 and the corresponding movies in the Supplementary Material section LABEL:Suppl.

IV DISCUSSION AND FUTURE DIRECTIONS

Although the role of heart anatomy and anisotropy in the origin and sustainability of cardiac arrhythmias has been appreciated for a long time, the experimental evidence capable to clarify the detail of the effects of the heart anatomy on the persistent cardiac arrhythmias and fibrillation are limited.

In particular, the theoretically plausible hypothesis that the anisotropic discontinuities in the heart might be a source of rise for cardiac re-entry due to the abrupt change in conduction velocity and wavefront curvatureFenton and Karma (1998); Spach (2001); Smaill et al. (2004) was in controversy with the observation that the transmural fiber arrangement, including the range of transmural change in fiber angle in ventricular wall, although varied between species(Hunter et al., 1997, p. 173), was consistent within a species. So that the question was that, if the pro-arrhythmic mechanism of cardiac re-entry initiation by the anisotropic discontinuities in a heartFenton and Karma (1998); Spach (2001); Smaill et al. (2004) was correct, what would then have been a reason for the consistent structure (Hunter et al., 1997, p. 173) of the anisotropic discontinuities in healthy mammalian hearts. The combination of the High Performance Computing with the high-resolution DT-MRI based anatomy models of the heart allows anatomically realistic in-silico testing of the effects of individual heart anatomy and anisotropy on the cardiac re-entry dynamics Kharche et al. (2015a, b); Antonioletti et al. (2017). In this paper, for the first time, we present the anatomy and myofiber structure realistic computer simulation study of the cardiac re-entry dynamics in the DT-MRI based model of the human foetal heart Pervolaraki et al. (2013).

The comparative isotropic vs anisotropic simulation of the otherwise homogeneous foetal heart shows that, in the 2D slice of the heart, the realistic fiber anisotropy might change the re-entry dynamics from pinning to the sharp end of the septum cuneiform opening, Fig. 2(a), into a fast anatomical re-entry around the opening, Fig. 2(b). Because of the 2D re-entry pinning to either the sharp end of the septum opening in the isotropic simulation in Fig. 2(a), or to the whole septum opening as an anatomical re-entry in the anisotropic simulation in Fig. 2(b), despite of the only basic segmentation of the MRI model into the tissue/not tissue points, and the ventricles not being isolated from the atria, the tip of the re-entry never got from the ventricles into the atria, Fig. 2. Although, from the cardiac physiology point of view, the

only basic segmentation of the raw DT-MRI data Pervolaraki et al. (2013) into the tissue/non-tissue pixels might be seen as a major limitation of the study, from the non-linear science point of view, the use of the raw MRI data as an example of a nature provided medium to study a re-entry dynamics gives an important insight into the pure anatomy induced drift in an otherwise homogeneous 2D medium, and into the possibility of pinning of the re-entry not to a major blood vessel but to a sharp end of an anatomical openning Biktasheva, Dierckx, and Biktashev (2015); and into that a real fiber anisotropy is capable to turn the pinned re-entry into an anatomical one. Importantly though, the 2D simulations in Fig. 2 are an important step to highlight the role and the necessity of the whole heart structure in the re-entry dynamics and self-termination.

In the 3D DT-MRI based isotropic model of the foetal heart, depending on the initial location/orientation of the organising filament of the excitation vortex, the geometry of the foetal heart might sustain perpetual cardiac re-entry even with a positive filament tension, Fig. 3. However, if the same positive filament tension vortex is initiated at the exactly same location/orientation and in the same anatomical environment in the full anisotropic 3D DT-MRI based model of the heart, the realistic fiber structure of the foetal heart facilitates fast self-termination of cardiac re-entry, Fig. 4.

From the respective comparison of the “isotropic vs anisotropic” simulations in FIG. 3 vs FIG. 4, and FIG. 7 vs FIG. 8, it can be seen that, whereas the re-entry filaments were capable to penetrate from the ventricles to atria in the isotropic simulations shown in FIG. 3 and FIG. 7, the abrupt change in the fiber angles between the atria and the ventricles, which can be seen in FIG. 1, did not allow the re-entry filaments to get from the ventricles to atria in the anisotropc simulations shown in FIG. 4 and FIG. 8, so that the ventricles’ anisotropy could complete the speedy elimination of the re-entry within its single rotation.

The comparison of the re-entry termination times in the raw DT-MRI data model whole heart simulations shown in Fig. 3, Fig. 4, Fig. 5, and Fig. 6, with the corresponding series of the “edited” MRI model whole heart simulations shown in Fig. 7, Fig. 8, Fig. 9, and Fig. 10, showed that, although the filament of the re-entry never got through into the small piece of excitable tissue accidentally adjacent to the apical region of the heart, the adjacent tissue served as a “capacitor” significantly prolongating the life time of the re-entry initiated at a particular location/orientation respective to the “capacitor’s” own location/orientation. See for the quantitative comparison of the re-entry termination times Fig. 11 and Fig. 12, where the bigger number and the total length of the filaments tend to correlate with the faster termination of re-entry, though these fail to identify the persistent re-entry in FIG. 3 simulation.

The “isotropic vs anisotropic” comparison of re-entry self-termination time in both the original raw DT-MRI simulations series, and in the “edited” MRI whole heart simulations, confirmed that the real anisotropy of the heart speeds up cardiac re-entry self-termination. The re-entry self-termination times provided for the summarised comparison in Fig. 11 and Fig. 12, show that, regardless of with or without the “leftover” piece adjacent to the apex, the anisotropy of the heart speeds up cardiac re-entry self-termination. Fig. 12 shows that anisotropy increases the transient number and the transient total length of the filaments. The bigger transient number and the total length of the filaments tend to correlate with the faster termination of re-entry. The biggest transient total length of the filaments was in case of the re-entry initiated along the $y$ axis, see panel (d) in Fig. 12, which ensured the re-entry fastest termination. It can be seen from FIG. 5, FIG. 6, FIG. 9, and FIG. 10, that the initial position of the filament along the $y$ axis allowed the filament to grow intramurally, thus maximally increasing the transient total length of the filaments, and speeding up their termination.

The simulations with the “edited” MRI image of thus completely isolated heart, Fig. 7, Fig. 8, Fig. 9, and Fig. 10, in comparison with the original DT-MRI model simulations, Fig. 3, Fig. 4, Fig. 5, Fig. 6, provide an important new biological insight into the problem of cardiac re-entry dynamics. Namely, that an excitable tissue accidentally adjacent to the heart might serve as a capacitor capable to prolongate time of cardiac re-entry self-termination, see for the respective comparison the simulations in Fig. 3 against Fig. 7, Fig. 4 against Fig. 8, Fig. 5 against Fig. 9, and Fig. 6 against Fig. 10, all also summarised in Fig. 11. The latter suggests a possible new mechanism for persistent cardiac re-entry. So that if, apart from the major blood vessels normally adjacent to the heart in vivo and affecting re-entry dynamics, there were also an accidental “touching” of the heart by an adjacent excitable tissue, for example, due to the change of posture in the night sleep, the “incidental capacitor” effect could prolongate the time of cardiac re-entry self-termination, or indeed failure to self-terminate, which could be an explanation to the elusive and difficult to reporduce but statistically salient data for longer episodes of arrhythmias reported in the night ECGs as opposed to the on average shorter arrhythmias in the day time ECGs. Although our simulations using the original raw DT-MRI data with the small piece of the foreign leftover tissue, could have been seen a limitation of the study, the real heart in vivo does not exist in complete isolation from the main blood vessels and other neighboring tissues. So, we believe that our “incidental” leftover tissue results only once more confirm the importance and the necessity of taking into account the real anatomical settings and surrounding of the heart for the full appreciation of cardiac re-entry dynamics.

The BeatBox DT-MRI based  in-silico model comparative study confirms the cardiac anatomy and anisotropy functional effect on cardiac re-entry sustainability as opposed to its self termination, the pinning of the re-entry to anatomical features, its transformation from pinned to anatomical re-entry, and the re-entry self-termination caused by the anisotropy of the tissue.

One of the limitations of the present study is the use of the simplified Fitz-Hugh-Nagumo Winfree (1991) excitation model Eq. (2). The simplified FHN model with the excitation kinetics parameters $\alpha=0.3$, $\beta=0.71$, $\gamma=0.5$, which, in an infinite homogeneous isotropic excitable medium, supports a rigidly rotating vortex with positive filament tension Biktashev, Holden, and Zhang (1994), was chosen for this study in order to eliminate the effects of realistic cell excitation kinetics, such as e.g. meander Winfree (1991), alternansKarma (1994), negative filament tensionBiktashev, Holden, and Zhang (1994), etc., in order to enhance and highlight the pure effects of the heart anatomy and anisotropy on the cardiac re-entry outcome. The realistic cell excitation models should be used in the future studies, in order to clarify the particular effects and interplay of the cell excitation kinetics with the heart anatomy and anisotropy.

As it can be seen from Fig. 1 (for the color-encoded fractional anisotropy (FA) and for the color-encoded all the three components of the fiber angles see Figure 4 in Pervolaraki et al (Pervolaraki et al., 2013, p. 5)), formation of the fiber structure at the epicardium and endocardium is not completed yet in the foetal heart, so that only the already formed intramural laminar structure of the fibers can affect the dynamics of cardiac re-entry. Although the use of the not fully formed foetal heart can be seen as a limitation of the study, on the other hand, it may be said that the chaotic epicardium and endocardium fiber orientation prevents the foetal heart re-entry from pinning to the fine anatomical features which were yet to be developed at the fully formed Pervolaraki et al. (2013) endocardium later on. The possible differences in the anatomy and fiber structure between the foetal heart used here and fully formed/adult hearts in general, could have seriously affected the simulations, such as in the case of e.g. reported pinning of cardiac re-entry to the junction of pectinate muscles with crystae terminalis in adult human atrium Wu et al. (1998); Yamazaki et al. (2012); Kharche et al. (2015a). That is, although it is possible to initiate a cardiac re-entry in the tiny $1.4g$ (at $143$ DGA) foetal heart Pervolaraki et al. (2013), the already formed intramural laminar fiber anisotropy of the foetal heart facilitates the re-entry self-termination, Fig. 11. With the hindsight of the present study, in a fully formed adult heart, because of the presence of the pinning opportunities provided by the endocardium anatomical features Wu et al. (1998); Yamazaki et al. (2012); Kharche et al. (2015a), there must exist additional mechanisms to facilitate cardiac re-entry self-termination Clayton et al. (1993).

The most serious limitation of the study is that only the basic segmentation of the raw DT-MRI data Pervolaraki et al. (2013) into the tissue/non-tissue pixels based on the MRI luminosity threshold, and only the primary eigenvalues of fibres orientation, were taken into account in the BeatBox Antonioletti et al. (2017) computer simulation of the cardiac re-entry dynamics. Further levels of the model segmentation, in order to take into account e.g. the heart collagen skeleton, isolation of ventricles from atria, etc., will inevitably change the outcome of the re-entry, by adding the electrically impermeable barriers to cardiac re-entry. Currently, this further segmentation is added into DT-MRI based models via complex rule based image post-processing Lombaert et al. (2012); Gahm, Kung, and Ennis (2013), which not only limits the available segmented DT-MRI cardiac anatomy data sets, but also inevitably brings in an artificial assumption/limitation element into these models. From the non-linear science point of view, which we have persued in this initial study, the rationale was to use the raw DT-MRI data as an example of a nature provided medium to study a re-entry dynamics. In the future, the multichannel computer tomography might offer an automatic tissue segmentation, so that the multi-level segmented DT-MRI based heart anatomy models might become more available, and be used in the BeatBox Antonioletti et al. (2017) anatomically and biophysically realistic simulations of cardiac re-entry dynamics.

Finally, we believe that a simple “mechanistic” explanation, although often craved for, might be rather inadequate/premature here, and will require better theoretical understanding of the demonstrated potential effect of the heart anisotropy on cardiac re-entry dynamics, for it is not a particular feature, or a sequence of anisotropy features, but rather the whole complex of the shape, anisotropy, and the exact heart position within the body surrounding, which affects the re-entry dynamics in a particular way, and which seems to have had evolved in order to ensure the fastest self-termination of cardiac re-entry. If our hypothesis is correct, it might explain the difficulties with reproducibility of the arrhythmia in vivo and in an isolated heart. The most important novel finding of the paper is that, contrary to what currently seems to be a commonly accepted view of the pro-arrhythmic nature of cardiac anisotropy, the point of view based on the mainly theoretical and simplified anatomy models studies, for the first time ever, and for the first time in a real whole heart DT-MRI based model, we have demonstrated that the real life heart anisotropy might have rather an anti-arrhythmic effect, as it facilitates the fastest self-termination of cardiac re-entry.

Acknowledgements.

We acknowledge the support of the UK Medical Research Council grant G1100357 for the human foetal heart DT-MRI data sets. We also wish to acknowledge the support of the BeatBox software development project by EPSRC (UK) grants EP/I029664 and EP/P008690/1. We thank all the developers of the BeatBox HPC Simulation Environment for Biophysically and Anatomically Realistic Cardiac Electrophysiology. We are grateful to Professor V.N.Biktashev for much appreciated advice and discussion.