March 10, 2011
15:57
WSPC/S0218-1274
02863
International Journal of Bifurcation and Chaos, Vol. 21, No. 2 (2011) 587–601 c World Scientific Publishing Company DOI: 10.1142/S0218127411028635
PROPAGATION AND SYNCHRONIZATION OF Ca2+ SPIRAL WAVES IN EXCITABLE MEDIA JUN MA Department of Physics, Central China Normal University, Wuhan 430079, P. R. China School of Science, Lanzhou University of Technology, Lanzhou 730050, P. R. China
[email protected] JUN TANG College of Science, China University of Mining and Technology, Xuzhou 221008, P. R. China
[email protected] CHUN-NI WANG School of Science, Lanzhou University of Technology, Lanzhou 730050, P. R. China
[email protected] YA JIA∗ Department of Physics, Central China Normal University, Wuhan 430079, P. R. China
[email protected] Received April 6, 2010 Ca2+ spiral wave is observed in a large number of cell types. Based on a Ca2+ model presented by Atri et al., the propagation and synchronization of the intracellular Ca2+ spiral waves are numerically studied and some interesting results are obtained. (i) The largest Lyapunov exponents versus bifurcation parameter is calculated to investigate the oscillation of Ca2+ and it is found that almost all of the largest Lyapunov exponents are negative except for some others that are very close to zero (with an order of magnitude about 10−4 to 10−5 ). (ii) The two controllable parameters — the concentration of inositol 1,4,5-trisphosphate (IP3 ) and the diffusion coefficient — play an important role in inducing and supporting the rotating spiral wave, and the critical thresholds for supporting Ca2+ spiral wave are identified by calculating the statistical factor of synchronization in two-dimensional space. (iii) The driving layer (layer-1) can activate the Ca2+ oscillations in the driven layer (layer-2), and Ca2+ spiral wave in the driven layer can also be developed to synchronize the Ca2+ spiral wave in the drive with appropriate coupling intensity no matter whether mono-directional or mutual coupling of Ca2+ is being used. (iv) An optimal synchronization is observed for an intermediate amount of noise in the presence of Gaussian white noise, and for this noise amount, the Ca2+ waves in the two layers reach complete synchronization most quickly, it is a resonance-like behavior due to noise. Keywords: Spiral wave; Ca2+ oscillation; synchronization; factor of synchronization. ∗
Author for correspondence 587
March 10, 2011
588
15:57
WSPC/S0218-1274
02863
J. Ma et al.
1. Introduction The spiral wave is observed in the chemical, physical and biological systems [Winfree, 1987; B¨ar et al., 1994; Wilkins & Sneyd, 1998; Mckenzie & Sneyd, 1998; Xin, 2000]. Many effective schemes have been proposed to remove the spiral wave in excitable and/or oscillatory media [Ouyang & Felesselles, 1996; Aliev & Panfilov, 1996; Holden, 1997; Mikhailov & Showalter, 2006; Ma et al., 2007; Liu & Wang, 2008] since some reliable evidence has confirmed that the spiral wave in cardiac tissue is responsible for certain harmful heart disease, like arrhythmia and it is thought that the breakup of spiral wave in the cardiac tissue can cause a rapid death of heart, which is called ventricular fibrillation. Most of the previous works [Zhang et al., 2006; Chen et al., 2007] focussed on eliminating the spiral wave, preventing breakup of spiral wave and studying the mechanism of instability of spiral wave in the media. It is also confirmed that local periodical forcing scheme, external polarized electric field and spatial parameter perturbation scheme are useful to remove the spiral wave in the media. On the other hand, spiral wave is also observed in neocortex in experiments [Huang et al., 2004; Schiff et al., 2007] and some authors suggested that the spiral wave in neocortex can be studied by using networks of neurons [Perc, 2007; Wang et al., 2008a, 2008b]. In most previous works on spiral wave in networks of neurons, special but appropriate initial values are selected for the membrane potentials and other variables of the neurons are selected from few sites, and the neurons in other sites are fixed at quiescent states. In this way, spiral wave seed is generated and a stable spiral wave is developed to occupy the entire networks of neurons within appropriate topologies and parameter regions. And propagation of spiral wave in the quiescent area of networks of neurons indicates that neuron signals can still be communicated to break through these quiescent domain where neurons are not active, the reason being that spiral wave is self-sustained. Oscillations in the concentration of intracellular Ca2+ are thought to play an important role in many cellular processes, for example, muscle contraction, etc., and the study on Ca2+ oscillation has been investigated extensively [Perc, 2007; Perc & Marhl, 2003a, 2003b; Marhl et al., 2006; Marhl et al., 2008; Wang et al., 2009]. The Ca2+ concentration is often observed to oscillate [Atri et al., 1993] with periods ranging from a few seconds to more than a
minute, and the oscillations seldom occur uniformly throughout the cell, but are organized into repetitive intracellular waves [Shuai & Jung, 2002]. Intracellular Ca2+ spiral wave [Lechleiter et al., 1991; Tang et al., 2008, 2009; Falcke et al., 1999; Falcke et al., 2003; Tang et al., 2008] is often observed in large cells such as Xenopus oocytes. Up to the best of our knowledge, the Ca2+ spiral wave describes the collective behavior of Ca2+ oscillation, and it is useful to explain the role of Ca2+ oscillation in signals communication, etc. [Sneyd et al., 1995; Cornell-Bell et al., 2004]. Up to now, the mechanism and initial steps of Ca2+ response is clear. An agonist activates a series of reactions that results in the production of the messenger called inositol 1,4,5-trisphosphate (IP3 ) [Putney & Bird, 1993] when it binds to its receptor. Furthermore, IP3 diffuses throughout the cell cytoplasm and binds to IP3 receptors in the endoplasmic reticulum (ER). As it is known, the concentration of Ca2+ in the ER is high and thus a steep concentration gradient of Ca2+ between the cytoplasm and ER is generated. IP3 receptors are a certain type of Ca2+ channels, these channels open and release large amount of Ca2+ when IP3 binds to the receptor. Then the released Ca2+ activates the IP3 receptors to release more Ca2+ , and release of Ca2+ through IP3 receptors is terminated when high cytoplasmic Ca2+ concentrations inactivate the receptors. That is to say, Ca2+ has both a positive and a negative feedback effect on the IP3 receptors. There are large amount Ca2+ in the cytoplasm and it is reasonable to study the collective behaviors of Ca2+ by investigating the intracellular Ca2+ spiral wave in the cell. In 1991, Lechleiter et al. reported that intracellular Ca2+ wave in immature Xenopus oocytes showed remarkable spatiotemporal organization. By loading the oocytes with a Ca2+ sensitive dye, releasing IP3 and observing Ca2+ release patterns with a confocal microscope, Lechleiter et al. observed spiral Ca2+ waves in vivo. The crucial feature of Xenopus oocytes that makes these observations possible is their large size. A typical Ca2+ wave often with a width of close to 100 µm cannot be observed in a small cell while it can be observed in Xenopus oocytes often have a diameter larger than 600 µm, thus it offers enough room for a spiral wave to form and develop. Therefore, Xenopus oocytes are often used to model the Ca2+ spiral waves in theoretical study. In this paper, our aim is to study the formation Ca2+ spiral wave and synchronization of Ca2+
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
spiral waves propagation between two layers in a cell due to Ca2+ coupling, and the statistical factor of synchronization is defined to analyze the statistical properties of spiral Ca2+ wave and critical threshold to support a stable rotating Ca2+ spiral wave. Particularly, we focus on the study of the formation and breakthrough of spiral Ca2+ wave in the quiescent area (certain parameters corresponding to nonexcitable area) by selecting some appropriate initial values of Ca2+ concentrations in few grids of the cell. In fact, the IP3 concentration is a changeable parameter in experiments and can be stimulated by binding extracellular agonists such as neurotranmitters to receptors in a membrane, hormones, or more accurately by photorelease of caged IP3 [Li et al., 1998; Sun et al., 1998]. For example, Sneyd et al. discovered that intercellular propagation of the Ca2+ wave is mediated by the movement of inositol 1,4,5-trisphosphate (IP3 ) through gap junctions [Sneyd et al., 1994; Sneyd et al., 1995]. Therefore, it is reasonable to study the transitions of Ca2+ waves when the controllable parameters IP3 and diffusion coefficient are fixed at different values. And the results on the propagation of Ca2+ spiral waves in a cell and synchronization of Ca2+ spiral waves between two layers can shed light on the understanding of how Ca2+ wave signal breaks through the quiescent domain by propagating Ca2+ spiral wave.
2. Theoretical Model In the following studies, the model of Atri et al. [1993], used to describe the intracellular Ca2+ spiral wave in Xenopus oocytes will be used to investigate the local dynamics, development and survival of the spiral wave, and synchronization of spiral waves between two layers (or two cells, each cell can be simplified as a layer). ∂t c = Dc ∇2 c + Jflux − Jpump ; γc ; kγ + c (1 − b)c h; = kf µ(p) b + (k1 + c)
Jpump = Jflux τh
k2 dh = 2 2 2 −h dt k2 + c
(1)
(2)
where c defines the Ca2+ concentration and p represents the concentration of IP3 . The variable h
589
is used to describe the fraction of IP3 receptors that have not been inactivated by Ca2+ . kf is a scaling parameter and it defines the Ca2+ flux when all IP3 receptors are open and activated, µ(p) describes the fraction of receptors that have bound IP3 . Jflux models the flux of Ca2+ through the IP3 receptor and depends on the fact that Ca2+ activates the IP3 receptor quickly and inactivates it with a slower time scale. Jpump denotes the amount of Ca2+ being pumped out of the cytoplasm back into endoplasmic reticulum (ER) or out through the plasma membrane. The binging sites for IP3 and Ca2+ are supposed to be independent and Ca2+ instantaneously activates the IP3 . More often, the parameters in Eq. (1) are selected with the following values. Parameters b = 0.111, k1 = k2 = 0.7 µM, γ = 2 µM/s, kµ = 0.7 µM, kf = 16.2 µM/s, kγ = 0.1 µM, τh = 2 s, and the biological meanings of these parameters can be found in [Atri et al., 1993] in detail.
3. Numerical Results and Discussion In this section, the local dynamics, formation and development of Ca2+ spiral wave, synchronization of spiral waves between two layers will be investigated, respectively. The media is discretized into 200×200 square grids in numerical simulation, time step is selected with 0.001, and no-flux boundary condition is used. In the case of synchronization of Ca2+ spiral waves of two layers, the same parameters are used in both layers, and the Ca2+ concentration coupling-induced synchronization of spiral waves is discussed in detail.
3.1. Local dynamics The parameter µ(p) or IP3 is changeable, and it is reasonable to study the local dynamics of the Ca2+ oscillation in the absence of Ca2+ diffusion. That is to say, parameter µ(p) can be regarded as a bifurcation parameter and it is confirmed that stable limit cycles appear via homoclinic bifurcation by increasing the value of µ(p). In Fig. 1 the Lyapunov exponent spectrums versus µ(p) at fixed time scale parameter τh = 2 and τh = 3 is plotted, respectively. The results in Fig. 1 show that almost all of the Lyapunov exponents are negative at certain fixed parameter µ(p) or IP3 which can induce a stable state, while few largest Lyapunov exponents vary close to zero which indicates a oscillating state.
March 10, 2011
590
15:57
WSPC/S0218-1274
02863
J. Ma et al.
(a) Fig. 1.
(b)
The largest Lyapunov exponent versus µ(p) at fixed time scale parameter (a) τh = 2.0, (b) τh = 3.0.
In our numerical studies, very few largest Lyapunov exponents are very close to zero with an order of magnitude about 10−4 ∼ 10−5 , it cannot be regarded as weak chaos because chaos cannot be observed in a two-variable autonomous system without time delay. The collective dynamics will be different from the local dynamics when the diffusion effect is considered.
and appropriate values for Ca2+ concentrations in a local area. According to local dynamics from the bifurcation and the curve for the largest Lyapunov exponents of Atri model, it is found that the bifurcation parameter µ(p) mainly decides the dynamics of this model. A statistical synchronization factor is defined to measure the statistical property of the spiral wave based on the mean field theory. The statistical factor is described by
3.2. Formation and survival of spiral wave In theoretical studies, spiral wave can be developed by two ways within appropriate parameter regions. The first one is to break a traveling wave to make a spiral seed, another way is to induce a spiral wave by selecting appropriate initial values for the activator and inhibitor of the system. Clearly, Ca2+ signals can be propagated when Ca2+ concentrations are in the oscillating regions. In fact, some Ca2+ concentrations in certain domain can be quiescent in cells and it is interesting to study how Ca2+ signals are communicated and waves are propagated in these domains and/or layers. As it is known, spiral wave is self-sustained and can propagate in the media without source forcing. Therefore, it can be a useful way to realize Ca2+ signals communication for those quiescent Ca2+ domain only when Ca2+ spiral wave is developed and occupies the whole quiescent region. In our numerical results, the Ca2+ concentrations of most sites are in the quiescent region except for local sites that are selected with special values. Up to the best of our knowledge, Ca2+ concentrations in the cell can be heterogeneous in the realistic cell, therefore, it is reasonable to select special
N
N
1 cij F = 2 N
(3)
j=1 i=1
R=
F 2 − F 2 N N 1 2 (cij − cij 2 ) N2
(4)
j=1 i=1
where the variable cij describes the Ca2+ concentration in site (i, j) and cs is the average value of Ca2+ concentration in all sites. N is constant and N 2 = 40 000 sites in two-dimensional space are used in this paper. The average value of Ca2+ concentrations F becomes invariable when the media becomes homogeneous or the concentrations in all sites are changed with the same value, and the factor of synchronization reaches certain constant close to zero. In fact, the calculation of the statistical factor of synchronization requires appropriate transient period to make a stable mean value, therefore, the factor of synchronization is the mean value of the Ca2+ concentrations in all sites with certain transient period, and it can be regarded as a spatiotemporal average value of Ca2+ concentrations. When each site is regarded as an oscillator,
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
Fig. 2. Cotnour line for factor of synchronization R in the phase space of parameter µ(p) versus diffusion coefficient D.
591
it indicates that all the oscillators reach synchronization when the factor of synchronization is very close to zero while other values can predicate generalized or phase synchronization. In the following numerical simulation studies, time step 0.001, the same initial values and no-flux boundary conditions are used. In the reaction–diffusion system, the diffusion coefficient D also plays an important role in supporting a spiral wave. Extensive numerical studies have been carried out to illustrate the factor of synchronization-measured transition from spiral wave to other states in the parameter µ(p) and diffusion coefficient D phase space. In Fig. 2, the abscissa denotes the diffusion coefficient D, longitudinal coordinates mark µ(p), and the factor of synchronization is shown by the contour line.
(a)
(b)
(c)
(d)
Fig. 3. (a) Evolution of average value of Ca2+ versus time, snapshots of Ca2+ at fixed µ(p) = 0.28 and D = 0.25 for (b) t = 150, (c) t = 450 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 0.65) and black minimal (about 0.05) values of Ca2+ , (d) relation between the factor of synchronization R and diffusion coefficient D at fixed µ(p) = 0.28.
March 10, 2011
592
15:57
WSPC/S0218-1274
02863
J. Ma et al.
Furthermore, the development and dynamical evolution of spiral wave of Ca2+ at certain parameter regions are simulated in detail. At first, the case for µ(p) = 0.26 is investigated, and it is found that no spiral seed can be generated at any fixed diffusion coefficient D = 0.15 to D = 1.05. Then the case for µ(p) = 0.28 is investigated, the evolution of average value of Ca2+ , the snapshots of Ca2+ at fixed time and the curve for the factor of synchronization versus diffusion coefficient are plotted in Fig. 3. The results in Fig. 3(a) show that the average Ca2+ concentration will reach a stable value after a transient period for most diffusion coefficients while long period oscillation is observed for the diffusion coefficient D = 0.35. The snapshots in Figs. 3(b)
and 3(c) show that sparse spiral wave is developed at fixed parameter D = 0.35 and µ(p) = 0.28. The curve in Fig. 3(d) shows that the sudden changing point at D = 0.35 indicates a sudden transition from spiral wave to other states. Furthermore, the case for µ(p) = 0.284 is also investigated, and the numerical results are plotted in Fig. 4. The results in Fig. 4(a) show that the average value of Ca2+ versus time are stable and thus no spiral wave can be developed when the diffusion coefficient D exceeds a certain threshold around D = 0.35. Figures 4(b) and 4(c) illustrate the development of spiral wave at fixed µ(p) = 0.284 and D = 0.25. The sudden changing point on the curve in Fig. 4(d) indicates the critical point for
µ
(a)
(b)
(c)
(d)
Fig. 4. (a) Evolution of average value of Ca2+ concentration versus time, snapshots of Ca2+ at fixed µ(p) = 0.284 and D = 0.25 for (b) t = 100, (c) t = 250 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 1) and black minimal (about 0.0) values of Ca2+ , (d) evolution of the factor of synchronization R versus diffusion coefficient D at fixed µ(p) = 0.284.
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
(a)
(b)
(c)
(d)
(e)
(f)
593
µ
(g) Fig. 5. Snapshots of Ca2+ spiral wave developed from the same initial values within a transient period about 450 time units at µ(p) = 0.30, for (a) D = 0.15, (b) D = 0.25, (c) D = 0.35, (d) D = 0.45, (e) D = 0.55, (f) D = 0.65. Gray scale coloring in all panels is linear, white depicting maximal value (about 1.2) and black minimal (about 0.0) values of Ca2+ , (g) evolution of the factor of synchronization R versus diffusion coefficient D at fixed µ(p) = 0.3.
the survival of the spiral wave. We also investigate the case for µ(p) = 0.30 with different diffusion coefficients, and the results are shown in Fig. 5. The results in Figs. 5(a)–5(f) show that the spiral wave can be easily developed to occupy the whole media within the shorter transient period, and the sudden changing point in the curve for the factor of synchronization versus diffusion coefficient indicates that the critical threshold of diffusion coefficient supports a rotating spiral wave. Furthermore, the case for µ(p) = 0.32 with different
diffusion coefficients is also investigated in extensive studies, and it is confirmed that the spiral wave seed is generated and can be developed within appropriate transient period to occupy the entire media with increasing diffusion coefficient D from 0.1 to 0.45. According to the results in Figs. 2–5, it is confirmed that the appropriate value of µ(p) to induce and develop a stable rotating spiral wave is close to µ(p) = 0.30, and sudden changing points in the curve for the factor of synchronization versus diffusion coefficient indicate that the critical parameter
March 10, 2011
594
15:57
WSPC/S0218-1274
02863
J. Ma et al.
supports a Ca2+ spiral wave in the media. As it is known, noise can change the excitable property of the media, and thus spiral wave propagation is also affected. The noise in intracellular Ca2+ system may originate from many complex sources. In this paper, a Gaussian white noise η(x, y, t) is used to simulate the noise sources by adding the Gaussian white noise to the right-hand side of Eq. (1). The statistical property of Gaussian white noise is often described by η(x, y, t) = 0 and η(x, y, t)η(x , y , t ) = 2D0 δ(x − x )δ(y − y )δ(t − t ), where D0 is the intensity of noise. In Fig. 6, the dynamical formation and development of Ca2+ spiral wave in the presence of weak noise is investigated. The results of Fig. 6 confirm that a sparse spiral wave can be generated in the presence of weak Gaussian white noise. And extensive studies show that no spiral wave seed can be generated when the intensity of Gaussian white noise exceeds a certain threshold. Above all, the dynamical formation and evolution of intracellular Ca2+ spiral wave in single cell are studied in detail, it is confirmed that Ca2+ spiral wave can be developed in the area where most Ca2+ concentrations are quiescent only by selecting appropriate initial values for Ca2+ in the few sites of a local area. And these special and appropriate initial values for Ca2+ concentrations are possible due to heterogeneity in realistic cells. In this way, Ca2+ wave is propagated from the region (domain) with
(a)
high Ca2+ concentration to quiescent Ca2+ concentrations regions. In the following subsection, how Ca2+ waves are propagated among two layers (or cells, each cell is thought as a layer for simplicity) will be discussed by measuring the synchronization of Ca2+ spiral waves in two layers (or cells) with Ca2+ concentration coupling.
3.3. Synchronization of Ca2+ spiral waves In this subsection, we will study the propagation of Ca2+ spiral waves between two layers due to complete Ca2+ coupling. The two layers (or cells, each cell can be simplified as a layer from the mathematical viewpoint) are marked by L-1 (drive system) and L-2 (driven system) and monodirectional and mutual coupling between L-1 and L-2 will be investigated, respectively. The same initial values will be used in L-1 (drive system) for generating spiral wave seed, we will discuss the problems on spiral wave propagation among the second layer due to different coupling when the second layer is selected with same initial values and/or random initial values, respectively. At first, monodirectional coupling is studied and the driven L-2 is described by ∂t c2 = Dc ∇2 c2 + Jflux − Jpump + k(c1 − c2 ); Jpump =
γc2 ; kγ + c2
(b)
Fig. 6. Evolution of average value of Ca2+ concentrations versus time, snapshots of Ca2+ concentrations at fixed µ(p) = 0.28, D = 0.25 and noise intensity D0 = 0.0002 for (a) t = 150, (b) t = 450 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 0.7) and black minimal (about 0.0) values of Ca2+ .
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
(a)
595
(b)
(c)
(d)
(e)
(f)
(g)
(h)
Fig. 7. (a) The evolution of errors between average values of Ca2+ concentrations in L-1 and L-2 versus time at different coupling intensity, (b) the evolution of factor of synchronization R versus coupling intensity k. The upper row snapshots illustrate evolution of Ca2+ spiral wave in L-1 for (c) t = 50, (d) 100, (e) 200 time units, and the bottom row snapshots illustrate evolution of Ca2+ spiral wave in L-2 for (f) t = 50, (g) 100, (h) 200 time units. Gray scale coloring in all panels is linear, white depicting maximal value (about 1.0) and black minimal (about 0.0) values of Ca2+ , all the Ca2+ concentrations in L-2 are quiescent. The monodirectional coupling intensity k = 0.1, µ(p) = 0.284 and diffusion coefficient D = 0.25 are used in both cells. The initial Ca2+ concentrations in L-2 are selected with quiescent values.
(1 − b)c2 h2 ; Jflux = kf µ(p2 ) b + (k1 + c2 ) τh
k2 dh2 = 2 2 2 − h2 dt k2 + c2 (5)
where the parameter k denotes the intensity of coupling, all other parameters are given the same values as above. The results of Fig. 7(a) show that the average Ca2+ in the Layer-2 becomes stable with increasing coupling coefficient to threshold k = 0.1, and the
March 10, 2011
596
15:57
WSPC/S0218-1274
02863
J. Ma et al.
mutual coupling of Ca2+ , and the same appropriate initial values for Ca2+ concentration are used in L-1 as above, the initial values of L-2 will be selected with quiescent state and random value, respectively. The results of Fig. 9(a) show that the factors of synchronization reach stable state on increasing the coupling coefficient k, the results of Fig. 9(b) confirm that the errors between the average Ca2+ concentrations of L-1 and L-2 decrease to zero versus time and the transient period for the errors to reach zero becomes shorter with stronger coupling coefficients being used. A new spiral wave is generated in L-2 and it begins to synchronize the spiral wave in L-1 due to Ca2+ concentration coupling, and the complete synchronization is reached when the errors between the average Ca2+ concentrations of L-1 and L-2 become zero. As a result, spiral wave is propagated to another domain by generating a new spiral wave in the cell with appropriate mutual coupling coefficient being used. Clearly, the snapshots of two layers and the error curves for average Ca2+ concentrations in Fig. 9(b) show that Ca2+ spiral wave can be developed completely and the two Ca2+ spiral waves reach complete synchronization within a
snapshots for dynamical evolution of spiral wave in L-2 show that Ca2+ spiral wave is generated and developed due to the monodirectional Ca2+ coupling. The curve in Fig. 7(b) shows that the factor of synchronization will reach a stable value with an increase in coupling coefficient. Extensive studies confirm that shorter transient period is taken to develop a stable rotating Ca2+ spiral wave in L-2 due to monodirectional Ca2+ coupling from L-1 to L-2. Furthermore, the case for smaller coupling coefficient is also investigated and it is found that no stable rotating Ca2+ spiral wave can be developed but a target-like wave in the L-2, and the results are plotted in Fig. 8 when the monodirectional Ca2+ coupling coefficient is selected with k = 0.05. The results in Fig. 8 show that target-like wave can be developed in the driven L-2 due to monodirectional Ca2+ coupling even if the coupling intensity is weak. In this way, Ca2+ signals are propagated from one cell to another one via driving from Ca2+ spiral wave. In fact, mutual coupling of Ca2+ between cells can be another possible way to propagate Ca2+ signals. In the following, we will investigate the case for
(a)
(b)
(c)
(d)
(e)
(f)
Fig. 8. The upper row snapshots illustrate dynamical evolution of Ca2+ spiral wave in L-1 for (a) t = 50, (b) 150, (c) 450 time units, and the bottom row snapshots illustrate evolution of Ca2+ spiral wave in L-2 for (d) t = 50, (e) 150, (f) 450 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 1.0) and black minimal (about 0.0) values of Ca2+ , all the Ca2+ concentrations in L-2 are quiescent. The coupling intensity k = 0.05 , µ(p) = 0.284 and diffusion coefficient D = 0.25 are used in both layers. The initial Ca2+ concentrations in L-2 are selected with quiescent values.
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
(a)
597
(b)
(c1)
(d1)
(e1)
(f1)
(c2)
(d2)
(e2)
(f2)
Fig. 9. (a) The evolution of the factor of synchronization R versus mutual coupling intensity k and (b) the evolution of errors between average value of Ca2+ concentrations in L-1 and L-2 versus time at fixed µ(p) = 0.284, diffusion coefficient D = 0.25. The upper row snapshots illustrate evolution of Ca2+ spiral wave in L-1 for (c1) t = 50, (d1) 200, (e1) 300 and (f1) 450 time units, and the bottom row snapshots illustrate evolution of Ca2+ spiral wave in L-2 for (c2) t = 50, (d2) 200, (e2) 300 and (f2) 450 time units, the mutual coupling begins to work at t = 200 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 1.0) and black minimal (about 0.0) values of Ca2+ , all the Ca2+ concentrations in L-2 are quiescent. The mutual coupling intensity k = 0.2, µ(p) = 0.284 and diffusion coefficient D = 0.25 are used in both layers. The initial Ca2+ concentrations in L-2 are selected with quiescent values.
transient period about 450 time units even if the initial state of Ca2+ concentrations of L-2 are quiescent. The snapshots in Figs. 9(c1)–9(g2) confirm that spiral waves in the L-1 and L-2 are developed and synchronization can be reached by using appropriate coupling intensity, and thus Ca2+ signal propagates and breaks through the quiescent domain. We also investigate the case when the initial state of Ca2+ concentrations in L-2 are random
values, and same appropriate initial values are selected for the Ca2+ concentrations in L-1. The factors of synchronization and errors of average Ca2+ concentrations will be calculated as well. The results in Fig. 10(a) show that the factors of synchronization decrease to a stable value with increasing mutual coupling coefficient to certain value and the curves in Fig. 10(b) illustrate that errors of average Ca2+ concentrations of two layers completely reaching zero versus time with
March 10, 2011
598
15:57
WSPC/S0218-1274
02863
J. Ma et al.
(a)
(b)
(c1)
(d1)
(e1)
(f1)
(c2)
(d2)
(e2)
(f2)
Fig. 10. (a) The evolution of factor of synchronization R versus mutual coupling intensity k and (b) the evolution of errors between average value of Ca2+ concentrations in L-1 and L-2 versus time at fixed µ(p) = 0.284, diffusion coefficient D = 0.25. The upper row snapshots illustrate the evolution of Ca2+ spiral wave in L-1 for (c1) t = 50, (d1) 200, (e1) 300 and (f1) 450 time units, and the bottom row snapshots illustrate evolution of Ca2+ spiral wave in L-2 for (c2) t = 50, (d2) 200, (e2) 300 and (f2) 450 time units, the mutual coupling begins to open at t = 200 time units. Grey scale coloring in all panels is linear, white depicting maximal value (about 1.0) and black minimal (about 0.0) values of Ca2+ , all the Ca2+ concentrations in L-2 are quiescent. The mutual coupling intensity k = 0.2, µ(p) = 0.284 and diffusion coefficient D = 0.25 are used in both cells. The initial Ca2+ concentrations in L-2 are selected with random values.
appropriate mutual coupling coefficients being used. The transient period for errors of average Ca2+ concentrations to reach zero is taken when stronger mutual coupling coefficients are used. The curves in Fig. 10(b) and snapshots of both layers show that Ca2+ spiral wave can be induced in L-2 and synchronization of Ca2+ spiral waves can be realized with appropriate mutual coupling coefficients being used even if random initial values are used in L-2. Comparing the results in Figs. 9 and 10, it is found that Ca2+ spiral waves can be induced and
developed in L-2 due to mutual coupling of Ca2+ concentrations and synchronization of Ca2+ spiral waves of two layers (L-1 and L-2) can be realized completely. And it is independent of the initial values of Ca2+ concentration of L-2. Comparing the results in the case of mutual coupling (Figs. 9 and 10) and the case of monodirectional coupling (Figs. 7 and 8), it is found that Ca2+ spiral wave can be easily developed in L-2 and the synchronization of Ca2+ spiral waves of the two layers (L-1 and L-2) can be easily reached due to mutual coupling.
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
In the case of monodirectional coupling, the Ca2+ spiral wave in L-1 is kept alive and provides continuous drive to activate the Ca2+ oscillations in L-2 and Ca2+ spiral wave is developed in L-2 , and synchronization of Ca2+ spiral waves is realized finally by appropriate coupling coefficients being used. In the case of mutual coupling, the developing Ca2+ spiral wave in L-1 can be affected by the driving from L-2, and no stable spiral wave can be developed in both cells. Finally, we discuss the case for synchronization of Ca2+ oscillation waves in the presence of Gaussian white noise. The factors of synchronization versus different intensities of noise at fixed coupling coefficient k = 0.2 are calculated and the results are plotted in Fig. 11. The results in Fig. 11 show that a maximal factor of synchronization is observed at fixed noise intensity D0 = 0.0003 and coupling coefficient k = 0.2, and it is also confirmed that a shortest transient period for synchronization of Ca2+ waves is taken with the intensity of noise corresponding to the peak value of the curve in Fig. 11, that is to say, an optimal synchronization of Ca2+ waves is found for an intermediate amount of noise, clearly it can be a resonance-like behavior due to noise. Above all, it extensively investigates the propagation and development of Ca2+ spiral wave in each cell and how Ca2+ spiral wave in one layer activates (or synchronizes) the quiescent Ca2+ in another layer to make Ca2+ waves break through the quiescent areas in the cell. It is found that Ca2+ waves can be propagated among layers (or cells when each cell is thought as a layer) due to the Ca2+ waves
Fig. 11. Correlation of factors of synchronization versus the intensity of Gaussian white noise. The initial Ca2+ concentrations in L-2 are selected with quiescent values.
599
coupling and synchronization. In this way, Ca2+ signals are communicated in the cells.
4. Conclusion In summary, propagation of Ca2+ spiral wave in each cell and between two layers (or two cells) are investigated to study the mechanism of Ca2+ signal propagating in the area where most of Ca2+ concentration are quiescent. At first, local dynamics analysis for Ca2+ oscillating is carried out by calculating the largest Lyapunov exponents versus parameter IP3 and it is found that almost all of the largest Lyapunov exponents are negative while very few are very close to zero which indicates possible oscillatory state for Ca2+ concentration. Thus an useful parameter region is provided to simulate the development of spiral wave in stable or quiescent domain. Then a statistical factor of synchronization is defined to measure the statistical property of the Ca2+ spiral wave in the cell. And the sudden changing points in the curve of the factor versus bifurcation parameters (IP3 and diffusion coefficient) indicate the critical threshold for supporting a stable rotating spiral wave. Appropriate initial Ca2+ concentration is used in some area and the Ca2+ concentration in the rest area is quiescent, and a stable rotating spiral wave can be developed for appropriate IP3 and diffusion coefficient. This indicates that Ca2+ signal can break through the quiescent area (parameter regions for stable Ca2+ concentration) in a cell by propagating Ca2+ spiral wave in the media. Furthermore, synchronization of Ca2+ spiral waves in two layers via layer coupling is investigated. In the case of monodirectional Ca2+ coupling, Ca2+ spiral wave in L-1 (layer-1 or system-1) can activate the quiescent Ca2+ and Ca2+ spiral wave can be developed in L-2 (layer-2 or system-2) within appropriate coupling intensity even when the Ca2+ concentrations in all sites of L-2 are quiescent. Target-like Ca2+ wave is observed in L-2 when the monodirectional coupling intensity is weak (k = 0.05) while synchronization of Ca2+ spiral waves is reached when appropriate coupling intensity is used. In the case of mutual Ca2+ coupling, a stable rotating Ca2+ spiral wave is developed in L-1, random initial values (and quiescent values) of Ca2+ concentrations are used in L-2. It is found that Ca2+ spiral wave can be generated in L-2 and the Ca2+ spiral waves reach complete synchronization with increasing mutual coupling
March 10, 2011
600
15:57
WSPC/S0218-1274
02863
J. Ma et al.
intensity to certain threshold. That is to say, a well developed Ca2+ spiral wave in L-1 can activate the Ca2+ oscillations in layer-2 subjected to mutual Ca2+ coupling. On the other hand, Ca2+ spiral wave cannot be developed (for appropriate Ca2+ concentrations in some areas within L-1) if the mutual Ca2+ coupling works at the beginning. Finally, we investigate the synchronization of Ca2+ spiral waves in the presence of Gaussian white noise, an optimal synchronization of Ca2+ spiral wave is observed for an intermediate amount of noise, and for this noise amount, the Ca2+ waves in two layers reach complete synchronization most quickly, it is a resonance-like behavior due to noise. The results in this paper can give some useful clues to understand the development of Ca2+ spiral wave in cell and propagation of spiral wave between different layers or domains.
Acknowledgments This work is partially supported by the National Natural Science Foundation of China under Grant No. 10747005, 10875049 and the Key Project of Chinese Ministry of Education under No. 108096 and the Natural Science of Lanzhou University of Technology under Grant Nos. Q200706.
References Aliev, R. R. & Panfilov, A. V. [1996] “A simple twovariable model of cardiac excitation,” Chaos Solit. Fract. 7, 293–301. Atri, A. et al. [1993] “A single-pool model for intracellular calcium oscillations and waves in the Xenopus laevis oocyte,” Biophys. J. 65, 1727–1739. B¨ar, M., Gottschalk, N. & Ertl, G. [1994] “Spiral waves in a surface reaction: Model calculations,” J. Phys. Chem. B 100, 1202–1204. Chen, W. Q., Zhang, H. & Ying, H. P. [2007] “Transition from Turing stripe patterns to hexagonal patterns induced by polarized electric fields,” J. Chem. Phys. 127, 154708. Cornell-Bell, A. H., Jun, P. & Trinkaus-Randall, V. [2004] “Decoding calcium wave signaling,” Adv. Molec. Cell Biol. 31, 661–687. Falcke, M. et al. [1999]“Impact of mitochondrial Ca2+ cycling on pattern formation and stability,” Biophys. J. 77, 37–44. Falcke, M., Li, Y. & Lechleiter, J. D. [2003] “Modeling the dependence of the period of intracellular Ca2+ waves on SERCA expression,” Biophys. J. 85, 1474– 1481. Holden, A. V. [1997] “The restless heart of a spiral,” Nature 387, 655–657.
Huang, X. Y., Troy, W. C., Yang, Q. et al. [2004] “Spiral waves in disinhibited mammalian neocortex,” J. Neuroscience 24, 9897–9902. Lechleiter, J. et al. [1991] “Spiral calcium wave propagation and annihilation in Xenopus laevis oocytes,” Science 252, 123–126. Li, W. H. et al. [1998] “Cell-permeant caged InsP3 ester shows that Ca2+ spike frequency can optimize gene expression,” Nature 392, 936–941. Liu, F. C. & Wang, X. F. [2008] “Elimination of antispiral waves by local inhomogeneity in oscillatory systems,” Chin. J. Chem. Phys. 21, 575–580. Ma, J., Jin, W. Y. & Li, Y. L. [2007] “Suppression of spiral waves by generating a self-exciting target wave,” Chin. J. Chem. Phys. 20, 53–58. Marhl, M., Perc, M. & Schuster, S. [2006]“A minimal model for decoding of time-limited Ca2+ oscillations,” Biophys. Chem. 120, 161–167. Marhl, M., Gosak, M. & Perc, M. [2008] “Spatiotemporal modelling explains the effect of reduced plasma membrane Ca2+ efflux on intracellular Ca2+ signalling in hepatocytes,” J. Theor. Biol. 252, 419– 426. Mckenzie, A. & Sneyd, J. [1998] “On the formation and breakup of spiral waves of calcium,” Int. J. Bifurcation and Chaos 8, 2003–2012. Mikhailov, A. S. & Showalter, K. [2006] “Control of waves, patterns and turbulence in chemical systems,” Phys. Rep. 425, 79–194. Ouyang, Q. & Felesselles, J. M. [1996] “Transition from spirals to defect turbulence driven by a convective instability,” Nature 379, 143–146. Perc, M. & Marhl, M. [2003a] “Sensitivity and flexibility of regular and chaotic calcium oscillations,” Biophys. Chem. 104, 509–522. Perc, M. & Marhl, M. [2003b] “Resonance effects determine the frequency of bursting Ca2+ oscillations,” Chem. Phys. Lett. 376, 432–437. Perc, M. [2007] “Effects of small-world connectivity on noise-induced temporal and spatial order in neural media,” Chaos Solit. Fract. 31, 280–291. Putney, J. & Bird, G. [1993] “The inositol phosphatecalcium signaling system in nonexcitable cells,” Endocr. Rev. 14, 610–631. Schiff, S. J. et al. [2007] “Dynamical evolution of spatiotemporal patterns in mammalian middle cortex,” Phys. Rev. Lett. 98, 178102. Shuai, J. W. & Jung, P. [2002] “Optimal intracellular calcium signaling,” Phys. Rev. Lett. 88, 68102. Sneyd, J., Charles, A. C. & Sanderson, M. J. [1994] “A model for the propagation of intercellular calcium waves,” Am. J. Physiol. Cell Physiol. 266, C293– C302. Sneyd, J., Keizer, J. & Sanderson, M. [1995] “Mechanisms of calcium oscillations and waves: A quantitative analysis,” FASEB J. 9, 1463–1472.
March 10, 2011
15:57
WSPC/S0218-1274
02863
Propagation and Synchronization of Ca2+ Spiral Waves in Excitable Media
Sneyd, J. et al. [1995] “Intercellular calcium waves mediated by diffusion of inositol trisphosphate: A twodimensional model,” Am. J. Physiol. Cell Physiol. 268, C1537–C1545. Sun, X. P. et al. [1998] “A continuum of IP3 mediated elementary Ca2+ signalling events in Xenopus oocytes,” J. Physiol. 509, 67–80. Tang, J. et al. [2008] “Numerical study of IP3 -dependent Ca2+ spiral waves in xenopus oocytes,” Europhys. Lett. 83, 68001. Tang, J. et al. [2008] “Numerical study of IP3 -induced Ca2+ spiral pattern evolution,” Chin. Phys. B 17, 4100–4106. Tang, J. et al. [2009] “Theoretical study on drift of Ca2+ spiral waves controlled by electric field,” Commun. Theor. Phys. 51, 941–946. Wang, Q. Y. et al. [2008a] “Synchronization transitions on small-world neuronal networks: Effects of
601
information transmission delay and rewiring probability,” Europhys. Lett. 83, 50008. Wang, Q. Y. et al. [2008b] “Delay-enhanced coherence of spiral waves in noisy Hodgkin–Huxley neuronal networks,” Phys. Lett. A 372, 5681–5687. Wang, B. H. et al. [2009] “Spatiotemporal multiple coherence resonances and calcium waves in a coupled hepatocyte system,” Chin. Phys. B 8, 872–880. Wilkins, M. & Sneyd, J. [1998] “Intercellular spiral waves of calcium,” J. Theor. Biol. 191, 299–308. Winfree, A. T. [1987] When Time Breaks Down (Princeton University Press, Princeton, NJ). Xin, H. W. [2000] “Theoretical study on stochastic resonancein chemical systems,” Chin. J. Chem. Phys. 13, 388–405. Zhang, H., Chen, J. X. & Li, Y. Q. [2006] “Control of spiral breakup by an alternating advective field,” J. Chem. Phys. 125, 204503.