Skip to main content
Advertisement
  • Loading metrics

Buffering and total calcium levels determine the presence of oscillatory regimes in cardiac cells

  • Miquel Marchena,

    Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Software, Writing – original draft, Writing – review & editing

    Affiliation Departament de Física, Universitat Politècnica de Catalunya-BarcelonaTech, Barcelona, Spain

  • Blas Echebarria,

    Roles Conceptualization, Formal analysis, Funding acquisition, Investigation, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

    Affiliation Departament de Física, Universitat Politècnica de Catalunya-BarcelonaTech, Barcelona, Spain

  • Yohannes Shiferaw,

    Roles Conceptualization, Formal analysis, Methodology, Validation, Writing – review & editing

    Affiliation Physics Department, California State University, Northridge, California 91330, USA

  • Enrique Alvarez-Lacalle

    Roles Conceptualization, Investigation, Methodology, Supervision, Validation, Writing – original draft, Writing – review & editing

    Affiliation Departament de Física, Universitat Politècnica de Catalunya-BarcelonaTech, Barcelona, Spain

Abstract

Calcium oscillations and waves induce depolarization in cardiac cells which are believed to cause life-threathening arrhythimas. In this work, we study the conditions for the appearance of calcium oscillations in both a detailed subcellular model of calcium dynamics and a minimal model that takes into account just the minimal ingredients of the calcium toolkit. To avoid the effects of homeostatic changes and the interaction with the action potential we consider the somewhat artificial condition of a cell without pacing and with no calcium exchange with the extracellular medium. Both the full subcellular model and the minimal model present the same scenarios depending on the calcium load: two stationary states, one with closed ryanodine receptors (RyR) and most calcium in the cell stored in the sarcoplasmic reticulum (SR), and another, with open RyRs and a depleted SR. In between, calcium oscillations may appear. The robustness of these oscillations is determined by the amount of calsequestrin (CSQ). The lack of this buffer in the SR enhances the appearance of oscillations. The minimal model allows us to relate the stability of the oscillating state to the nullcline structure of the system, and find that its range of existence is bounded by a homoclinic and a Hopf bifurcation, resulting in a sudden transition to the oscillatory regime as the cell calcium load is increased. Adding a small amount of noise to the RyR behavior increases the parameter region where oscillations appear and provides a gradual transition from the resting state to the oscillatory regime, as observed in the subcellular model and experimentally.

Author summary

In cardiac cells, calcium plays a very important role. An increase in calcium levels is the trigger used by the cell to initiate contraction. Besides, calcium modulates several transmembrane currents, affecting the cell transmembrane potential. Thus, dysregulations in calcium handling have been associated with the appearance of arrhythmias. Often, this dysregulation results in the appearance of periodic calcium waves or global oscillations, providing a pro-arrhythmic substrate. In this paper, we study the onset of calcium oscillations in cardiac cells using both a detailed subcellular model of calcium dynamics and a minimal model that takes into account the essential ingredients of the calcium toolkit. Both reproduce the main experimental results and link this behavior with the presence of different steady-state solutions and bifurcations that depend on the total amount of calcium in the cell and in the level of buffering present. We expect that this work will help to clarify the conditions under which calcium oscillations appear in cardiac myocytes and, therefore, will represent a step further in the understanding of the origin of cardiac arrhythmias.

Introduction

Cardiovascular diseases represent one of the main causes of death worldwide [1]. Often, mortality is related to the appearance of rapid cardiac rhythms, such as tachycardia and fibrillation, that result in contractibility loss, reducing cardiac output and eventually leading to sudden cardiac death [2]. Although the onset of rapid arrhythmias can be due to a large variety of factors [3], including changes in the properties of cardiac tissue [4], often arrhythmias are triggered by spontaneous intracellular calcium releases [5, 6]. In cardiac cells, calcium is responsible for regulating cell contraction, but it also modulates several currents that affect the action potential. Thus, spontaneous calcium release in the interbeat interval, during diastole, may elicit extra action potential depolarizations and excitation waves, potentially disrupting normal wave propagation. This sometimes leads to the formation of rotors (functional reentry) and eventually a disordered electrical state characteristic of fibrillation [79].

Often, this focal activity is due to the presence of periodic calcium waves, that result in calcium oscillations [1015]. In paced cardiac cells, oscillations necessarily compete with the external pacing frequency and they may be behind occurrences of spontaneous calcium release events during diastole [16]. Calcium oscillations arise typically due to a malfunction of the Ryanodine Receptor (RyR) [1619], a ligand-gated calcium channel [20] that controls the amplitude of the intracellular calcium transient, by regulating the release of calcium stored at the sarcoplasmic reticulum (SR). Since calcium dynamics in cardiac cells is regulated by the release of calcium at several tens of thousands of RyR clusters (termed calcium release units, CaRUs), global oscillations must appear as a result of an oscillatory regime at the local cluster level that can later be coordinated by diffusion of free calcium. Alternatively, when synchronization is not complete, oscillations at the local level can give rise to periodic calcium waves, providing a pro-arrhythmic substrate [21, 22]. Calcium oscillations have been observed to appear in ventricular myocytes under elevated values of cytosolic calcium [23], due to periodic opening and closing of the RyRs. An increase in cytosolic calcium concentration results in a higher frequency of the oscillations until, at larger values, the SR is depleted because the RyR becomes permanently open [23]. A similar transition has also been studied in models under conditions of SR calcium overload [2426].

A key element in the regulation of intracellular calcium dynamics is Calsequestrin (CSQ). This is one of the major Ca2+-binding proteins in the SR, with high capacity and low affinity [27]. CSQ is preferentially anchored close to the RyR channels [28, 29]. Nowadays it is known that CSQ is an important regulator of RyR gating [3032] through the calcium bound to CSQ. Therefore, CSQ plays two roles in cardiac cells: it can be considered a Ca2+ reservoir for release from the luminal space, but it is also thought to act as a modulator of RyR channel gating. So that, in addition to simply storing Ca2+, CSQ influences the release process by controlling free Ca2+ dynamics near the luminal regulatory sites of the RyR complex [3335]. Dysregulations in CSQ dynamics have been associated with heart diseases [36, 37]. In particular, it has been shown that the absence of CSQ enhances the RyR channel opening and, thus, promotes premature spontaneous SR Ca2+ macro events. The ablation of CSQ causes a form of tachycardia, namely, catecholaminergic polymorphic ventricular tachycardia (CPVT) [36].

In this paper, we use a detailed subcellular calcium model [38] to show the appearance of periodic calcium waves and then analyze this phenomenon using a deterministic model of calcium in a cardiac cell (or in a CaRU). Within this model, we study the existence and stability of different solutions and their dependence on CSQ levels and function. We show that oscillations typically appear at high global calcium concentration and/or high RyR open probability. Their appearance depend on a delicate balance between the total calcium level in the cell and the level of buffering of calcium available. For instance, at high values of CSQ, the system presents a transition from a low concentration, excitable state, to a high concentration state. Such a transition has been proposed to be the basis of complex states, such as long-lasting sparks [39]. At low concentrations of CSQ, in between these two stable states, oscillations appear. We study this transition using a minimal model, that includes the concentration of dyadic and SR calcium and the open probability of the RyR and show that it suffices to explain the appearance of oscillations. A further reduction to a minimal two-dimensional model allows us to explain the transition to the oscillatory regime in terms of the nullcline structure of the system.

Materials and methods

