- Email: [email protected]
Spiral wave breakup in excitable media with an inhomogeneity of conduction anisotropy
Computers in Biology and Medicine 40 (2010) 775–780
Contents lists available at ScienceDirect
Computers in Biology and Medicine journal homepage: www.elsevier.com/locate/cbm
Spiral wave breakup in excitable media with an inhomogeneity of conduction anisotropy P. Kuklik a,c,n, L. Szumowski b, P. Sanders c, J.J. Z˙ebrowski a a
Physics of Complex Systems, Faculty of Physics, Warsaw University of Technology, ul. Koszykowa 75, Warszawa, Poland Institute of Cardiology, ul. Alpejska 42, Warsaw, Poland c Cardiovascular Research Center, Department of Cardiology, Royal Adelaide Hospital, and the Disciplines of Medicine and Physiology, University of Adelaide, Adelaide, Australia b
a r t i c l e in fo
abstract
Article history: Received 31 March 2010 Accepted 22 July 2010
Many conditions remodel the heart muscle such that it results in a perturbation of cells coupling. The effect of this perturbation on the stability of the spiral waves of electrochemical activity is not clear. We used the FitzHugh–Nagumo model of an excitable medium to model the conduction of the activation waves in a two-dimensional system with inhomogeneous anisotropy level. Inhomogeneity of the anisotropy level was modeled by adding Gaussian noise to diffusion coefficients corresponding with lateral coupling of the cells. Low noise levels resulted in a stable propagation of the spiral wave. For large noise level conduction was not possible due to insufficient coupling in direction perpendicular to fibers. For intermediate noise intensities, the initial wave broke up into several independent spiral waves or waves circulating around conduction obstacles. At an optimal noise intensity, the number of wavelets was maximized—a form of anti-coherent resonance was obtained. Our results suggest that the inhomogeneity of conduction anisotropy may promote wave breakup and hence play an important role in the initiation and perpetuation of the cardiac arrhythmias. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Inhomogeneous excitable media Spiral wave dynamics Spiral wave breakup Atrial fibrillation Heart arrhythmias
1. Introduction The contraction of the heart muscle is driven by a propagating electrochemical wave of activation. In normal conditions, the activation wave spreads from the sino-atrial node through the heart in a precise, complex manner and ends propagation when there is no more excitable muscle ahead of it [1]. One of pathological forms of wave activity is a spiral wave—a spiral shaped activation wave reentering the same fragment of tissue it has already propagated through [2,3]. Such a perturbation interferes with normal heart rhythm since the spiral wave circulates at a very high rate disturbing the mechanical function of the heart. There are several mechanisms of spiral wave initiation and breakup leading to creation of multiple wavelets [2–11]. In this work, we show a wave breakup caused by the inhomogeneity of the conduction anisotropy. Atrial fibrillation is the most common heart rhythm disorder to affect humans [12]. It has been postulated that one of the mechanisms that leads to the initiation and maintenance of this arrhythmia is the anchoring of spiral waves to myocardial tissue.
n Corresponding author at: Cardiovascular Research Center, Department of Cardiology, L5 McEwin Building, Royal Adelaide Hospital, Adelaide, SA 5000, Australia. Tel.: + 61 8 822 22723; fax: + 61 8 8313 1166. E-mail address: [email protected] (P. Kuklik).
0010-4825/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.compbiomed.2010.07.005
To date, only high density electrical and optical mapping techniques were able to provide some insights into dynamics of spiral waves during atrial fibrillation [13]. However, emerging picture is complex and full description of spiral waves dynamics during atrial fibrillation remains still unknown. A deeper understanding of the spiral wave dynamics may facilitate the design of newer treatment strategies. In this work, we concentrated on spiral wave initiation due to a wave breakup. Wave breakup may occur due to a structural inhomogeneity (e.g. inhomogeneity of the diffusion coefficient [4]) or due the dynamical properties of the system [5]. Several mechanisms of wave breakup due to structural inhomogeneity were recently investigated. Among them: random removal of small fragments of the system [6], perturbation of the position of cells in a cellular automaton model [7], local decrease of the conduction velocity [8], periodic modulation of the conduction velocity [9], inhomogeneity of the density of the system elements [10], periodic modulation of the kinetic parameters [11], inhomogeneity of the conduction velocity [4,14], mechanical deformation of the system [15]. Initiation of the spiral wave due to the inhomogeneity of the tissue fibers orientation was investigated, but only in the three dimensional case [16,17]. In this study, we investigate a role of the inhomogeneous anisotropic uncoupling which would be a model of lateral uncoupling of the cardiac cells due to pathological muscle remodeling.
776
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Conduction of the electrochemical wave of activity in atrial myocardium is determined by several factors including the tissue architecture and the extent and distribution of cell-to-cell connections by gap junction proteins (connexins). Strong coupling between cells enables the rapid transport of the ions between cells and therefore an increase of the wave velocity. Similarly, poor coupling may result in a slowing down of the wave and even in propagation failure. The nature of this electrical transit of conduction makes it sensitive to the distance between the cells and to the number and distribution of connexin proteins. Furthermore, many of the conditions that predispose to the development of arrhythmia are known to remodel the heart muscle leading to a perturbation of wave propagation [18,19]. There are several studies focusing on the effect of the cell coupling on wave conduction [20–24]. Interesting results were obtained by the group of Jacqes de Bakker in an experiment focusing on the mechanisms of the slow conduction in the heart muscle after infarction [20]. A detailed, microscopic study of the patterns of activation revealed a pattern named by the authors a ‘‘zig-zag pattern of activation’’ and characterized by a dissociation of the wave front into independent wavelets propagating through separate muscle strands. This dissociation was possible due to a poor coupling between muscle cells and was caused by an excessive expression of collagen in the infarcted muscle. Such a pattern of propagation was suggested as a possible explanation of the arrhythmogeneity of the post-infarct heart muscle. In this study, we used the FitzHugh–Nagumo model of an excitable medium to investigate the effect of the random perturbation of cell coupling on the stability of a spiral wave. We obtained a form of spatial anti-coherent resonance. Note that in the usual coherent resonance there is an optimal intensity of the external noise at which noise-induced temporal oscillations in the system are most coherent [25]. In our case, the optimal level of the spatial noise added to the cell coupling resulted in an increase of the spatiotemporal complexity of the system. Therefore, we used a term ‘‘anti-coherent resonance’’ to describe the result found in this study.
constant. For lower values of b (b o0.64) the system became oscillatory due to a very low inhibition. For b 40.78 and m 40.36 the size of the meandering pattern of the spiral wave exceeded the size of the sample leading to the termination of the wave. Integration of the FHN equations was carried out in a two dimensional geometry based on a 50 50 rectangular grid using zero-flux boundary conditions. Calculation of the number of the wavelets was conducted on a part of the sample, defined by circle subtended in the center of the sample with radius equal to half of the sample size. The wheel geometry was used to minimize the effect of the sample corners on counting the number of separate wavelets (corners provide a possibility for a wave to be separated into distinct parts and introduce artifacts to waves counting). The Euler integration scheme was used with a time step Dt ¼0.07 and the space step Dx¼ Dy ¼0.4. Since we do not aim to quantitatively reconstruct physiological phenomena, we used dimensionless units in our study. We introduced two diffusion coupling coefficients: for the i,j horizontal direction (Di,j x ) and for the vertical direction (Dy ). As i,j the initial condition, we set Di,j ¼ 0:5 and D ¼ 0:125. Those values x y were set in order to model the anisotropic conduction in the heart tissue. Due to a lack of consensus in the literature (published coupling coefficients ratios range from 1:1 to 9:1 [28,29]), we set the ratio of the diffusion coefficients in the horizontal and perpendicular direction to 1:4. An inhomogeneity of the cells coupling was introduced by addition of spatial, Gaussian noise to Di,j y which would correspond with lateral decoupling of the cells: i,j Di,j y -maxð0,Dy ZÞ
ð3Þ
where Z is given by a Gaussian distribution with a fixed standard deviation s. An example of the distribution of the diffusion coefficients is shown in Fig. 1. With increase in s the amplitude of the coupling in vertical direction (Di,j y ) is decreasing, reaching at
2. Methods We used simplified FitzHugh–Nagumo (FHN) model of an excitable medium [26]. The system was modeled as a twodimensional lattice of the diffusively coupled FHN elements. For an in depth discussion of model properties see [27].The equations of the model were the following: dvi,j vi þ 1,j þvi1,j 2vi,j 1 1 ¼ vi,j v3i,j ui,j þ Di,j x 2 3 dt m Dx i,j vi,j þ 1 þvi,j þ 1 2vi,j ð1Þ þ Dy Dy2
dui,j ¼ mðvi,j þ bgui,j Þ dt
ð2Þ
where vi,j is the transmembrane voltage of the element {i,j}, ui,j is the total slow current, g ¼ 0.5, b ¼0.7, m ¼ 0.3. In order to investigate phenomena in different dynamical states of the FHN model, we varied the parameters b and m in a range in which the propagating spiral wave are sustained. The parameter b responsible for the magnitude of the inhibition current was kept between 0.65 and 0.75 so as to keep m ¼0.3 constant. Alternatively, the parameter m responsible for relative time scales of v and u variables was set between 0.22 and 0.32 to keep b ¼ 0.7
Fig. 1. Example of the spatial distribution of the diffusion coefficients after the addition of Gaussian noise with s ¼0.1. The length of a horizontal line corresponds to the value of Dxi,j and the length of a vertical line with the value of Dyi,j . Introduction of the noise at each grid point perturbs the initial structure of the system.
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
some point a level at which there is a conduction block—a failure of propagation in vertical direction. For in depth discussion of propagation failure in anisotropic media see [30]. Spiral waves were initiated in the system by setting the spatial distribution of v and u to the distribution corresponding with a propagating spiral wave: for each element of the simulation grid, its phase was calculated as
y¼
rkj2 pkn 2 pk
ð4Þ
where r and j are polar coordinates of a given point with respect to the center of the spiral wave, k is its winding number and n is an integer set so as to obtain a phase between 0 and 1. Values of v and u for a given element of the grid were set depending on the phase. y ¼0 corresponded to the resting state and y ¼1 corresponded to the wave front depolarization. All intermediate values of the variables were copied from the stationary distribution of the variables within the action potential cycle. Such prepared distribution of v and u variables (see Fig. 5, t¼5 for example of such distribution) was an initial condition from which spiral wave easily evolved. The number of wavelets was calculated by assessing the number of distinct depolarized areas (i.e. such with v40).
3. Results A low magnitude of the spatial noise (s o0.1) added to the diffusion coefficients did not significantly affect the stability of the spiral wave. Only in few cases the wave broke up (see average number of wavelets for s o0.1 in Figs. 3 and 4). Large magnitude of the spatial noise added to the diffusion coefficients resulted in most cases in the termination of the spiral wave: the wave was unable to rotate due to the high anisotropy level (see average number of wavelets for s 40.5 in Figs. 3 and 4). For intermediate noise intensity, the initial wave broke up into several independent spiral waves or waves circulating around conduction obstacles. Due to the stochastic nature of the inhomogeneity generation, we conducted a series of iterations for a given value of the standard deviation of the noise, s. The number of distinct wavelets as a function of the time is shown in Fig. 2a. The average number of the wavelets calculated in a sliding window
777
20 t.u. (time units) in size saturates after 500 t.u. (we obtained a similar result for all s). Consequently, we fixed end time for the simulation at Tend ¼ 1000 t.u. We performed calculations for s between 0 and 0.7 in increments of 0.02. Next, for a given value of s, we repeated the simulation calculating the average number of wavelets after time Tend. An example of the result for s ¼0.16 is shown in Fig. 2b. For further calculations, we assumed 100 repetitions as sufficient to obtain average wavelet number for given s. For each s, we calculated the average number of distinct wavelets and their standard deviation. The example result is shown in Fig. 3. In this example, for s ¼0.24 the number of wavelets (8.1171.75) attained a maximum—a form of spatial anti-coherent resonance was obtained. Distributions of the number of wavelets for different values of parameters b and m had similar unimodal characteristics (Figs. 3b and 4a, b). Increase of the current inhibition (parameter b) caused a decrease of the mean number of wavelets and a shift of the peak position toward smaller s (Fig. 3b). An increase of the parameter m corresponded with a shift of the peak position from s ¼0.32 (m ¼0.2) to s ¼0.1 (m ¼ 0.32). The maximum number of wavelets was greatest for intermediate values of m (8.23 for m ¼ 0.26 and 8.3 for m ¼0.28) (Fig. 4) and lower for smaller and large values of m. An example of the evolution of the system for s ¼0.2 is shown in Fig. 5. The spiral wave (t ¼5) is unable to rotate initially due to the high anisotropy level resulting in the sliding of the tip to the left and a block of lower part of the arm (t ¼10). After some time, the sliding tip was able to pass the excitation to a row of the elements above, which resulted in the creation of a wavelet propagating backwards (upper arrow, t ¼20). A similar event happened for a separated end of the lower part of the wave (lower arrow, t ¼20). The two additional wavelets evolved further into a more complex organization shown on the remaining images. Wavelets observed in the system fall into two categories: wavelets propagating in loops or small spiral waves rotating in areas able to support wave rotation. An example of a loop is shown in Fig. 6. A loop consists of at least two sets of disconnected rows of elements (with a small value of Dy preventing propagation) ended with a kind of passage consisting of elements with diffusion coefficients allowing propagation in the horizontal and vertical directions. The wavelet in such a loop propagates along one of disconnected rows jumping to other one
Fig. 2. Number of distinct wavelets in the system as a function of the time for s ¼ 0.16, b ¼0.7 and m ¼ 0.3 (a). The average number saturates after t¼ 500 t.u. (also for other values of s). Average number of the distinct wavelets after t¼ 500 t.u. as a function of the number of simulation repetitions for s ¼ 0.16 (b).
778
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Fig. 3. Mean number of wavelets in the system as a function of the standard deviation of noise for b ¼ 0.675 and m ¼ 0.3. For each value of s, the simulations were repeated 100 times and the result averaged. The vertical bars denote the standard deviation of the number of wavelets (a). Mean number of wavelets in the system as a function of the standard deviation of noise for different values of parameter b (with constant m ¼ 0.3) (b).
Fig. 4. Mean number of wavelets in the system as a function of the standard deviation of noise for different values of parameter m (with constant b ¼0.7). Plots for m ¼ 0.2, 0.22, 0.24 and 0.26 (a). Plots for m ¼ 0.28, 0.3 and 0.32 (b).
Fig. 5. Example of the spiral wave breakup for s ¼ 0.15, b ¼ 0.7, m ¼ 0.3. A part of the initial spiral wave arm propagating perpendicularly to the direction of maximal diffusion (vertical) is blocked (t ¼20) leading ultimately to the formation of new spiral or circulating waves.
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Fig. 6. Example of a loop created by addition of noise to the distribution of the diffusion coefficients (explanation in the text). Grey arrows denote path of the circulating wave of activity.
through the passages at the ends and propagating backward to reenter the first row and repeat the cycle.
4. Discussion The results presented in this study show the relationship between the inhomogeneity of the diffusion anisotropy and the number of the distinct wavelets created after introduction of the spiral wave into system. Our main finding is that the local perturbation of the diffusion anisotropy level promotes spiral wave breakup. We found that there is an optimal perturbation level, which maximizes the number of distinct wavelets. For low and high perturbation levels, we obtained, respectively, a stable propagation and the termination and of the wave. Such behavior may be classified as an anti-coherent resonance due to the inhomogeneity introduced as a spatial noise of the anisotropy level. Termination of the wave in case of high level of spatial noise intensity is to be expected, due to conduction block in vertical direction for too low values of Dy. Such coupling results in a system equivalent to a set of 1D uncoupled cables in which rotation of spiral wave is not possible. A decrease of the mean number of wavelets with b increase (Fig. 3b) may be explained by a less likely propagation of small wavelets along the loops caused by large inhibition current. Increasing inhibition results with total propagation block for b ¼0.75 (Fig. 3b). For lower values of b, the number of the wavelets increases due to lower excitation threshold for lower b-loops unable to sustain wavelets for base parameters (b ¼0.7) due to propagation failure caused by low coupling, are now able to conduct due to increased excitability of the system. Thus, the increase of a number of available loops for propagation results in an increase of wavelets number. A decrease of the mean number of wavelets for smaller m (m o0.24, Fig. 4a) and for larger m (m 40.28, Fig. 4b) was related with different mechanisms. Decrease for large m was related with propagation failure of small wavelets. The onset of the inhibiting current is more rapid in this case causing, in effect, an inhibition of the voltage. In effect, small wavelets are more prone to propagation failure due to an insufficient depolarization area in front of the wave thus decreasing the mean number of wavelets. However, for small values of m the decrease in a number of wavelets was related with a less likely breakup of the wave. Late offset of the inhibition current increased the ability of the wave to propagate through areas of the low coupling thus resulting in more rare incident of the wave breakup and formation of the wavelets. A cause of the spiral wave breakup by conduction block may be either of dynamical origin (e.g. inability to excite due to short interval between consecutive activations) [5] or may be related with a structure of the system (e.g. poor cells coupling) [6–11,14,16,17]. Wave breakup analyzed in this study falls into structural category, related with a modification of the cells
779
coupling—a breakup of the activation wave was caused by propagation failure related to low diffusive coupling. There are several factors resulting in a propagation failure [30] and diffusive coupling used in our system may be regarded as an effective coupling factor affecting wave conduction. Low coupling in vertical direction resulted in a conduction block that either stopped a whole wave (if wavelet size was similar to area of low coupling) or caused annihilation of the part of the wavelet. The mechanism of the wave breakup presented here may be one of the mechanisms of an arrhythmia in diseased heart muscle. Remodeling of the diseased heart muscle leads to an inhomogeneous perturbation of cells coupling [19]. Waves of activity entering such a region may break up and form complex patterns of activation containing spiral waves and wavelets circulating around loops caused by a local separation of the muscle fibers, as shown in de Bakker study reporting ‘‘zig-zag’’ type of propagation [20]. Such type of propagation would correspond with a wavelet in a loop propagating along one of disconnected rows and ‘‘jumping’’ to other one through the connections between rows (example of such loop shown in Fig. 6). Remodeling of the tissue structure and function resulting in ‘‘zig-zag’’ type of activation was caused by infarction in de Bakker paper, but it is possible that other disease related factors may also result in similar type of electrical remodeling [1]. In our previous work, we showed the formation of the spiral waves and waves circulating around non-conducting areas in an isotropic system with an inhomogeneous depression of the diffusion coefficient [4]. The areas of the depressed diffusion formed non-conducting obstacles enabling waves to circulate around such loops. In the system presented in here, we obtained similar propagation properties; however, the loops were formed by a specific inhomogeneity of the lateral coupling while keeping the coupling along the fiber constant. In a heart tissue, such scenario would correspond to an inhomogeneous fibrosis (decrease of the lateral cell coupling with a preserved ability to propagate along the fiber orientation) whereas our previous work would correspond to total uncoupling from the neighboring cells. The frequency of the waves propagating in the loops is a function of the loop length. Shorter loops result in a higher circulation frequency. In a case of the coexistence of several loops (typical for our system), waves emanating from the shortest loop would interfere with the waves circulating around the slower loops. This mechanism was described in [31] where the authors showed for inhomogeneous media that the spiral wave with the highest frequency sweeps away slower spiral waves. However, in our case, the slower activity is not due to the kinetic properties of the model, but rather due to the structural inhomogeneity. In our system, the fastest source of activity will not eliminate slower ones but rather override them with a more frequent train of waves. Thus, after termination of the fastest activity, the wavelets may resume circulation along loops. This property may be relevant for the understanding of the mechanism of atrial fibrillation and radiofrequency ablation procedure. During the ablation procedure, atrial tissue is progressively destroyed. Clinical studies showed, that the frequency of the atrial electrical activity in course of the ablation tends to decrease (ultimately resulting in the restoration of the sinus rhythm) [32]. The mechanism of this phenomenon is still unknown. Our study suggests, that this phenomenon may be explained as progressive removal of the loops formed by fibrotic areas sustaining circulating wavelets. However, the subject requires further study. Decrease of the wavelets number for a large level of noise may lead to a conclusion that noise may have certain anti-arrhythmic role, preventing progression of the arrhythmia into more dangerous and persistent form. Such conclusion has to be formulated very carefully due to overwhelming evidence of a pro-arrhythmic role of
780
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
the heart tissue heterogeneity and fibrosis. Similar conclusions (a warning against oversimplification of role of the heterogeneity in suppressing arrhythmias) were formulated in work of Panfilov et al. [33,34] in which authors studied break up of spiral wave in AlievPanfilov model of excitable media and shown, that addition of non-excitable cells prevents spiral wave from breakup. In our model, the increase of the noise intensity promoted the spiral wave breakup but at some point it caused a propagation termination due to too low lateral coupling to sustain propagation. Since this conclusion contradicts to clinical experience (e.g. [20]), this result should be carefully interpreted. We used the following simplifications in our model:
we used FitzHugh–Nagumo model of the electrical activity of
the heart tissue which is a simplified with respect to the complex behavior of ionic channels; simplified geometry of the system – a two dimensional flat rectangle – which neglects three dimensional nature of the atria, architecture of the muscle fibers, pectinate structure and heterogeneity of the kinetic characteristics of the cells.
The essential features of the mechanism presented in this study comprise:
propagation failure due to low coupling in the direction of propagation;
ability of the system to sustain conduction along rows of laterally uncoupled cells; and
ability of the system to sustain spiral wave. In order to assess whether our simplifications have an effect on the obtained results, all above mechanism features should be independent on the model simplifications. Propagation failure due to low coupling is a generic feature of every excitable system irrespectively on system geometry and dimension (for every configuration of the activation there is a critical value of the coupling that will induce a propagation failure). Ability of the system to sustain conduction along rows of cells (equivalent to 1D cable) is another generic feature of the excitable media reasonably expected to be present in all models of the heart electric activity and irrespectively of model geometry (which is reduced to a fragment of 1D cable by lateral uncoupling). Ability of the system to sustain spiral waves is more complex, since it was shown that both kinetics of the model and system geometry play a role in sustainability of the spiral wave [35]. Therefore, our result is valid only for systems able to sustain spiral wave.
Conflict of interest statement None declared.
Acknowledgement Dr. Sanders is supported by the National Heart Foundation of Australia. References [1] D.P. Zipes, J. Jalife, Cardiac Electrophysiology from Cell to Bedside, W.B. Saunders Company, 2004. [2] E.M. Cherry, F.H. Fenton, Visualization of spiral and scroll waves in simulated and experimental cardiac tissue, New J. Phys. 10 (2008) 125016.
[3] J.M. Davidenko, A.M. Pertsov, R. Salomonsz, W.T. Baxter, J. Jalife, Stationary and drifting spiral waves of excitation in isolated cardiac muscle, Nature 355 (1992) 349–351. [4] P. Kuklik, J.J. Z˙ebrowski, Reentry wave formation in excitable media with stochastically generated inhomogeneities, Chaos 15 (2005) 033301. [5] F.H. Fenton, E.M. Cherry, H.M. Hastings, S.J. Evans, Multiple mechanisms of spiral wave breakup in a model of cardiac electrical activity, Chaos 12 (2002) 852–892. [6] G. Bub, A. Shrier, Propagation through heterogeneous substrates in simple excitable media models, Chaos 12 (2002) 747–753. [7] M. Markus, B. Hess, Isotropic cellular automaton for modelling excitable media, Nature 347 (1990) 56–58. [8] T.J. Lewis, J.P. Keener, Wave-block in excitable media due to regions of depressed excitability, SIAM J. Appl. Math. 61 (1) (2000) 293–316. [9] I. Schebesch, H. Engel, Wave propagation in heterogeneous excitable media, Phys. Rev. E. 57 (4) (1998) 3905–3910. [10] G. Bub, A. Shrier, L. Glass, Spiral wave generation in heterogeneous excitable media, Phys. Rev. Lett. 88 (2002) 5. [11] J. Kneller, R. Zou, E.J. Vigmond, Z. Wang, L.J. Leon, S. Nattel, Cholinergic atrial fibrillation in a computer model of a two-dimensional sheet of canine atrial cells with realistic ionic properties, Circ. Res. 90 (2002) e73–e87. [12] F.D. Murgatroyd, A.J. Camm, Atrial arrhythmiasm, Lancet 341 (1993) 1317–1322. [13] D.S. Rosenbaum, J. Jalife, in: Optical Mapping of Cardiac Excitation and Arrhythmias, Futura Publishing Company, New York, 2001 Armonk. [14] B.E. Steinberg, L. Glass, A. Shrier, et al., The role of heterogeneities and intercellular coupling in wave propagation in cardiac tissue, Philos. Trans. Roy. Soc. A 364 (2006) 1299–1311. [15] H. Zhang, X.S. Ruan, B. Hu, Q. Ouyang, Spiral breakup due to mechanical deformation in excitable media, Phys. Rev. E 70 (2004) 016212. [16] A.V. Panfilov, J.P. Keener, Reentry in three-dimensional Fitzhugh–Nagumo medium with rotational anisotropy, Physica D 84 (1995) 545–552. [17] F. Fenton, A. Karma, Fiber-rotation-induced vortex turbulence in thick myocardium, Phys. Rev. Lett. 81 (1998) 2. [18] P. Kohl, Heterogeneous cell coupling in the heart: an electrophysiological role for fibroblasts, Circ. Res. 93 (2003) 381. [19] M.S. Spach, J.F. Heidlage, P.C. Dolber, R.C. Barr, Electrophysiological effects of remodeling cardiac gap junctions and cell size: experimental and model studies of normal cardiac growth, Circ. Res. 86 (2000) 302. [20] J.M.T. De Bakker, F.J.L. Van Capelle, M.J. Janse, et al., Slow conduction in the infarcted human heart :‘‘zigzag’’ course of activation, Circulation 88 (1993) 915–926. [21] J.E. Saffitz, H.L. Kanter, K.G. Green, T.K. Tolley, E.C. Beyer, Tissue-specific determinants of anisotropic conduction velocity in canine atrial and ventricular myocardium, Circ. Res. 74 (1994) 1065–1070. [22] N.S. Peters, J. Coromilas, N.J. Severs, A.L. Wit, Disturbed connexin43 gap junction distribution correlates with the location of reentrant circuits in the epicardial border zone of healing canine infarcts that cause ventricular tachycardia, Circulation 95 (1997) 988–996. [23] M.S. Spach, Anisotropic structural complexities in the genesis of reentrant arrhythmias, Circulation 84 (1991) 1447–1450. [24] C.Y Chung, H. Bien, E. Entcheva, The role of cardiac tissue alignment in modulating electrical function, J. Cardiovasc. Electrophysiol. 18 (12) (2007) 1323–1329. [25] B. Lindner, J. Garcia-Ojalvo, A. Neimand, L. Schimansky-Geier, Effects of noise in excitable systems, Phys. Rep. 392 (2004) 321–424. [26] R.A. FitzHugh, Impulses and physiological states in theoretical models of nerve membrane, Biophys. J. 1 (1961) 445–466. [27] A.T. Winfree, Varieties of spiral wave behavior: an experimentalist’s approach to the theory of excitable media, Chaos 1 (3) (1991) 303–334. [28] C. Tobon, C. Ruiz, et al., Reentrant mechanisms triggered by ectopic activity in a three-dimensional realistic model of human atrium. a computer simulation study, Comput. Cardiol. 1 (2) (2008) 629–632 2008. [29] G. Seemann, C. Hoper, et al., Heterogeneous three-dimensional anatomical and electrophysiological model of human atria, Philos. Trans. Roy. Soc. A—Math. Phys. Eng. Sci. 364 (1843) (2006) 1465–1481. [30] J.P. Keener, The effects of discrete gap junction coupling on propagation in myocardium, J. Theor. Biol. 148 (1) (1991) 49–82. [31] F. Xie, Z. Qu, J.N. Weiss, A. Garfinkel, Interactions between stable spiral waves with different frequencies in cardiac tissue, Phys. Rev. E. 59 (1999) 2. [32] M. Haissaguerre, P. Sanders, M. Hocini, L.F. Hsu, D.C. Shah, C. Scavee, Y. Takahashi, M. Rotter, J.L. Pasquie, S. Garrigue, J. Clementy, P. Jais, Changes in atrial fibrillation cycle length and inducibility during catheter ablation and their relation to outcome, Circulation 109 (2004) 3007–3013. [33] A.V. Panfilov, Spiral breakup in an array of coupled cells: the role of the intercellular conductance, Phys. Rev. Lett. 88 (2002) 11. [34] K.H.W.J. ten Tusscher, A.V. Panfilov, Influence of nonexcitable cells on spiral breakup in two-dimensional and three-dimensional excitable media, Phys Rev E 68 (2003) 062902. [35] Z. Qu, Critical mass hypothesis revisited: role of dynamical wave stability in spontaneous termination of cardiac fibrillation, Am. J. Physiol. Heart Circ. Physiol. 290 (2006) H255–H263.
Contents lists available at ScienceDirect
Computers in Biology and Medicine journal homepage: www.elsevier.com/locate/cbm
Spiral wave breakup in excitable media with an inhomogeneity of conduction anisotropy P. Kuklik a,c,n, L. Szumowski b, P. Sanders c, J.J. Z˙ebrowski a a
Physics of Complex Systems, Faculty of Physics, Warsaw University of Technology, ul. Koszykowa 75, Warszawa, Poland Institute of Cardiology, ul. Alpejska 42, Warsaw, Poland c Cardiovascular Research Center, Department of Cardiology, Royal Adelaide Hospital, and the Disciplines of Medicine and Physiology, University of Adelaide, Adelaide, Australia b
a r t i c l e in fo
abstract
Article history: Received 31 March 2010 Accepted 22 July 2010
Many conditions remodel the heart muscle such that it results in a perturbation of cells coupling. The effect of this perturbation on the stability of the spiral waves of electrochemical activity is not clear. We used the FitzHugh–Nagumo model of an excitable medium to model the conduction of the activation waves in a two-dimensional system with inhomogeneous anisotropy level. Inhomogeneity of the anisotropy level was modeled by adding Gaussian noise to diffusion coefficients corresponding with lateral coupling of the cells. Low noise levels resulted in a stable propagation of the spiral wave. For large noise level conduction was not possible due to insufficient coupling in direction perpendicular to fibers. For intermediate noise intensities, the initial wave broke up into several independent spiral waves or waves circulating around conduction obstacles. At an optimal noise intensity, the number of wavelets was maximized—a form of anti-coherent resonance was obtained. Our results suggest that the inhomogeneity of conduction anisotropy may promote wave breakup and hence play an important role in the initiation and perpetuation of the cardiac arrhythmias. & 2010 Elsevier Ltd. All rights reserved.
Keywords: Inhomogeneous excitable media Spiral wave dynamics Spiral wave breakup Atrial fibrillation Heart arrhythmias
1. Introduction The contraction of the heart muscle is driven by a propagating electrochemical wave of activation. In normal conditions, the activation wave spreads from the sino-atrial node through the heart in a precise, complex manner and ends propagation when there is no more excitable muscle ahead of it [1]. One of pathological forms of wave activity is a spiral wave—a spiral shaped activation wave reentering the same fragment of tissue it has already propagated through [2,3]. Such a perturbation interferes with normal heart rhythm since the spiral wave circulates at a very high rate disturbing the mechanical function of the heart. There are several mechanisms of spiral wave initiation and breakup leading to creation of multiple wavelets [2–11]. In this work, we show a wave breakup caused by the inhomogeneity of the conduction anisotropy. Atrial fibrillation is the most common heart rhythm disorder to affect humans [12]. It has been postulated that one of the mechanisms that leads to the initiation and maintenance of this arrhythmia is the anchoring of spiral waves to myocardial tissue.
n Corresponding author at: Cardiovascular Research Center, Department of Cardiology, L5 McEwin Building, Royal Adelaide Hospital, Adelaide, SA 5000, Australia. Tel.: + 61 8 822 22723; fax: + 61 8 8313 1166. E-mail address: [email protected] (P. Kuklik).
0010-4825/$ - see front matter & 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.compbiomed.2010.07.005
To date, only high density electrical and optical mapping techniques were able to provide some insights into dynamics of spiral waves during atrial fibrillation [13]. However, emerging picture is complex and full description of spiral waves dynamics during atrial fibrillation remains still unknown. A deeper understanding of the spiral wave dynamics may facilitate the design of newer treatment strategies. In this work, we concentrated on spiral wave initiation due to a wave breakup. Wave breakup may occur due to a structural inhomogeneity (e.g. inhomogeneity of the diffusion coefficient [4]) or due the dynamical properties of the system [5]. Several mechanisms of wave breakup due to structural inhomogeneity were recently investigated. Among them: random removal of small fragments of the system [6], perturbation of the position of cells in a cellular automaton model [7], local decrease of the conduction velocity [8], periodic modulation of the conduction velocity [9], inhomogeneity of the density of the system elements [10], periodic modulation of the kinetic parameters [11], inhomogeneity of the conduction velocity [4,14], mechanical deformation of the system [15]. Initiation of the spiral wave due to the inhomogeneity of the tissue fibers orientation was investigated, but only in the three dimensional case [16,17]. In this study, we investigate a role of the inhomogeneous anisotropic uncoupling which would be a model of lateral uncoupling of the cardiac cells due to pathological muscle remodeling.
776
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Conduction of the electrochemical wave of activity in atrial myocardium is determined by several factors including the tissue architecture and the extent and distribution of cell-to-cell connections by gap junction proteins (connexins). Strong coupling between cells enables the rapid transport of the ions between cells and therefore an increase of the wave velocity. Similarly, poor coupling may result in a slowing down of the wave and even in propagation failure. The nature of this electrical transit of conduction makes it sensitive to the distance between the cells and to the number and distribution of connexin proteins. Furthermore, many of the conditions that predispose to the development of arrhythmia are known to remodel the heart muscle leading to a perturbation of wave propagation [18,19]. There are several studies focusing on the effect of the cell coupling on wave conduction [20–24]. Interesting results were obtained by the group of Jacqes de Bakker in an experiment focusing on the mechanisms of the slow conduction in the heart muscle after infarction [20]. A detailed, microscopic study of the patterns of activation revealed a pattern named by the authors a ‘‘zig-zag pattern of activation’’ and characterized by a dissociation of the wave front into independent wavelets propagating through separate muscle strands. This dissociation was possible due to a poor coupling between muscle cells and was caused by an excessive expression of collagen in the infarcted muscle. Such a pattern of propagation was suggested as a possible explanation of the arrhythmogeneity of the post-infarct heart muscle. In this study, we used the FitzHugh–Nagumo model of an excitable medium to investigate the effect of the random perturbation of cell coupling on the stability of a spiral wave. We obtained a form of spatial anti-coherent resonance. Note that in the usual coherent resonance there is an optimal intensity of the external noise at which noise-induced temporal oscillations in the system are most coherent [25]. In our case, the optimal level of the spatial noise added to the cell coupling resulted in an increase of the spatiotemporal complexity of the system. Therefore, we used a term ‘‘anti-coherent resonance’’ to describe the result found in this study.
constant. For lower values of b (b o0.64) the system became oscillatory due to a very low inhibition. For b 40.78 and m 40.36 the size of the meandering pattern of the spiral wave exceeded the size of the sample leading to the termination of the wave. Integration of the FHN equations was carried out in a two dimensional geometry based on a 50 50 rectangular grid using zero-flux boundary conditions. Calculation of the number of the wavelets was conducted on a part of the sample, defined by circle subtended in the center of the sample with radius equal to half of the sample size. The wheel geometry was used to minimize the effect of the sample corners on counting the number of separate wavelets (corners provide a possibility for a wave to be separated into distinct parts and introduce artifacts to waves counting). The Euler integration scheme was used with a time step Dt ¼0.07 and the space step Dx¼ Dy ¼0.4. Since we do not aim to quantitatively reconstruct physiological phenomena, we used dimensionless units in our study. We introduced two diffusion coupling coefficients: for the i,j horizontal direction (Di,j x ) and for the vertical direction (Dy ). As i,j the initial condition, we set Di,j ¼ 0:5 and D ¼ 0:125. Those values x y were set in order to model the anisotropic conduction in the heart tissue. Due to a lack of consensus in the literature (published coupling coefficients ratios range from 1:1 to 9:1 [28,29]), we set the ratio of the diffusion coefficients in the horizontal and perpendicular direction to 1:4. An inhomogeneity of the cells coupling was introduced by addition of spatial, Gaussian noise to Di,j y which would correspond with lateral decoupling of the cells: i,j Di,j y -maxð0,Dy ZÞ
ð3Þ
where Z is given by a Gaussian distribution with a fixed standard deviation s. An example of the distribution of the diffusion coefficients is shown in Fig. 1. With increase in s the amplitude of the coupling in vertical direction (Di,j y ) is decreasing, reaching at
2. Methods We used simplified FitzHugh–Nagumo (FHN) model of an excitable medium [26]. The system was modeled as a twodimensional lattice of the diffusively coupled FHN elements. For an in depth discussion of model properties see [27].The equations of the model were the following: dvi,j vi þ 1,j þvi1,j 2vi,j 1 1 ¼ vi,j v3i,j ui,j þ Di,j x 2 3 dt m Dx i,j vi,j þ 1 þvi,j þ 1 2vi,j ð1Þ þ Dy Dy2
dui,j ¼ mðvi,j þ bgui,j Þ dt
ð2Þ
where vi,j is the transmembrane voltage of the element {i,j}, ui,j is the total slow current, g ¼ 0.5, b ¼0.7, m ¼ 0.3. In order to investigate phenomena in different dynamical states of the FHN model, we varied the parameters b and m in a range in which the propagating spiral wave are sustained. The parameter b responsible for the magnitude of the inhibition current was kept between 0.65 and 0.75 so as to keep m ¼0.3 constant. Alternatively, the parameter m responsible for relative time scales of v and u variables was set between 0.22 and 0.32 to keep b ¼ 0.7
Fig. 1. Example of the spatial distribution of the diffusion coefficients after the addition of Gaussian noise with s ¼0.1. The length of a horizontal line corresponds to the value of Dxi,j and the length of a vertical line with the value of Dyi,j . Introduction of the noise at each grid point perturbs the initial structure of the system.
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
some point a level at which there is a conduction block—a failure of propagation in vertical direction. For in depth discussion of propagation failure in anisotropic media see [30]. Spiral waves were initiated in the system by setting the spatial distribution of v and u to the distribution corresponding with a propagating spiral wave: for each element of the simulation grid, its phase was calculated as
y¼
rkj2 pkn 2 pk
ð4Þ
where r and j are polar coordinates of a given point with respect to the center of the spiral wave, k is its winding number and n is an integer set so as to obtain a phase between 0 and 1. Values of v and u for a given element of the grid were set depending on the phase. y ¼0 corresponded to the resting state and y ¼1 corresponded to the wave front depolarization. All intermediate values of the variables were copied from the stationary distribution of the variables within the action potential cycle. Such prepared distribution of v and u variables (see Fig. 5, t¼5 for example of such distribution) was an initial condition from which spiral wave easily evolved. The number of wavelets was calculated by assessing the number of distinct depolarized areas (i.e. such with v40).
3. Results A low magnitude of the spatial noise (s o0.1) added to the diffusion coefficients did not significantly affect the stability of the spiral wave. Only in few cases the wave broke up (see average number of wavelets for s o0.1 in Figs. 3 and 4). Large magnitude of the spatial noise added to the diffusion coefficients resulted in most cases in the termination of the spiral wave: the wave was unable to rotate due to the high anisotropy level (see average number of wavelets for s 40.5 in Figs. 3 and 4). For intermediate noise intensity, the initial wave broke up into several independent spiral waves or waves circulating around conduction obstacles. Due to the stochastic nature of the inhomogeneity generation, we conducted a series of iterations for a given value of the standard deviation of the noise, s. The number of distinct wavelets as a function of the time is shown in Fig. 2a. The average number of the wavelets calculated in a sliding window
777
20 t.u. (time units) in size saturates after 500 t.u. (we obtained a similar result for all s). Consequently, we fixed end time for the simulation at Tend ¼ 1000 t.u. We performed calculations for s between 0 and 0.7 in increments of 0.02. Next, for a given value of s, we repeated the simulation calculating the average number of wavelets after time Tend. An example of the result for s ¼0.16 is shown in Fig. 2b. For further calculations, we assumed 100 repetitions as sufficient to obtain average wavelet number for given s. For each s, we calculated the average number of distinct wavelets and their standard deviation. The example result is shown in Fig. 3. In this example, for s ¼0.24 the number of wavelets (8.1171.75) attained a maximum—a form of spatial anti-coherent resonance was obtained. Distributions of the number of wavelets for different values of parameters b and m had similar unimodal characteristics (Figs. 3b and 4a, b). Increase of the current inhibition (parameter b) caused a decrease of the mean number of wavelets and a shift of the peak position toward smaller s (Fig. 3b). An increase of the parameter m corresponded with a shift of the peak position from s ¼0.32 (m ¼0.2) to s ¼0.1 (m ¼ 0.32). The maximum number of wavelets was greatest for intermediate values of m (8.23 for m ¼ 0.26 and 8.3 for m ¼0.28) (Fig. 4) and lower for smaller and large values of m. An example of the evolution of the system for s ¼0.2 is shown in Fig. 5. The spiral wave (t ¼5) is unable to rotate initially due to the high anisotropy level resulting in the sliding of the tip to the left and a block of lower part of the arm (t ¼10). After some time, the sliding tip was able to pass the excitation to a row of the elements above, which resulted in the creation of a wavelet propagating backwards (upper arrow, t ¼20). A similar event happened for a separated end of the lower part of the wave (lower arrow, t ¼20). The two additional wavelets evolved further into a more complex organization shown on the remaining images. Wavelets observed in the system fall into two categories: wavelets propagating in loops or small spiral waves rotating in areas able to support wave rotation. An example of a loop is shown in Fig. 6. A loop consists of at least two sets of disconnected rows of elements (with a small value of Dy preventing propagation) ended with a kind of passage consisting of elements with diffusion coefficients allowing propagation in the horizontal and vertical directions. The wavelet in such a loop propagates along one of disconnected rows jumping to other one
Fig. 2. Number of distinct wavelets in the system as a function of the time for s ¼ 0.16, b ¼0.7 and m ¼ 0.3 (a). The average number saturates after t¼ 500 t.u. (also for other values of s). Average number of the distinct wavelets after t¼ 500 t.u. as a function of the number of simulation repetitions for s ¼ 0.16 (b).
778
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Fig. 3. Mean number of wavelets in the system as a function of the standard deviation of noise for b ¼ 0.675 and m ¼ 0.3. For each value of s, the simulations were repeated 100 times and the result averaged. The vertical bars denote the standard deviation of the number of wavelets (a). Mean number of wavelets in the system as a function of the standard deviation of noise for different values of parameter b (with constant m ¼ 0.3) (b).
Fig. 4. Mean number of wavelets in the system as a function of the standard deviation of noise for different values of parameter m (with constant b ¼0.7). Plots for m ¼ 0.2, 0.22, 0.24 and 0.26 (a). Plots for m ¼ 0.28, 0.3 and 0.32 (b).
Fig. 5. Example of the spiral wave breakup for s ¼ 0.15, b ¼ 0.7, m ¼ 0.3. A part of the initial spiral wave arm propagating perpendicularly to the direction of maximal diffusion (vertical) is blocked (t ¼20) leading ultimately to the formation of new spiral or circulating waves.
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
Fig. 6. Example of a loop created by addition of noise to the distribution of the diffusion coefficients (explanation in the text). Grey arrows denote path of the circulating wave of activity.
through the passages at the ends and propagating backward to reenter the first row and repeat the cycle.
4. Discussion The results presented in this study show the relationship between the inhomogeneity of the diffusion anisotropy and the number of the distinct wavelets created after introduction of the spiral wave into system. Our main finding is that the local perturbation of the diffusion anisotropy level promotes spiral wave breakup. We found that there is an optimal perturbation level, which maximizes the number of distinct wavelets. For low and high perturbation levels, we obtained, respectively, a stable propagation and the termination and of the wave. Such behavior may be classified as an anti-coherent resonance due to the inhomogeneity introduced as a spatial noise of the anisotropy level. Termination of the wave in case of high level of spatial noise intensity is to be expected, due to conduction block in vertical direction for too low values of Dy. Such coupling results in a system equivalent to a set of 1D uncoupled cables in which rotation of spiral wave is not possible. A decrease of the mean number of wavelets with b increase (Fig. 3b) may be explained by a less likely propagation of small wavelets along the loops caused by large inhibition current. Increasing inhibition results with total propagation block for b ¼0.75 (Fig. 3b). For lower values of b, the number of the wavelets increases due to lower excitation threshold for lower b-loops unable to sustain wavelets for base parameters (b ¼0.7) due to propagation failure caused by low coupling, are now able to conduct due to increased excitability of the system. Thus, the increase of a number of available loops for propagation results in an increase of wavelets number. A decrease of the mean number of wavelets for smaller m (m o0.24, Fig. 4a) and for larger m (m 40.28, Fig. 4b) was related with different mechanisms. Decrease for large m was related with propagation failure of small wavelets. The onset of the inhibiting current is more rapid in this case causing, in effect, an inhibition of the voltage. In effect, small wavelets are more prone to propagation failure due to an insufficient depolarization area in front of the wave thus decreasing the mean number of wavelets. However, for small values of m the decrease in a number of wavelets was related with a less likely breakup of the wave. Late offset of the inhibition current increased the ability of the wave to propagate through areas of the low coupling thus resulting in more rare incident of the wave breakup and formation of the wavelets. A cause of the spiral wave breakup by conduction block may be either of dynamical origin (e.g. inability to excite due to short interval between consecutive activations) [5] or may be related with a structure of the system (e.g. poor cells coupling) [6–11,14,16,17]. Wave breakup analyzed in this study falls into structural category, related with a modification of the cells
779
coupling—a breakup of the activation wave was caused by propagation failure related to low diffusive coupling. There are several factors resulting in a propagation failure [30] and diffusive coupling used in our system may be regarded as an effective coupling factor affecting wave conduction. Low coupling in vertical direction resulted in a conduction block that either stopped a whole wave (if wavelet size was similar to area of low coupling) or caused annihilation of the part of the wavelet. The mechanism of the wave breakup presented here may be one of the mechanisms of an arrhythmia in diseased heart muscle. Remodeling of the diseased heart muscle leads to an inhomogeneous perturbation of cells coupling [19]. Waves of activity entering such a region may break up and form complex patterns of activation containing spiral waves and wavelets circulating around loops caused by a local separation of the muscle fibers, as shown in de Bakker study reporting ‘‘zig-zag’’ type of propagation [20]. Such type of propagation would correspond with a wavelet in a loop propagating along one of disconnected rows and ‘‘jumping’’ to other one through the connections between rows (example of such loop shown in Fig. 6). Remodeling of the tissue structure and function resulting in ‘‘zig-zag’’ type of activation was caused by infarction in de Bakker paper, but it is possible that other disease related factors may also result in similar type of electrical remodeling [1]. In our previous work, we showed the formation of the spiral waves and waves circulating around non-conducting areas in an isotropic system with an inhomogeneous depression of the diffusion coefficient [4]. The areas of the depressed diffusion formed non-conducting obstacles enabling waves to circulate around such loops. In the system presented in here, we obtained similar propagation properties; however, the loops were formed by a specific inhomogeneity of the lateral coupling while keeping the coupling along the fiber constant. In a heart tissue, such scenario would correspond to an inhomogeneous fibrosis (decrease of the lateral cell coupling with a preserved ability to propagate along the fiber orientation) whereas our previous work would correspond to total uncoupling from the neighboring cells. The frequency of the waves propagating in the loops is a function of the loop length. Shorter loops result in a higher circulation frequency. In a case of the coexistence of several loops (typical for our system), waves emanating from the shortest loop would interfere with the waves circulating around the slower loops. This mechanism was described in [31] where the authors showed for inhomogeneous media that the spiral wave with the highest frequency sweeps away slower spiral waves. However, in our case, the slower activity is not due to the kinetic properties of the model, but rather due to the structural inhomogeneity. In our system, the fastest source of activity will not eliminate slower ones but rather override them with a more frequent train of waves. Thus, after termination of the fastest activity, the wavelets may resume circulation along loops. This property may be relevant for the understanding of the mechanism of atrial fibrillation and radiofrequency ablation procedure. During the ablation procedure, atrial tissue is progressively destroyed. Clinical studies showed, that the frequency of the atrial electrical activity in course of the ablation tends to decrease (ultimately resulting in the restoration of the sinus rhythm) [32]. The mechanism of this phenomenon is still unknown. Our study suggests, that this phenomenon may be explained as progressive removal of the loops formed by fibrotic areas sustaining circulating wavelets. However, the subject requires further study. Decrease of the wavelets number for a large level of noise may lead to a conclusion that noise may have certain anti-arrhythmic role, preventing progression of the arrhythmia into more dangerous and persistent form. Such conclusion has to be formulated very carefully due to overwhelming evidence of a pro-arrhythmic role of
780
P. Kuklik et al. / Computers in Biology and Medicine 40 (2010) 775–780
the heart tissue heterogeneity and fibrosis. Similar conclusions (a warning against oversimplification of role of the heterogeneity in suppressing arrhythmias) were formulated in work of Panfilov et al. [33,34] in which authors studied break up of spiral wave in AlievPanfilov model of excitable media and shown, that addition of non-excitable cells prevents spiral wave from breakup. In our model, the increase of the noise intensity promoted the spiral wave breakup but at some point it caused a propagation termination due to too low lateral coupling to sustain propagation. Since this conclusion contradicts to clinical experience (e.g. [20]), this result should be carefully interpreted. We used the following simplifications in our model:
we used FitzHugh–Nagumo model of the electrical activity of
the heart tissue which is a simplified with respect to the complex behavior of ionic channels; simplified geometry of the system – a two dimensional flat rectangle – which neglects three dimensional nature of the atria, architecture of the muscle fibers, pectinate structure and heterogeneity of the kinetic characteristics of the cells.
The essential features of the mechanism presented in this study comprise:
propagation failure due to low coupling in the direction of propagation;
ability of the system to sustain conduction along rows of laterally uncoupled cells; and
ability of the system to sustain spiral wave. In order to assess whether our simplifications have an effect on the obtained results, all above mechanism features should be independent on the model simplifications. Propagation failure due to low coupling is a generic feature of every excitable system irrespectively on system geometry and dimension (for every configuration of the activation there is a critical value of the coupling that will induce a propagation failure). Ability of the system to sustain conduction along rows of cells (equivalent to 1D cable) is another generic feature of the excitable media reasonably expected to be present in all models of the heart electric activity and irrespectively of model geometry (which is reduced to a fragment of 1D cable by lateral uncoupling). Ability of the system to sustain spiral waves is more complex, since it was shown that both kinetics of the model and system geometry play a role in sustainability of the spiral wave [35]. Therefore, our result is valid only for systems able to sustain spiral wave.
Conflict of interest statement None declared.
Acknowledgement Dr. Sanders is supported by the National Heart Foundation of Australia. References [1] D.P. Zipes, J. Jalife, Cardiac Electrophysiology from Cell to Bedside, W.B. Saunders Company, 2004. [2] E.M. Cherry, F.H. Fenton, Visualization of spiral and scroll waves in simulated and experimental cardiac tissue, New J. Phys. 10 (2008) 125016.
[3] J.M. Davidenko, A.M. Pertsov, R. Salomonsz, W.T. Baxter, J. Jalife, Stationary and drifting spiral waves of excitation in isolated cardiac muscle, Nature 355 (1992) 349–351. [4] P. Kuklik, J.J. Z˙ebrowski, Reentry wave formation in excitable media with stochastically generated inhomogeneities, Chaos 15 (2005) 033301. [5] F.H. Fenton, E.M. Cherry, H.M. Hastings, S.J. Evans, Multiple mechanisms of spiral wave breakup in a model of cardiac electrical activity, Chaos 12 (2002) 852–892. [6] G. Bub, A. Shrier, Propagation through heterogeneous substrates in simple excitable media models, Chaos 12 (2002) 747–753. [7] M. Markus, B. Hess, Isotropic cellular automaton for modelling excitable media, Nature 347 (1990) 56–58. [8] T.J. Lewis, J.P. Keener, Wave-block in excitable media due to regions of depressed excitability, SIAM J. Appl. Math. 61 (1) (2000) 293–316. [9] I. Schebesch, H. Engel, Wave propagation in heterogeneous excitable media, Phys. Rev. E. 57 (4) (1998) 3905–3910. [10] G. Bub, A. Shrier, L. Glass, Spiral wave generation in heterogeneous excitable media, Phys. Rev. Lett. 88 (2002) 5. [11] J. Kneller, R. Zou, E.J. Vigmond, Z. Wang, L.J. Leon, S. Nattel, Cholinergic atrial fibrillation in a computer model of a two-dimensional sheet of canine atrial cells with realistic ionic properties, Circ. Res. 90 (2002) e73–e87. [12] F.D. Murgatroyd, A.J. Camm, Atrial arrhythmiasm, Lancet 341 (1993) 1317–1322. [13] D.S. Rosenbaum, J. Jalife, in: Optical Mapping of Cardiac Excitation and Arrhythmias, Futura Publishing Company, New York, 2001 Armonk. [14] B.E. Steinberg, L. Glass, A. Shrier, et al., The role of heterogeneities and intercellular coupling in wave propagation in cardiac tissue, Philos. Trans. Roy. Soc. A 364 (2006) 1299–1311. [15] H. Zhang, X.S. Ruan, B. Hu, Q. Ouyang, Spiral breakup due to mechanical deformation in excitable media, Phys. Rev. E 70 (2004) 016212. [16] A.V. Panfilov, J.P. Keener, Reentry in three-dimensional Fitzhugh–Nagumo medium with rotational anisotropy, Physica D 84 (1995) 545–552. [17] F. Fenton, A. Karma, Fiber-rotation-induced vortex turbulence in thick myocardium, Phys. Rev. Lett. 81 (1998) 2. [18] P. Kohl, Heterogeneous cell coupling in the heart: an electrophysiological role for fibroblasts, Circ. Res. 93 (2003) 381. [19] M.S. Spach, J.F. Heidlage, P.C. Dolber, R.C. Barr, Electrophysiological effects of remodeling cardiac gap junctions and cell size: experimental and model studies of normal cardiac growth, Circ. Res. 86 (2000) 302. [20] J.M.T. De Bakker, F.J.L. Van Capelle, M.J. Janse, et al., Slow conduction in the infarcted human heart :‘‘zigzag’’ course of activation, Circulation 88 (1993) 915–926. [21] J.E. Saffitz, H.L. Kanter, K.G. Green, T.K. Tolley, E.C. Beyer, Tissue-specific determinants of anisotropic conduction velocity in canine atrial and ventricular myocardium, Circ. Res. 74 (1994) 1065–1070. [22] N.S. Peters, J. Coromilas, N.J. Severs, A.L. Wit, Disturbed connexin43 gap junction distribution correlates with the location of reentrant circuits in the epicardial border zone of healing canine infarcts that cause ventricular tachycardia, Circulation 95 (1997) 988–996. [23] M.S. Spach, Anisotropic structural complexities in the genesis of reentrant arrhythmias, Circulation 84 (1991) 1447–1450. [24] C.Y Chung, H. Bien, E. Entcheva, The role of cardiac tissue alignment in modulating electrical function, J. Cardiovasc. Electrophysiol. 18 (12) (2007) 1323–1329. [25] B. Lindner, J. Garcia-Ojalvo, A. Neimand, L. Schimansky-Geier, Effects of noise in excitable systems, Phys. Rep. 392 (2004) 321–424. [26] R.A. FitzHugh, Impulses and physiological states in theoretical models of nerve membrane, Biophys. J. 1 (1961) 445–466. [27] A.T. Winfree, Varieties of spiral wave behavior: an experimentalist’s approach to the theory of excitable media, Chaos 1 (3) (1991) 303–334. [28] C. Tobon, C. Ruiz, et al., Reentrant mechanisms triggered by ectopic activity in a three-dimensional realistic model of human atrium. a computer simulation study, Comput. Cardiol. 1 (2) (2008) 629–632 2008. [29] G. Seemann, C. Hoper, et al., Heterogeneous three-dimensional anatomical and electrophysiological model of human atria, Philos. Trans. Roy. Soc. A—Math. Phys. Eng. Sci. 364 (1843) (2006) 1465–1481. [30] J.P. Keener, The effects of discrete gap junction coupling on propagation in myocardium, J. Theor. Biol. 148 (1) (1991) 49–82. [31] F. Xie, Z. Qu, J.N. Weiss, A. Garfinkel, Interactions between stable spiral waves with different frequencies in cardiac tissue, Phys. Rev. E. 59 (1999) 2. [32] M. Haissaguerre, P. Sanders, M. Hocini, L.F. Hsu, D.C. Shah, C. Scavee, Y. Takahashi, M. Rotter, J.L. Pasquie, S. Garrigue, J. Clementy, P. Jais, Changes in atrial fibrillation cycle length and inducibility during catheter ablation and their relation to outcome, Circulation 109 (2004) 3007–3013. [33] A.V. Panfilov, Spiral breakup in an array of coupled cells: the role of the intercellular conductance, Phys. Rev. Lett. 88 (2002) 11. [34] K.H.W.J. ten Tusscher, A.V. Panfilov, Influence of nonexcitable cells on spiral breakup in two-dimensional and three-dimensional excitable media, Phys Rev E 68 (2003) 062902. [35] Z. Qu, Critical mass hypothesis revisited: role of dynamical wave stability in spontaneous termination of cardiac fibrillation, Am. J. Physiol. Heart Circ. Physiol. 290 (2006) H255–H263.