In this paper we introduce two different approaches to understand Ca2+ oscillations in cardiac cells. First, we use a fully detailed subcellular stochastic model of calcium handling to report numerical results showing calcium oscillations. We analyze under which conditions oscillations appear in a controlled scenario where no external pacing is present, and there are no calcium fluxes with the extracellular medium. Later, to gain insight regarding the origin of the oscillations that we observe in the full model, we construct a minimal deterministic model for the local dynamics of calcium at the level of the Calcium Release Unit. The numerical and mathematical analysis of this model allows us to analyze the underlying mechanism of Ca2+ oscillations disregarding the coordination effects of the full model.

Detailed subcellular calcium model

We model the spatial structure of the cell as in a previous model of a cardiomyocyte presented in Marchena and Echebarria [38], which has been modified to add the effects of calsequestrin. A full description of the model is given in the S1 File. We provide now its key features. The equations of the model read: (1) (2) (3) where ci is the calcium concentration in the cytosol, the total calcium concentration in the SR, and cbi represents the concentration of a given buffer in the cytosol. Besides, Jrel and Jup are the release flux from the SR and the uptake by SERCA, respectively, and Jbi represents the binding of free calcium to the different buffers in the cytosol (TnC, SR binding buffer and CaM). These currents are given by: (4) (5) (6)

All the details of the spatial model structure and the values of the parameters can be found in Marchena and Echebarria [38]. The spatial structure of the model includes cytoplasmic and SR spaces, with a spatial discretization of 100 nm. The volume fraction between cytosolic and SR spaces, vi/vsr, is considered to vary spatially, with different values whether the point is close to the z-line or in the inter z-line space. The release flux Jrel carries Ca2+ ions from the SR to the cytoplasm through the RyRs. In our model, the RyR channels, indicated by a yellow dot in Fig 1a, are distributed over the cell along the z-lines, that we identify with periodic narrow strips (0.3 μm width) with a predefined period (Tx). Experimental data shows that the SR domain coincides with these z-lines [40]. A collection of grid points presenting RyRs forms a cluster, i.e., a CaRU. In cardiac cells, CaRUs are arranged periodically in the longitudinal and transversal directions, with some—seemingly Gaussian—dispersion [41]. We place the centers of the clusters on the perimeter following an exact periodic distribution with a period Tx = Ty = 0.5μm. Inside the cell, CaRUs are placed following a Gaussian distribution centered at the z-lines and with a fixed dispersion σ = 0.4μm. The average distance between CaRUs is Tx = 1.6μm and Ty = 0.5μm. We consider that a CaRU contains 36 RyRs, divided equally among 4 grid points, each one containing 9 RyRs. Each RyR can be in one of four states: open (O), close (C) and two inactivated states (I1 and I2) as it is shown in Fig 1b. The transitions among these states is considered to be stochastic. In the release flux, the variable ORyR is the fraction of RyRs that are in the open state and is calculated for all grid points that have a group of RyRs.

thumbnail
Fig 1.

a) RyR distribution in the cell. Each CaRU is formed by four simulation voxels, each one containing 9 RyRs. Thus, all CaRUs are formed by 36 RyRs. The CaRUs are distributed over the cell along the z-lines with a Gaussian distribution in both transversal and longitudinal axes. b) Each RyR follows a four state model, with stochastic transitions among the different states.

https://doi.org/10.1371/journal.pcbi.1007728.g001

We have considered two types of RyR regulation. First, we have considered an open rate and termination dependent on both cytosolic and luminal calcium. The dependence on luminal calcium is due to the presence of CSQ, which mediates its opening, as described in the introduction. However, the core of the paper uses a RyR which is not affected by luminal calcium. One of the main goals of this paper is to test which are the effects of buffering on the appearance of oscillations purely in terms of its binding properties to calcium and not its regulatory effects. Introducing a CSQ-mediated RyR opening would mix different possible mechanisms. For comparison and analysis, we think it is better to leave this dependence initially out. Once the mechanism of CSQ as a buffer is understood we have used section 2 in the S1 File. to show that a full detailed model of RyR with its open rate dependent on SR calcium content gives exactly the same transition to oscillations as described in the main text. If anything, this regulatory mechanism just makes them more robust.

Besides the presence of CSQ in the SR, we also consider in the model the concentrations of other buffers. In particular, TnC, CaM and SR-bound buffers in the cytosol [38]. Due to the addition of CSQ in the model, two parameters have been adjusted from the parameters published in Marchena and Echebarria [38]: the opening rate parameter, now ka = 2.1 ⋅ 10−3 μM−2ms−1, and the dependence of the open probability of the RyR on luminal calcium, now EC50−SR = 450 μM. Contrary to the buffers in the cytosol, the dynamics of CSQ is considered to be faster [4244] than release with a time scale hundreds of time shorter. If we denote by cbSQ the calcium concentration bound to CSQ in the SR, then, the amount of bound calcium is given by: (7)

Assuming fast binding, the stationary condition for cbSQ () is: (8) where KSQ = koffSQ/konSQ is the dissociation constant. From this, the concentration of free calcium can be obtained solving Eq (2) for the total amount of calcium in the SR, (9)

Solving this equation we obtain the value of free calcium in the SR, (10)

The advantage of using this formulation of the rapid buffer approximation over the more usual, for example in [45], is that it conserves mass exactly.

Under physiological conditions, the total amount of calcium in the cell at steady state is fixed by calcium homeostasis, i.e. the complex interaction of LCC, exchanger, and pumps, which affect the steady state level at which the calcium entering the cell balances the calcium extruding. In this work, we are interested in studying instabilities in calcium cycling, under constant cell calcium content. This allows us to focus the analysis on the conditions for the appearance of calcium oscillations under different possible calcium homeostatic levels. Thus, we neglect calcium exchange with the extracellular medium, setting the conductances of the L-type calcium channels and the NCX equal to zero. Then, the total amount of calcium in the cell, QT, is given by: (11)

For a better comparison with the results from a reduced calcium model, described later, we will consider as our control parameter the average calcium content of the cell . Thus, in our simulations, is a constant value that is determined by the initial conditions for cytosolic and luminal calcium (free and bound to buffers).

Reduced calcium model

The minimal model for the local dynamics of calcium is based on the schematics shown in Fig 2. We consider a simplified description of the system, with dynamics of the total calcium concentration in the SR, , and in the cytosolic space close to the RyR2, or dyadic space, cd, and of the open probability of the RyR, Po, (12) (13) (14) with the currents given by (15)

thumbnail
Fig 2. Sketch of the different compartments considered in the simplified model, with the internal variables and the equations of the respective calcium fluxes.

https://doi.org/10.1371/journal.pcbi.1007728.g002

A detailed derivation of these equations and their range of validity can be found in Section 3 of the S1 File. Notice that, as in the full model, we consider a situation where no external pacing is imposed. In this sense, neither external intake from the LCC is considered, nor any extrusion via the sodium-calcium exchanger.

For simplicity, we consider a SERCA pump without an equilibrium condition, that always pumps calcium from the cytosol to the SR (see Eq (15). This gives a trivial solution at ci = cd = 0, instead of the physiological value of ∼ 100nM. However, at these values of cd and ci the luminal calcium is roughly 1mM and, thus, the trivial solution is a reasonable simplification. As in the detailed subcellular model, we assume the approximation of rapid CSQ buffer, so we can compute the amount of free luminal calcium csr from the total luminal calcium from Eq (9).

To close the system we should introduce an extra equation for calcium concentration in the cytosol ci. However, as we assume that the total calcium content in the cell is constant, then we have a conservation equation. Therefore, we can compute ci solving the following quadratic equation for the conservation of (16) where v is the unit volume defined as vvi + vsr + vd, and Bb is the concentration of a generic buffer in the cytosol.

To simplify the analysis, we proceed to work with the assumption that the dynamics of opening of the RyRs is faster than that of calcium concentration (), obtaining then a minimal two-variable model. This will be our baseline minimal model. However, we will later also consider an alternative model with fast dynamics for the dyadic calcium concentration (). Clearly, under the assumption that both processes are fast, as considered often in the literature [46], in the model oscillations would dissappear.

A third theoretical possibility would be to consider fast luminal dynamics. Although analytically possible, this approach does not have much physiological sense, as SERCA is typically slow compared to release or diffusion from the dyadic space.

Fast RyR dynamics.

In this case, we assume that the open and close dynamics of the RyR are fast, so we can assume that it is in a quasi-steady state (). Then, from Eq (14), we obtain: (17) where the parameter is the ratio of the open and close rates of the RyR.

Then, with these assumptions, the simplified model becomes (18) (19)

Fast dyadic calcium dynamics.

In order to test the robustness of the analysis, we also consider a simplified model given by Eqs (12)–(14), in the limit of fast dynamics in the dyadic space and take . Then, from Eq (12): (20)

Substituting this expression in Eqs (13) and (14), we obtain another minimal model, given by (21) (22) where again, ci must be computed solving the quadratic equation for the conservation of mass [Eq (16)]. For simplicity we will consider the case when no calsequestrin is present BSQ = 0.

Results

We first present the results of the numerical simulations of both the full detailed model and the minimal model of calcium cycling. Both produce the same basic scenarios for intracellular calcium dynamics, with three different dynamical behaviors, which we then proceed to analyze. The goal of the development of the minimal model is, precisely, to be able to perform this analytical treatment and check how the behavior depends on total calcium and buffering levels.

Subcellular model

The full detailed model allows us to investigate the different behaviors present in cardiomyocyte calcium cycling when there is no external pacing. We should point out that, under these conditions, the average calcium concentration in the cell is conserved since the total amount of calcium QT in the cell is constant. We produce simulations with different levels of average calcium concentration and observe very different behaviors (Fig 3a). For the lowest value of , the RyR remains almost closed, and most of the calcium content is stored in the SR. Despite the stochasticity of the system, the average values obtained are reproduced reliably with only the presence of local sparks as fluctuations of this global state. This state corresponds to an excitable state, which is the expected behavior of the cell if it has to react properly to external excitation. We call this general state a global shutdown state.

thumbnail
Fig 3.

A: Calcium traces obtained with the full subcellular model and three different values of the average calcium concentration, . B: Line-scans at different values of the load. Increasing the load, the system undergoes a transition from a low cytosolic calcium state (at ), where RyRs remain in the closed state, to spontaneous oscillations, giving rise to calcium waves (). Finally, at high calcium loads () oscillations give rise to a high cytosolic calcium state, where the RyRs remain open, resulting in SR calcium depletion.

https://doi.org/10.1371/journal.pcbi.1007728.g003

When the calcium load increases, the system starts to spontaneously show calcium waves. These waves persist in time with different shapes and durations, giving rise to a nearly periodic oscillation in the global calcium signal. Roughly, we observe one calcium wave per second (Fig 4). Waves are normally initiated at different sites each time but they appear systematically indicating a strong oscillation at the substrate level that we will address in the discussion.

thumbnail
Fig 4. The average period of oscillations at different values of the average calcium concentration , for a concentration of CSQ of BSQ = 2mM (green dots), and in the absence of CSQ (blue dots).

https://doi.org/10.1371/journal.pcbi.1007728.g004

Finally, at large values of oscillations disappear, giving place to a stable state with a high level of calcium in the cytosolic space and a depleted SR. In this state, the RyRs are generally open allowing for the depletion of the SR and the increase of cytosolic calcium. Except for local fluctuations this state is globally stable and we can call it the open state. This state would not respond to external pacing. However, it would produce the activation of the NCX exchanger, which would slowly decrease the average concentration model. As we pointed out previously, the elimination of the calcium intake and extrusion in the model allows us to focus on the general behavior of the cell under different homeostatic scenarios. Numerical simulations indicate that as the calcium level is increased, the cell goes from a shut-down and ready-to-respond state to an oscillatory regime to a global open state where the cell does not respond. Including RyR regulation by SR Ca does not change this scenario, but rather, it actually increases the region where oscillations appear (Fig. 2 in S1 File).

A similar trend has been observed experimentally by Stevens et al [23], even if in the experimental preparation the control parameter was the amount of cytosolic calcium, and not total calcium, as in our simulations. Oscillations appear as the amount of calcium in the cell increases, giving rise to a state with depleted SR calcium (and RyRs in the open state), at high calcium concentrations. Furthermore, experimentally it has been shown that changes in buffering levels can have also important effects on this transition. More specifically, Stevens et. al [23] have shown that the reduction of CSQ in the SR bulk enhances the appearance of oscillations. We have checked if this situation is also present in our simulation and found this to be the case. As shown in Fig 4, when we reduce the CSQ concentration, the oscillations appear at lower values of , they have a higher frequency and the range of oscillations in terms of becomes broader.

Minimal model without calsequestrin

We proceed to explain the results obtained in the minimal model where we find the same qualitative behavior as in the results obtained with the full subcellular stochastic model. In this respect, the minimal three variable model, Eqs (12)–(14), reproduces the appearance of oscillations (S1 Fig). However, since we are interested in using these models to obtain a qualitative understanding of the instability, we have rather considered the reduction to two variables. We have performed simulations in the approximation of fast RyR open dynamics at different values of the cell average calcium concentration . We consider first the case when no calsequestrin is present, BSQ = 0. As we observed in the full subcellular model, at low values of the system remains in a low concentration steady state (Fig 5). In this state the system is excitable, so the fixed point is locally stable, but a large enough perturbation produces an increase in calcium concentration, resulting in a calcium transient typical from CICR. On the other hand, at high calcium loads, there is a new fixed point with high cytosolic calcium concentration, that has been related to the appearance of long-lasting sparks [39]. As in the subcellular model, in between the low concentration fixed point and the permanently open state, the system presents oscillations (Fig 5), that are stable for a quite broad range of loads.

thumbnail
Fig 5. Time traces of the different calcium concentrations for different values of total calcium concentration, , and calsequestrin concentration set to zero (BSQ = 0).

After a transient, the system ends up in either a steady state which is excitatory at low levels of total calcium in the cell with observed low levels of calcium in the cytosol, in an oscillatory state with intermediate levels of total calcium in the cell, or in a state of high total levels of calcium in the cell with observed high cytosolic calcium levels.

https://doi.org/10.1371/journal.pcbi.1007728.g005

Indeed, the number of stationary solutions changes with the calcium concentration . The fixed points of the system can be found from: (23) (24) together with: (25)

Eqs (23)–(25) represent three algebraic equations that give the concentrations ci, cd and csr as a function of total calcium concentration in the cell . When no calsequestrin is present, BSQ = 0, it is easy to obtain that (26) (27)

Introducing these expressions into Eq (23) we obtain an equation of the form . For each value of we can obtain the values of ci that solve the equation. For instance, for a global average calcium concentration of we have only one solution, as shown in Fig 6a. This solution, given by cd = ci = 0 and , exists for all values of . At high values of the average concentration, , another two solutions appear, as depicted in Fig 6b and 6c.

thumbnail
Fig 6. Plot of the function f(ci) for different values of the concentration.

a) At low concentrations there is a single fixed point. b) At higher concentrations two extra unstable fixed points appear. c) At high concentrations the upper fixed point becomes stable.

https://doi.org/10.1371/journal.pcbi.1007728.g006

To calculate the stability of the stationary solutions, we have computed the value of the eigenvalues of the Jacobian matrix, corresponding to Eqs (18) and (19). We find that, while the lower branch is always stable, the other branch of solutions is unstable for a large range of parameters (Fig 7), due to the appearance of oscillations. The stability of the corresponding periodic orbit has been calculated using XPPAUT [47] (Fig 8). We obtain that, as is increased, a limit cycle appears in a global homoclinic bifurcation, with zero frequency (Fig 8c). Increasing , this limit cycle finally disappears in a Hopf bifurcation, at which the upper fixed point becomes stable (Fig 8b).

thumbnail
Fig 7. Solutions for cytosolic calcium concentration, ci, as a function of total calcium concentration, .

Discontinuous lines represent unstable solutions while continuous lines stable ones.

https://doi.org/10.1371/journal.pcbi.1007728.g007

thumbnail
Fig 8.

a) Solutions for cytosolic calcium concentration, ci, as a function of total calcium concentration, . Discontinuous lines represent unstable solutions and continuous lines stable ones. A closer look at the transitions is shown in b). When reducing the total concentration, at , a limit cycle emerges in a Hopf bifurcation, from the upper state, that then becomes unstable. The red lines represent the lower and upper values of the limit cycle. At, , the intermediate unstable fixed point collides with the limit cycle, that disappears in a homoclinic bifurcation. Below , the RyR close state is the only solution. c) Oscillation periods as a function of .

https://doi.org/10.1371/journal.pcbi.1007728.g008

Minimal model with calsequestrin

We present now the main results of the simulation in the fast RyR minimal model when calsequestrin is present. We use Eqs (18) and (19), with Eq (9) that relates free with CSQ-bound SR calcium concentrations. Typically KSQ = 650μM and BSQ somewhere between zero and 20 mM. We have then calculated the different fixed points as a function of the total concentration for different values of BSQ (Fig 9). When BSQ ≠ 0, increasing , the appearance of two extra solutions occurs at larger values of , meaning that the close state solution is stable for a wider range of . Besides, the oscillatory range becomes narrower when BSQ increases, until at certain point oscillations disappear. The disappearance of oscillations is due to the transformation of the bifurcation at which the upper fixed value gains stability, from a Hopf bifurcation into a saddle-node bifurcation.

thumbnail
Fig 9.

a) Number of fixed points and stability of those fixed points as a function of total calcium in the unit, for different values of calsequestrin concentration. In b) and c) details of the transition are shown, for selected values of BSQ.

https://doi.org/10.1371/journal.pcbi.1007728.g009

In addition, at this point, the system presents five fixed points where only the lowest and uppermost are stable. It is important to note that this shows that CSQ enhances the elimination of oscillations. Finally, for high values of , the system has again three fixed points.

Analytical results

The mathematical tractability of the minimal model allows us to get a better understanding of the transition to the upper state via oscillations. A first insight can be obtained by plotting the nullclines of the system (Fig 10), which can help us understand the main mechanisms behind the transition to the oscillatory state. In particular, we obtain the critical average calcium concentration for the onset of oscillations, which depends on buffering levels, and the conditions for the appearance of the upper state.

thumbnail
Fig 10. Structure of the nullclines at different values of indicated in the title of each panel.

The black line indicates the first nullcline , while the orange lines corresponds to . Dots indicate the fixed points. Filled dot: stable fixed point and unfilled dot: unstable fixed point.

https://doi.org/10.1371/journal.pcbi.1007728.g010

Nullclines and stability of solutions.

Besides the transition from one to three solutions (Fig 10a and 10b), nullclines present a clear restructuring of their branches well before the upper state becomes stable. Increasing the total concentration, there is a sudden pinch-off in the cd-nullcline (Fig 10b and 10c). Before this change in nullcline topology, the lower state (with all the calcium in the luminal space and none in the cytosol) is the only stable attractor. Once the cd-nullclines split, oscillations may appear around the upper unstable fixed point.

We can understand the effect of the pinch-off, plotting the trajectories for values of close to the transition. Below the pinch-off, the trajectory follows the fast dynamics of the cd-nullcline (the black line in Figs 10 and 11), until it reaches the fixed point. As the load is increased, the branches of the cd-nullclines come closer until at a certain point the break-up occurs (Fig 11c). Due to the emergence of the pinch-off, the system dynamics follows the lower nullcline up to the tip, at the largest value of csr, where, due to the fast dynamics in the cd direction, it jumps to follow again the nullcline, at larger values of ci. Since there is no stable point, the trajectory starts a persistent oscillation around the unstable fixed point. We can state that, for this problem, the nullcline break-up is the necessary and sufficient condition to obtain oscillations.

thumbnail
Fig 11. Structure of the nullclines at different values of indicated in the title of each panel.

The black line indicates the nullcline , while the orange line corresponds to . The red curve is a trajectory with a direction indicated by the red arrows.

https://doi.org/10.1371/journal.pcbi.1007728.g011

Onset of oscillations.

Using this observation, we can calculate analytically the critical value of beyond which the system oscillates. At the pinch-off of the nullclines, the cubic solution of (black nullcline) loses two of the three solutions at a given csr. To calculate this point analytically, first we observe that the pinch-off occurs at values of cd < Ko (Ko = 15μM). To simplify the calculations, let us make the approximations that, at the pinch off, the cd satisfies cdKo and cdcsr. Being this the case, then Eq (18) reduces to (28)

Furthermore, from Eq (25), we can write ci in terms of csr (assuming BSQ = 0, cdcsr) (29)

Assuming that the concentration of bound cytosolic calcium is much larger than that of free cytosolic calcium, we obtain (30)

From the polynomial equation for cd, Eq (28), solutions of cd are lost at values of csr given by , with . Expanding, this gives the critical value of as: (31)

The same equation can be written as (32)

Once the pinch-off has been produced (Fig 10c), there are two values of , corresponding to the lowest and highest values of the upper and lower nullclines, respectively. Just at the pinch-off these two points merge. This allows to establish the critical value at which the oscillation appears as the one that makes zero the discriminant of the second order polynomial of . After some algebra, the condition for the critical becomes (33) which can be written as (34)

This gives (35)

Using the parameters in Table 2 in S1 File, this expression gives a critical value for the onset of oscillations at around , that given the approximations considered, agrees quite well with the value obtain from the simulations of . Thus, when the total calcium content exceeds this critical value , corresponding to a calcium concentration in the SR of (in the lower state), the system starts to oscillate, at a value of diastolic SR calcium load, given by: (36) which gives a value of mM. At the onset of oscillations, there is thus a sudden decrease in basal SR calcium concentration, to less than half its previous value before the oscillations.

An increase in the quantity of cytosolic buffers (higher Bb) results in a delay in the onset of oscillations, that would occur at higher calcium load (Fig 12a). A higher calcium affinity (lower Kb) on the contrary, would result in oscillations at lower loads (Fig 12b). It is interesting to notice, too, that the strength of SERCA does not affect the onset of oscillations.

thumbnail
Fig 12. Dependency of the onset of oscillations with buffer parameters.

The filled dots represent the control values Bb = 80μM, Kb = 0.5μM, given in Table 2 in S1 File.

https://doi.org/10.1371/journal.pcbi.1007728.g012

Transition to the upper state.

The oscillations disappear at a Hopf bifurcation when the upper state becomes stable. It is possible to relate this transition to the structure of the nullclines in Fig 10. For that, let us recover the definition of the Jacobian matrix J: (37)

A fixed point will be stable provided f1 + g2 < 0. When f1g2f2 g1 > 0 the stable fixed point corresponds to a node and if f1g2f2g1 < 0 to a stable spiral. We can use this to relate the slope of the nullclines to the stability of the upper fixed point. Let us denote and the time derivatives of the independent variables.

At large values of the total concentration (see Fig 13a, with ), the slope of the black nullcline (α = 0) at the fixed point is positive, while the slope of the orange nullcline (β = 0) is negative. Then, increasing cd at constant csr, α goes from being positive to negative. This means that f1 ≡ ∂α/∂cd < 0. Using the same argument, it is easy to check that also g2 ≡ ∂β/∂csr < 0, and therefore the fixed point is stable. At lower values of the total concentration (see, for instance, Fig 13a, for ) the slope of the black nullcline becomes negative. Thus, in this case, while g2 is still negative, f1 becomes positive, and it is not possible to determine the stability of the fixed point. It will depend on the speed of rate of cd and csr close to the fixed point. If the dynamics of cd is faster, then one expects this point to be unstable, if it is csr that varies fast, then stable. One would then expect that buffers that change the dynamics either in the cytosol or in the SR would effect the stability of the fixed point and, therefore, the range of existence of the limit cycle.

thumbnail
Fig 13. Structure of the nullclines at two different values of a) , b) .

A black line corresponds to the nullcline , while the orange line indicates . The functions f1 and g2 are the elements of the diagonal of the Jacobian matrix defined as and .

https://doi.org/10.1371/journal.pcbi.1007728.g013

Robustness of the results

Fast dyadic calcium dynamics.

It is useful to check that the basic points of our discussion hold when different possible approximations are applied to the minimal model. Namely, if the time scale of RyR opening is not as fast as the time scale of calcium diffusion near the dyadic space we should analyze the fast dyadic calcium approximation and not the fast RyR opening approximation to obtain information from the nullcline analysis. To this end, we have performed simulations of Eqs (21) and (22) at different values of the total calcium concentration (Fig 14) and test that we find the same basic behavior: at low values of the system remains in a low concentration steady state (Fig 15), but oscillations appear for a range of , up to a limit where an upper state becomes stable.

thumbnail
Fig 14. Solutions as a function of total calcium concentration , with calsequestrin concentration set to zero when the fast dyadic approximation is used.

We obtain the same type of structure as expected. The system can be in a monostable state, which is excitatory (low load), in an oscillatory state (intermediate load), or in a bistable state (high loads), where it usually ends in a state of open RyR and depleted SR calcium concentration.

https://doi.org/10.1371/journal.pcbi.1007728.g014

thumbnail
Fig 15. Structure of the nullclines at different values of indicated in the title of each panel.

The black line indicates the nullcline , while the orange line corresponds to . The red curve is a trajectory with a direction indicated with the red arrows.

https://doi.org/10.1371/journal.pcbi.1007728.g015

More importantly, the basic structure of the nullclines determines the possible solutions and again, oscillations appear when pinch off is produced (Fig 14). The fixed points in this case are the same as in Fig 7, since they correspond to fixed points of Eqs (12)–(14). However, different slaving conditions may change the stability of the fixed points, that now are analyzed in the plane (csr, Po). Similarly to what we found in the previous analysis, calculating the stability, we find that the intermediate state is always unstable while the upper branch is stable above certain value of . Both states appear at around , indicating the robustness of our analysis.

Stochastic effects.

In real systems, the opening and closing of the RyR2 channels presents intrinsic stochastic dynamics. This is also the case in the full subcellular model. When considering average variables overall the cell, most of this stochasticity is averaged out. However, even if small, stochastic effects may play an important role [48] and, in particular, affect the properties of oscillations [49, 50]. We have thus studied these effects by including stochasticity in the minimal model (18) and (19), adding to the open probability (Po) a small random fluctuation, such that: (38) where σ is the strength of the noise and U is a random number between 0 and 1. The magnitude of σ has been adjusted to be much smaller than Po when the system is in the open state or oscillates but large enough to be able to affect the dynamics close to the homoclinic bifurcation. In Fig 16a the calcium traces for both domain (cytosol and SR) are plotted. For high values of the load, the system oscillates in the same way as in the deterministic limit (). In the deterministic model, the homoclinic bifurcation appears at , so it is reasonable to expect the same behavior for loads above this threshold. Below it, the oscillations in the stochastic system persist and the period increases as the load is decreased, as is typical for homoclinic bifurcations in the presence of noise [51]. For low values of () the period has changed an order of magnitude (∼ 1s) and now the calcium traces resemble closer those observed in the subcellular model, or in experiments [23]. The dependence of the period on the total load is shown in Fig 16b.

thumbnail
Fig 16.

a) Traces of the cytosol and SR calcium concentrations for three different values of total calcium concentration . b) Oscillation periods as a function of . The strength of the noise is σ = 2 ⋅ 10−3.

https://doi.org/10.1371/journal.pcbi.1007728.g016

Discussion and conclusion

Calcium oscillations play an important role in cardiac cells, from the regulation of growth in human cardiac progenitor cells [52], to the control of the pacemaker rhythm in both early embryonic heart cells [53, 54] and in sinoatrial nodal pacemaker cells (SANCs) [5557]. Calcium oscillations have been observed under conditions of high cytosolic calcium concentration [23] or SR calcium overload. High levels of cytosolic calcium affect the opening probability of the RyR, which may result in oscillations or in a permanently open state [23, 58]. Calcium overload can be obtained, for instance, by inhibition of the Na+-K+ pump current INaK, that results in [Na+]i overload. The consequent build-up of [Na+]i reduces the effectiveness of the Na+-Ca2+ exchanger at removing calcium from the cell and intracellular calcium concentrations become elevated. A similar effect is observed in models of hypercalcemia [59]. The effect of elevated [Na+]i has been studied in computational models, finding calcium oscillations [60], that, depending on the model, appear via a supercritical Hopf [24] or a homoclinic bifurcation [25, 26].

The instability of the diastolic resting state is thus well-known [6163] and it plays an important role in the initiation of various cardiac arrhythmias [22]. However, the precise mechanistic relationship between Ca and dangerous AP repolarization is still not well understood. One possibility that has been suggested in previous studies is that Ca cycling induces voltage deflections during the AP which can induce a substrate for wave break and reentry [64]. These voltage deflections are referred to as early-after-depolarizations (EADs), which can be caused by a variety of mechanisms, such as reactivation of the L-type Ca current during AP repolarization. An alternative possibility is that EADs may be caused by Ca oscillations that are regulated by Ca buffers such as CSQ. In this picture, calcium oscillations can perturb the AP via voltage sensitive currents such as the L-type Ca current and sodium-calcium exchange. This idea is consistent with several studies in the literature showing that mutations in CSQ cause abnormal Ca cycling rhythms that are directly linked to sudden cardiac death [65, 66]. Thus, a natural extension of our work would be to study the effect of calcium oscillations on the AP, for which we should extend our model to include detailed, and coupled, calcium and membrane dynamics.

We have shown, using a full subcellular model, that under global Ca overload, SR oscillations appear, leading finally to a state with permanently open RyR and depleted SR. In these simulations, oscillations give rise to periodic calcium waves, that propagate along the myocyte. To obtain a better understanding of the origin and parameter dependence of these transitions, we have studied them in a simplified model of the calcium dynamics. Despite not including all the physiological details (or because of that), minimal models are often useful to gain a better understanding of the origin of complex calcium rhythms [67], as oscillations [68]. We have thus analyzed the different transitions within a minimal model, that takes into account the RyR dynamics, as well as fluxes among different compartments. This allows us to explain the origin of oscillations using a nullcline analysis, as well as to give analytic expressions for the different transitions. A key conclusion is that all buffers can affect heavily the dynamics. The effect of buffers on oscillations has been studied previously [45], and experimentally a change in CSQ levels has been observed to alter the range of values of cytosolic calcium at which oscillations are observed [23]. Here we find that an increase in the levels of CSQ prevents altogether the oscillations, obtaining a direct transition to an open state, that also occurs at larger values of total calcium content having, as mentioned, very important consequences for the study of EADs. We have shown that this picture is robust and independent on which is faster, whether the opening of the RyR2 or the diffusion of calcium near it. Although in this paper we have focused on the effect of CSQ, an interesting extension of this work would be to study in more detail the effect of changes in cytosolic buffer dynamics, allowing, for instance, for mobile buffers. This would require to modify the minimal model to account for difussion of the buffer to nearby calcium release units, and we defer it for future work.

The fact that oscillations appear both in a fully detailed stochastic model and in a simple deterministic model seems another proof of the robustness of the behavior. Our results with noise added to the minimal model seem to suggest that allowing the RyRs to behave stochastically could actually increase the parameter region where oscillations appear. This effect has been observed previously, where some amount of noise coming from the stochasticity of the ion channels may sustain oscillations by the process of coherence resonance [69, 70]. Thus, contrary to the results from the deterministic minimal model, where there is a sudden decrease in basal SR calcium load at the onset of oscillation, the addition of noise provides a gradual transition from the resting state to the oscillatory regime. This would agree with the results obtained by Stevens et al [23], where oscillations appeared gradually as cytosolic calcium content was increased, together with a continuous decrease in the level of basal SR. The addition of noise also gives a longer oscillation period, as observed in the subcellular model or in experiments. Even with this limitation, the analysis of the deterministic minimal model is still useful since it determines the minimal period of oscillation and the possibility of having periodic calcium waves. In a cell, when each calcium release unit presents or is very close to the onset of oscillations, periodic calcium waves appear. Under these circumstances, all units have a tendency to fire periodically, with noise making some units have a slightly early release. Diffusion of calcium to neighbors triggers the release of other units that were already close to open again. Waves are thus the coordinating mechanism of the oscillatory substrate of each unit.

Contrary to other scenarios studied for the onset of calcium waves [61], in our case, the rest state is not unstable, but rather it coexists with the oscillatory state, although, as discussed above, small fluctuations may favor the latter. Oscillations appear when the SR Ca flux through the RyRs becomes large, either because of a large luminal calcium load or to sensitization of the RyR due to an increase in diastolic calcium. These conditions can be met from a lack of CSQ, under an increase in global cell calcium content, or upon phosphorylation of the RyR, for instance. Under these circumstances, there is a strong release flux from the SR that results in SR depletion. As the SR is refilled, the SR release flux again increases, resulting in periodic oscillations. Although RyR inactivation is present in the subcellular model, in our simulations release refractoriness is due mostly to SR depletion, and its recovery is dictated by the slow refilling of the SR since this time scale is slower than the RyRs recovery from inactivation. Thus, the latter modulates the conditions of oscillations but does not determine them. We cannot discard, however, that under conditions of fast SR refilling, RyR refractoriness could not become the determining factor of the oscillations. This is reminiscent of the situation observed in calcium alternans, where both a slow refilling of the SR and RyR inactivation have been deemed as possible mechanisms for the onset of that rhythm [71, 72].

Finally, in this paper, we consider a cell which has achieved calcium balance where intrusion and extrusion match, and have neglected calcium fluxes across the cell membrane to focus on the internal calcium dynamics in order to decouple both processes. Under normal pacing, extracellular calcium fluxes typically represent about 10-20% of the total calcium fluxes, so it is not unreasonable to consider that the total calcium content remains constant once at a steady state. Under these conditions, cytosolic and SR calcium concentrations are not independent but linked, and a clear control parameter is the total calcium concentration. Here we show that it uniquely determines the state of the system. Of course, in the presence of transmembrane fluxes, calcium oscillations or waves, or a permanently open RyR, would result in an extrusion of calcium out of the cell and an eventual decrease in the total calcium load of the cell, that would transition back to the quiescent state (in the absence of external pacing). It seems interesting to study in the future the effect of oscillations and waves in the action potential, as well as a paced cell at different periods and the interaction with the pacing period. Observing how the time scale related to oscillations interacts with the time scale needed to extrude the calcium in the cytosol if the open (upper state) is reached may lead to new interesting phenomena.

Supporting information

S1 File. Supporting information file.

In this supplemental information file, we give a full description of the subcellular model, with a table of the different parameters used, as well as the effect of including inactivation in the RyR2. We also show how to obtain the minimal model from a reduction of a compartmental model of a CaRU.

https://doi.org/10.1371/journal.pcbi.1007728.s001

(PDF)

S1 Fig. Time traces of the different calcium concentrations for different values of total calcium concentration, and calsequestrin concentration set to zero (BSQ = 0), obtained simulating the three variable model (12)-(14).

After a transient, the system ends up in either a steady state which is excitatory at low levels of total calcium in the cell with observed low levels of calcium in the cytosol, in an oscillatory state with intermediate levels of total calcium in the cell, or in a state of high total levels of calcium in the cell with observed high cytosolic calcium levels.

https://doi.org/10.1371/journal.pcbi.1007728.s002

(TIF)

S2 Fig. Periodic calcium waves in the subcellular model with (a) and without (b) inactivation.

In the lower panel we show the fraction of RyRs in the different states shown in the schematics of Fig 1b. For the simulations in (b) we have set all the inactivation rates equal to zero.

https://doi.org/10.1371/journal.pcbi.1007728.s003

(TIF)

Acknowledgments

We want to thank I.R.Cantalapiedra, A.Peñaranda and L.Hove-Madsen from fruitful discussions.

References

  1. 1. Mendis S, Puska P, Norrving B, Organization WH, et al. Global atlas on cardiovascular disease prevention and control. Geneva: World Health Organization; 2011.
  2. 2. Huikuri HV, Castellanos A, Myerburg RJ. Sudden death due to cardiac arrhythmias. New England Journal of Medicine. 2001;345(20):1473–1482. pmid:11794197
  3. 3. Qu Z, Weiss JN. Mechanisms of ventricular arrhythmias: from molecular fluctuations to electrical turbulence. Annual review of physiology. 2015;77:29–55. pmid:25340965
  4. 4. Alonso S, Bär M, Echebarria B. Nonlinear physics of electrical wave propagation in the heart: a review. Reports on Progress in Physics. 2016;79(9):096601. pmid:27517161
  5. 5. Song Z, Qu Z, Karma A. Stochastic initiation and termination of calcium-mediated triggered activity in cardiac myocytes. Proceedings of the National Academy of Sciences. 2017;114(3):E270–E279. pmid:28049836
  6. 6. Colman MA. Arrhythmia mechanisms and spontaneous calcium release: Bi-directional coupling between re-entrant and focal excitation. PLoS computational biology. 2019;15(8). pmid:31393876
  7. 7. Hoeker GS, Katra RP, Wilson LD, Plummer BN, Laurita KR. Spontaneous calcium release in tissue from the failing canine heart. American Journal of Physiology-Heart and Circulatory Physiology. 2009;297(4):H1235–H1242. pmid:19648256
  8. 8. Pogwizd SM, McKenzie JP, Cain ME. Mechanisms underlying spontaneous and induced ventricular arrhythmias in patients with idiopathic dilated cardiomyopathy. Circulation. 1998;98(22):2404–2414. pmid:9832485
  9. 9. Zhou S, Chang CM, Wu TJ, Miyauchi Y, Okuyama Y, Park AM, et al. Nonreentrant focal activations in pulmonary veins in canine model of sustained atrial fibrillation. American Journal of Physiology-Heart and Circulatory Physiology. 2002;283(3):H1244–H1252. pmid:12181156
  10. 10. Orchard C, Eisner D, Allen D. Oscillations of intracellular Ca2+ in mammalian cardiac muscle. Nature. 1983;304(5928):735. pmid:6888540
  11. 11. Kort AA, Lakatta EG. Calcium-dependent mechanical oscillations occur spontaneously in unstimulated mammalian cardiac tissues. Circulation research. 1984;54(4):396–404. pmid:6713605
  12. 12. Cheng H, Lederer M, Lederer W, Cannell M. Calcium sparks and [Ca2+] i waves in cardiac myocytes. American Journal of Physiology-Cell Physiology. 1996;270(1):C148–C159.
  13. 13. Marchant JS, Parker I. Role of elementary Ca 2+ puffs in generating repetitive Ca 2+ oscillations. The EMBO Journal. 2001;20(1-2):65–76. pmid:11226156
  14. 14. Uhlén P, Fritz N. Biochemistry of calcium oscillations. Biochemical and biophysical research communications. 2010;396(1):28–32. pmid:20494106
  15. 15. Dupont G, Combettes L, Bird GS, Putney JW. Calcium oscillations. Cold Spring Harbor perspectives in biology. 2011;3(3):a004226. pmid:21421924
  16. 16. Plummer BN, Cutler MJ, Wan X, Laurita KR. Spontaneous calcium oscillations during diastole in the whole heart: the influence of ryanodine reception function and gap junction coupling. American Journal of Physiology-Heart and Circulatory Physiology. 2011;300(5):H1822–H1828. pmid:21378143
  17. 17. Katra RP, Oya T, Hoeker GS, Laurita KR. Ryanodine receptor dysfunction and triggered activity in the heart. American Journal of Physiology-Heart and Circulatory Physiology. 2007;292(5):H2144–H2151. pmid:17189349
  18. 18. Priori SG, Napolitano C, Tiso N, Memmi M, Vignati G, Bloise R, et al. Mutations in the cardiac ryanodine receptor gene (hRyR2) underlie catecholaminergic polymorphic ventricular tachycardia. Circulation. 2001;103(2):196–200. pmid:11208676
  19. 19. Lou Q, Belevych AE, Radwański PB, Liu B, Kalyanasundaram A, Knollmann BC, et al. Alternating membrane potential/calcium interplay underlies repetitive focal activity in a genetic model of calcium-dependent atrial arrhythmias. The Journal of physiology. 2015;593(6):1443–1458. pmid:25384790
  20. 20. Van Petegem F. Ryanodine receptors: structure and function. Journal of Biological Chemistry. 2012;287(38):31624–31632. pmid:22822064
  21. 21. Xie LH, Weiss JN. Arrhythmogenic consequences of intracellular calcium waves. American Journal of Physiology-Heart and Circulatory Physiology. 2009;297(3):H997–H1002. pmid:19561309
  22. 22. Shiferaw Y, Aistrup GL, Wasserstrom JA. Intracellular Ca2+ waves, afterdepolarizations, and triggered arrhythmias; 2012.
  23. 23. Stevens SC, Terentyev D, Kalyanasundaram A, Periasamy M, Györke S. Intra-sarcoplasmic reticulum Ca2+ oscillations are driven by dynamic regulation of ryanodine receptor function by luminal Ca2+ in cardiomyocytes. The Journal of physiology. 2009;587(20):4863–4872. pmid:19703963
  24. 24. Varghese A, Winslow RL. Dynamics of the calcium subsystem in cardiac Purkinje fibers. Physica D: Nonlinear Phenomena. 1993;68(3-4):364–386.
  25. 25. Winslow RL, Varghese A, Noble D, Adlakha C, Hoythya A. Generation and propagation of ectopic beats induced by spatially localized Na–K pump inhibition in atrial network models. Proceedings of the Royal Society of London Series B: Biological Sciences. 1993;254(1339):55–61. pmid:8265676
  26. 26. Benson AP, Holden AV. Calcium oscillations and ectopic beats in virtual ventricular myocytes and tissues: bifurcations, autorhythmicity and propagation. In: International Workshop on Functional Imaging and Modeling of the Heart. Springer; 2005. p. 304–313.
  27. 27. Mitchell RD, Simmerman H, Jones L. Ca2+ binding effects on protein conformation and protein interactions of canine cardiac calsequestrin. Journal of Biological Chemistry. 1988;263(3):1376–1381. pmid:3335548
  28. 28. Jorgensen AO, Campbell KP. Evidence for the presence of calsequestrin in two structurally different regions of myocardial sarcoplasmic reticulum. The Journal of cell biology. 1984;98(4):1597–1602. pmid:6371026
  29. 29. Franzini-Armstrong C, Kenney LJ, Varriano-Marston E. The structure of calsequestrin in triads of vertebrate skeletal muscle: a deep-etch study. The Journal of cell biology. 1987;105(1):49–56. pmid:3497158
  30. 30. Györke I, Hester N, Jones LR, Györke S. The role of calsequestrin, triadin, and junctin in conferring cardiac ryanodine receptor responsiveness to luminal calcium. Biophysical journal. 2004;86(4):2121–2128. pmid:15041652
  31. 31. Knollmann BC, Chopra N, Hlaing T, Akin B, Yang T, Ettensohn K, et al. Casq2 deletion causes sarcoplasmic reticulum volume increase, premature Ca 2+ release, and catecholaminergic polymorphic ventricular tachycardia. The Journal of clinical investigation. 2006;116(9):2510–2520. pmid:16932808
  32. 32. Song L, Alcalai R, Arad M, Wolf CM, Toka O, Conner DA, et al. Calsequestrin 2 (CASQ2) mutations increase expression of calreticulin and ryanodine receptors, causing catecholaminergic polymorphic ventricular tachycardia. The Journal of clinical investigation. 2007;117(7):1814–1823. pmid:17607358
  33. 33. Beard NA, Sakowska MM, Dulhunty AF, Laver DR. Calsequestrin is an inhibitor of skeletal muscle ryanodine receptor calcium release channels. Biophysical journal. 2002;82(1):310–320. pmid:11751318
  34. 34. Beard N, Laver DR, Dulhunty AF. Calsequestrin and the calcium release channel of skeletal and cardiac muscle. Progress in biophysics and molecular biology. 2004;85(1):33–69. pmid:15050380
  35. 35. Beard NA, Casarotto MG, Wei L, Varsányi M, Laver DR, Dulhunty AF. Regulation of ryanodine receptors by calsequestrin: effect of high luminal Ca2+ and phosphorylation. Biophysical journal. 2005;88(5):3444–3454. pmid:15731387
  36. 36. Dirksen WP, Lacombe VA, Chi M, Kalyanasundaram A, Viatchenko-Karpinski S, Terentyev D, et al. A mutation in calsequestrin, CASQ2D307H, impairs Sarcoplasmic Reticulum Ca2+ handling and causes complex ventricular arrhythmias in mice. Cardiovascular research. 2007;75(1):69–78. pmid:17449018
  37. 37. Rizzi N, Liu N, Napolitano C, Nori A, Turcato F, Colombi B, et al. Unexpected structural and functional consequences of the R33Q homozygous mutation in cardiac calsequestrin: a complex arrhythmogenic cascade in a knock in mouse model. Circulation research. 2008;103(3):298–306. pmid:18583715
  38. 38. Marchena M, Echebarria B. Computational model of calcium signaling in cardiac atrial cells at the submicron scale. Frontiers in physiology. 2018;9:1760. pmid:30618786
  39. 39. Song Z, Karma A, Weiss JN, Qu Z. Long-lasting sparks: Multi-metastability and release competition in the calcium release unit network. PLoS computational biology. 2016;12(1):e1004671. pmid:26730593
  40. 40. Soeller C, Crossman D, Gilbert R, Cannell MB. Analysis of ryanodine receptor clusters in rat and human cardiac myocytes. Proceedings of the National Academy of Sciences. 2007;104(38):14958–14963. pmid:17848521
  41. 41. Chen-Izu Y, McCulle SL, Ward CW, Soeller C, Allen BM, Rabang C, et al. Three-dimensional distribution of ryanodine receptor clusters in cardiac myocytes. Biophysical journal. 2006;91(1):1–13. pmid:16603500
  42. 42. Robertson S, Johnson JD, Potter J. The time-course of Ca2+ exchange with calmodulin, troponin, parvalbumin, and myosin in response to transient increases in Ca2+. Biophysical journal. 1981;34(3):559–569. pmid:7195747
  43. 43. Cannell M, Allen D. Model of calcium movements during activation in the sarcomere of frog skeletal muscle. Biophysical Journal. 1984;45(5):913–925. pmid:6733242
  44. 44. Donoso P, Prieto H, Hidalgo C. Luminal calcium regulates calcium release in triads isolated from frog and rabbit skeletal muscle. Biophysical journal. 1995;68(2):507–515. pmid:7696504
  45. 45. Wagner J, Keizer J. Effects of rapid buffers on Ca2+ diffusion and Ca2+ oscillations. Biophysical Journal. 1994;67(1):447–456. pmid:7919018
  46. 46. Hinch R, Greenstein J, Tanskanen A, Xu L, Winslow R. A simplified local control model of calcium-induced calcium release in cardiac ventricular myocytes. Biophysical journal. 2004;87(6):3723–3736. pmid:15465866
  47. 47. Ermentrout B. Simulating, analyzing, and animating dynamical systems: a guide to XPPAUT for researchers and students. vol. 14. Siam; 2002.
  48. 48. Alvarez-Lacalle E, Echebarria B, Spalding J, Shiferaw Y. Calcium alternans is due to an order-disorder phase transition in cardiac cells. Physical review letters. 2015;114(10):108101.
  49. 49. Kummer U, Krajnc B, Pahle J, Green AK, Dixon CJ, Marhl M. Transition from stochastic to deterministic behavior in calcium oscillations. Biophysical journal. 2005;89(3):1603–1611. pmid:15994893
  50. 50. Dupont G, Croisier H. Spatiotemporal organization of Ca2+ dynamics: A modeling-based approach. HFSP journal. 2010;4(2):43–51.
  51. 51. Santos GJE, Rivera M, Eiswirth M, Parmananda P. Effects of noise near a homoclinic bifurcation in an electrochemical system. Physical Review E. 2004;70(2):021103.
  52. 52. Ferreira-Martins J, Rondon-Clavo C, Tugal D, Korn JA, Rizzi R, Padin-Iruegas ME, et al. Spontaneous calcium oscillations regulate human cardiac progenitor cell growth. Circulation research. 2009;105(8):764–774. pmid:19745162
  53. 53. Viatchenko-Karpinski S, Fleischmann B, Liu Q, Sauer H, Gryshchenko O, Ji G, et al. Intracellular Ca2+ oscillations drive spontaneous contractions in cardiomyocytes during early development. Proceedings of the National Academy of Sciences. 1999;96(14):8259–8264. pmid:10393982
  54. 54. Sasse P, Zhang J, Cleemann L, Morad M, Hescheler J, Fleischmann BK. Intracellular Ca 2+ oscillations, a potential pacemaking mechanism in early embryonic heart cells. The Journal of general physiology. 2007;130(2):133–144. pmid:17664344
  55. 55. Lakatta EG, Maltsev VA, Vinogradova TM. A coupled SYSTEM of intracellular Ca2+ clocks and surface membrane voltage clocks controls the timekeeping mechanism of the heart’s pacemaker. Circulation research. 2010;106(4):659–673. pmid:20203315
  56. 56. Vinogradova TM, Zhou YY, Maltsev V, Lyashkov A, Stern M, Lakatta EG. Rhythmic ryanodine receptor Ca2+ releases during diastolic depolarization of sinoatrial pacemaker cells do not require membrane depolarization. Circulation research. 2004;94(6):802–809.
  57. 57. Vinogradova TM, Lyashkov AE, Zhu W, Ruknudin AM, Sirenko S, Yang D, et al. High basal protein kinase A–dependent phosphorylation drives rhythmic internal Ca2+ store oscillations and spontaneous beating of cardiac pacemaker cells. Circulation research. 2006;98(4):505–514. pmid:16424365
  58. 58. Saftenku E, Williams AJ, Sitsapesan R. Markovian models of low and high activity levels of cardiac ryanodine receptors. Biophysical journal. 2001;80(6):2727–2741. pmid:11371448
  59. 59. Alvarez-Lacalle E, Penaranda A, Cantalapiedra IR, Hove-Madsen L, Echebarria B. Effect of RyR2 refractoriness and hypercalcemia on calcium overload, spontaneous release, and calcium alternans. In: Computing in Cardiology 2013. IEEE; 2013. p. 683–686.
  60. 60. Varghese A, Winslow RL. Dynamics of abnormal pacemaking activity in cardiac Purkinje fibers. Journal of theoretical biology. 1994;168(4):407–420. pmid:8072299
  61. 61. Tveito A, Lines GT, Hake J, Edwards AG. Instabilities of the resting state in a mathematical model of calcium handling in cardiac myocytes. Mathematical biosciences. 2012;236(2):97–107.
  62. 62. Hernandez-Hernandez G, Alvarez-Lacalle E, Shiferaw Y. Role of connectivity and fluctuations in the nucleation of calcium waves in cardiac cells. Physical Review E. 2015;92(5):052715.
  63. 63. Shiferaw Y. Nonlinear onset of calcium wave propagation in cardiac cells. Physical Review E. 2016;94(3):032405. pmid:27739779
  64. 64. Weiss JN, Garfinkel A, Karagueuzian HS, Chen PS, Qu Z. Early afterdepolarizations and cardiac arrhythmias. Heart rhythm. 2010;7(12):1891–1899. pmid:20868774
  65. 65. Viatchenko-Karpinski S, Terentyev D, Györke I, Terentyeva R, Volpe P, Priori SG, et al. Abnormal calcium signaling and sudden cardiac death associated with mutation of calsequestrin. Circulation research. 2004;94(4):471–477. pmid:14715535
  66. 66. Terentyev D, Nori A, Santoro M, Viatchenko-Karpinski S, Kubalova Z, Gyorke I, et al. Abnormal interactions of calsequestrin with the ryanodine receptor calcium release channel complex linked to exercise-induced sudden cardiac death. Circulation research. 2006;98(9):1151–1158. pmid:16601229
  67. 67. Cantalapiedra IR, Alvarez-Lacalle E, Peñaranda A, Echebarria B. Minimal model for calcium alternans due to SR release refractoriness. Chaos: An Interdisciplinary Journal of Nonlinear Science. 2017;27(9):093928. pmid:28964152
  68. 68. Sneyd J, Han JM, Wang L, Chen J, Yang X, Tanimura A, et al. On the dynamical structure of calcium oscillations. Proceedings of the National Academy of Sciences. 2017;114(7):1456–1461. pmid:28154146
  69. 69. Schmid G, Goychuk I, Hänggi P. Stochastic resonance as a collective property of ion channel assemblies. EPL (Europhysics Letters). 2001;56(1):22.
  70. 70. Lindner B, Garcıa-Ojalvo J, Neiman A, Schimansky-Geier L. Effects of noise in excitable systems. Physics reports. 2004;392(6):321–424.
  71. 71. Echebarria B, Alvarez-Lacalle E, Cantalapiedra IR, Peñaranda A. Mechanisms Underlying Electro-Mechanical Cardiac Alternans. In: Nonlinear Dynamics in Biological Systems. Springer; 2016. p. 113–128.
  72. 72. Alvarez-Lacalle E, Cantalapiedra IR, Penaranda A, Cinca J, Hove-Madsen L, Echebarria B. Dependency of calcium alternans on ryanodine receptor refractoriness. PloS one. 2013;8(2). pmid:23390511