Introduction
============

Sudden cardiac death is primarily caused by ventricular fibrillation, a disorganized and irregular electrical rhythm that arises as a consequence of electrical reentry and wavebreaks.[@b1-cmc-suppl.1-2016-001] Cardiac electrical dynamics are governed primarily at the cellular level by the interaction between action potential properties, including action potential duration (APD) restitution, and intracellular calcium (Ca) handling, and at the tissue level, by action potential propagation properties, including conduction velocity (CV) restitution and cell--cell coupling.[@b2-cmc-suppl.1-2016-001]--[@b8-cmc-suppl.1-2016-001] Irregular dynamics at the cellular level can lead to instabilities in cardiac tissue, including alternans, a beat-to-beat alternation in APD and/or the intracellular Ca transient.[@b9-cmc-suppl.1-2016-001],[@b10-cmc-suppl.1-2016-001] Complex patterns in the alternations of the APD and Ca transient can arise at multiple spatial scales.  In the cardiac myocyte, alternans can be electromechanically concordant (EMC), in which the long--short APD alternation accompanies a large--small alternation in the Ca transient, or electromechanically discordant (EMD), in which long--short APD alternation accompanies a small--large alternation in the Ca transient.[@b9-cmc-suppl.1-2016-001],[@b11-cmc-suppl.1-2016-001] In cardiac tissue, alternans can be spatially concordant (SC), in which the long--short alternation of the APD is in phase throughout the tissue on each beat, or spatially discordant (SD), in which tissue regions of long--short alternations are adjacent to regions of short--long alternations.[@b9-cmc-suppl.1-2016-001],[@b12-cmc-suppl.1-2016-001],[@b13-cmc-suppl.1-2016-001] Critically, SD APD alternans is highly arrhythmogenic, promoting a steep gradient in repolarization, which can lead to conduction block and reentry initiation.[@b14-cmc-suppl.1-2016-001]

Coupled maps are a valuable tool for investigating the complex interactions between subcellular-, cellular-, and tissue-level cardiac dynamics.  Early studies focused on the relationship between APD and preceding diastolic interval (DI), the APD restitution curve, and identified APD restitution slope as a critical parameter governing the onset of alternans.[@b15-cmc-suppl.1-2016-001],[@b16-cmc-suppl.1-2016-001] Subsequent experimental and computational work has shown that APD may depend on more than just the previous DI, referred to as short-term memory.[@b17-cmc-suppl.1-2016-001],[@b18-cmc-suppl.1-2016-001] Short-term memory effects may collectively refer to ion channel recovery dynamics, intracellular Ca cycling, and other ionic concentration changes that occur on time scales longer than a single beat.[@b19-cmc-suppl.1-2016-001] Coupled maps have been expanded to account for the influence of memory on APD.[@b20-cmc-suppl.1-2016-001]--[@b24-cmc-suppl.1-2016-001] However, it is difficult to make mechanistic predictions using these models, since memory effects are typically collectively encompassed by one or two variables and a few parameters that govern the time scale and relative influence of memory.  Although it is clear that alternans can arise at the cellular level via instabilities in both electrical and Ca signaling, experimental evidence in recent years has suggested that Ca-driven instabilities may be the primary mechanism underlying alternans formation,[@b4-cmc-suppl.1-2016-001],[@b10-cmc-suppl.1-2016-001],[@b14-cmc-suppl.1-2016-001] demonstrating the need to specifically account for intracellular Ca signaling, in particular sarcoplasmic reticulum (SR) Ca handling, in computational models.  In 2007, Qu et al developed a detailed coupled map model, coupling APD and Ca cycling dynamics, ie, excitation--contraction coupling, including bidirectional APD--Ca coupling and SR Ca uptake and release dynamics, and investigated bifurcations that lead to cellular instabilities, including alternans and chaos.[@b25-cmc-suppl.1-2016-001]

By incorporating coupled maps of cellular dynamics into a model of action potential wave propagation, one-dimensional cable coupled map models have been developed to investigate tissue-level properties that have shown, for example, how a steep CV restitution relationship, ie, the relationship between CV and the preceding DI, can lead to a transition from SC to SD alternans.[@b5-cmc-suppl.1-2016-001],[@b12-cmc-suppl.1-2016-001],[@b26-cmc-suppl.1-2016-001] This arrhythmogenic phenomenon can arise purely via restitution dynamics, in the absence of spatial heterogeneity: if APD alternans arises, via instabilities in calcium or voltage, at an upstream location in the cardiac tissue, a short DI will precede a short APD and a long DI will precede a long DI.  For steep CV restitution, the shorter DI will result in slower electrical wave propagation, lengthening the DI at locations downstream.  For sufficiently steep CV restitution, the lengthened DI will facilitate a long APD (via APD restitution).  Similarly, the upstream long DI results in faster propagation, leading to a shorter DI at downstream locations, and thus a short APD downstream.  Thus, SD alternans arises in the tissue on a beat-to-beat basis. 

Cell coupling effects on repolarization have also been incorporated in one-dimensional cable coupled map models and shown to influence the onset of alternans, even for steep APD restitution.[@b22-cmc-suppl.1-2016-001],[@b27-cmc-suppl.1-2016-001] In nearly all prior studies utilizing one-dimensional cable-coupled maps, APD dynamics have been minimally defined by the APD restitution relationship, with the addition, in some studies, of a short-term memory variable.  In this study, we incorporate a single cell coupled map, with detailed Ca handling and bidirectional APD-Ca coupling, into a one-dimensional cable model of cardiac tissue.  It is well-established that SR Ca handling is perturbed in pathological settings, in particular in heart failure.[@b28-cmc-suppl.1-2016-001]--[@b34-cmc-suppl.1-2016-001] Our analysis focuses on predicting how altered SR Ca handling leads to arrhythmogenic instabilities at the cellular level, ie, EMD alternans, and at the tissue level, ie, SD alternans. 

Methods
=======

Single cell excitation--contraction coupling map model
------------------------------------------------------

In this section, we describe in brief the single myocyte-coupled map, incorporating excitation--contraction coupling, previously described by Qu et al.[@b25-cmc-suppl.1-2016-001] and illustrated in [Figure 1A](#f1-cmc-suppl.1-2016-001){ref-type="fig"}.  As mentioned in the "Introduction" section, APD is often assumed to be a function of the preceding DI via APD restitution.[@b16-cmc-suppl.1-2016-001] In this model, APD is governed by both voltage-dependent recovery kinetics and Ca-dependent processes, ie, $$\begin{matrix}
{a_{n + 1} = f(d_{n}) + \gamma c_{n + 1}^{p}a_{n} + 1,} \\
{= \frac{1}{1 - \gamma c_{n + 1}^{p}}f(d_{n}),} \\
\end{matrix}$$where *a~n~*~+1~ is the APD of the (*n* + 1)th beat, $c_{n + 1}^{p}$ is the peak intracellular Ca on the (*n* + 1)th beat, γ governs Ca--APD coupling, *f*(*d~n~*) is the APD restitution function ([Eq.  A1](#fd23-cmc-suppl.1-2016-001){ref-type="disp-formula"}), *d~n~* is the DI of the *n*th beat, and *an* and *d~n~* are related by *t~n~* = *a~n~* + *d~n~*, where *t~n~* is the interbeat interval (IBI) of the *n*th beat. 

Peak intracellular Ca $\left( c_{n + 1}^{p} \right)$ is given by the sum of the diastolic intracellular Ca concentration (*c~n~*) and the total SR Ca released (*r~n~*~+1~) during the (*n* + 1)th beat, ie, $$c_{n + 1}^{p} = c_{n} + r_{n + 1},$$where *c~n~* is the diastolic Ca in the cytoplasm at the end of the *n*th beat and is given by the difference of the total cellular Ca concentration, *b~n~*, and the SR Ca load, *l~n~*, at the end of the *n*th beat, ie, $$c_{n} = b_{n} - l_{n}.$$

The SR Ca load on the (*n* + 1)th beat (*l~n~*~+1~) is given by the SR Ca load on the *n*th beat (*l~n~*), less the SR Ca release (*r~n~*~+1~), and plus the SR Ca uptake (*u~n~*~+1~), ie, $$l_{n + 1} = l_{n} - r_{n + 1} + u_{n + 1}.$$

SR Ca release on the (*n* + 1)th beat is given by $$r_{n + 1} = \lambda q(d_{n})g(l_{n}),$$where *q*(*d~n~*) relates the restitution properties of SR Ca release ([Eq.  A2](#fd24-cmc-suppl.1-2016-001){ref-type="disp-formula"}), depending on the previous DI, *g*(*l~n~*) relates the dependency of SR Ca release on the current SR Ca load ([Eq.  A3](#fd25-cmc-suppl.1-2016-001){ref-type="disp-formula"}), and *λ* ∈ (0,1\] describes the fraction of SR Ca released.  In the setting of pathological impairment, *λ* \< 1 may occur, eg, due to reduced ryanodine receptor functional expression, open channel gating, or reduced cluster density. 

SR Ca uptake via SERCA pumps on the (*n* + 1)th beat is given by $$u_{n + 1} = \nu u(t_{n + 1})h(c_{n + 1}^{p}),$$where *u*(*t~n~*~+1~) ([Eq.  A4](#fd26-cmc-suppl.1-2016-001){ref-type="disp-formula"}) and $h(c_{n + 1}^{p})$ ([Eq.  A5](#fd27-cmc-suppl.1-2016-001){ref-type="disp-formula"}) relate the dependence of SR Ca uptake on the stimulus period and peak intracellular Ca levels, respectively, and *ν* is the rate of SR Ca uptake. 

Finally, the total cellular Ca on the (*n* + 1)th beat (*b~n~*~+1~) is given by the total Ca on the *n*th beat (*b~n~*), plus the net Ca influx or efflux during the (*n* + 1)th beat, $$b_{n + 1} = b_{n} - \kappa\lbrack c_{n} - c(t_{n + 1})\rbrack + \eta(a_{n + 1} - a_{n}),$$where *c*(*t~n~*~+1~) is the steady-state intracellular Ca for stimulus period *t~n~*~+1~ ([Eq.  A6](#fd28-cmc-suppl.1-2016-001){ref-type="disp-formula"}), *κ* is a rate of Ca accumulation, and *η* is the APD--Ca coupling.  Total Ca increases at shorter stimulus periods.  Positive or negative APD--Ca coupling results in a net Ca influx or efflux, respectively, during an alternating or chaotic state, ie, when *a~n~*~+1~ *≠ a~n~*.  Note that the dynamics of Ca concentrations, ie, *b~n~*, *c~n~*, $c_{n}^{p}$, and *l~n~* do not account for cellular or subcellular compartment volumes nor Ca buffering, and thus are more representative of a cellular or subcellu-lar compartment Ca content.  As a consequence, concentration values presented somewhat differ from typical experimental measurements.  In particular, accounting for the SR--cytoplasm volume ratio would result in larger SR Ca concentrations subsequently presented.  However, for comparison with prior work using the model,[@b25-cmc-suppl.1-2016-001] we have used the original model formulation, with Ca dynamics and concentrations presented as previously described. 

Combining [Eqs.  1](#fd1-cmc-suppl.1-2016-001){ref-type="disp-formula"}--[7](#fd7-cmc-suppl.1-2016-001){ref-type="disp-formula"}, the single myocyte map model can be written succinctly as $$a_{n + 1} = F_{a}(a_{n},c_{n},l_{n},t_{n};\lambda)$$ $$b_{n + 1} = F_{b}(a_{n},b_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda)$$ $$c_{n + 1} = F_{c}(a_{n},b_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda,\nu)$$ $$l_{n + 1} = F_{l}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda,\nu)$$

In the single myocyte, IBI is equal to the stimulus period and assumed to be constant, ie, *t~n~* = *t~n~*~+1~ = *T*.  All model equations and parameters are given in the [Appendix](#app1-cmc-suppl.1-2016-001){ref-type="app"} and [Table A1](#tA1-cmc-suppl.1-2016-001){ref-type="table"}, respectively. 

One-dimensional myocyte cable map model
---------------------------------------

We next extend the single cell-coupled map to a one-dimensional cable model that accounts for APD and CV restitution, excitation--contraction coupling, and cell--cell coupling.  APD, DI, IBI, total Ca, diastolic intracellular Ca, and SR Ca load on the *n*th beat, given by *a~n~*(*x*), *d~n~*(*x*), *t~n~*(*x*), *b~n~*(*x*), *c~n~*(*x*), and *l~n~*(*x*) as before, are now functions of spatial position *x* ∈ \[0, *L*\], where *L* = 10 cm is the cable length ([Fig.  1B](#f1-cmc-suppl.1-2016-001){ref-type="fig"}).  If the cable is stimulated at *x* = 0 and action potentials on the *n*th beat propagate down the cable length with velocity given by the CV restitution curve, *v*(*d~n~*) ([Eq.  A7](#fd29-cmc-suppl.1-2016-001){ref-type="disp-formula"}), ie, the relationship between CV and the preceding DI, then the IBI at position *x* is governed by the following relationship: $$\frac{dt_{n + 1}(x)}{dx} = \frac{1}{v(d_{n + 1}(x))} - \frac{1}{v(d_{n}(x))}.$$

We next augment the equation for APD in the single cell ([Eq.  8a](#fd8-cmc-suppl.1-2016-001){ref-type="disp-formula"}) with additional terms to account for cell coupling effects on repolarization[@b26-cmc-suppl.1-2016-001],[@b27-cmc-suppl.1-2016-001],[@b35-cmc-suppl.1-2016-001]--[@b39-cmc-suppl.1-2016-001]: $$a_{n + 1}(x) = F_{a}(a_{n},c_{n},l_{n},t_{n};\lambda) + \xi^{2}\frac{d^{2}a_{n + 1}}{dx^{2}} - \omega\frac{da_{n + 1}}{dx}$$

The first spatial term, proportional to $\frac{d^{2}a_{n + 1}}{dx^{2}}$ accounts for diffusion of the membrane potential during repolarization, and the second spatial term, proportional to $\frac{da_{n + 1}}{dx}$, accounts for the asymmetry that arises due to pacing one end of the cable.  Finally, [Eqs.  8b](#fd9-cmc-suppl.1-2016-001){ref-type="disp-formula"}--[8d](#fd11-cmc-suppl.1-2016-001){ref-type="disp-formula"}, which do not contain any additional coupling terms since we assume that intercellular Ca coupling is minimal, are similarly defined at each location *x*: $$b_{n + 1}(x) = F_{b}\left( {a_{n}(x),b_{n}(x),c_{n}(x),l_{n}(x),t_{n}(x),t_{n + 1}(x);\lambda)} \right)$$ $$c_{n + 1}(x) = F_{c}\left( {a_{n}(x),b_{n}(x),c_{n}(x),l_{n}(x),t_{n}(x),t_{n + 1}(x);\lambda\nu)} \right)$$ $$l_{n + 1}(x) = F_{l}\left( {a_{n}(x),c_{n}(x),l_{n}(x),t_{n}(x),t_{n + 1}(x);\lambda\nu)} \right)$$

Numerical integration
---------------------

We can numerically integrate [Eq.  9](#fd12-cmc-suppl.1-2016-001){ref-type="disp-formula"} by first defining spatial location *x~j~* = *j*Δ*x*, where Δ*x* = 0.1 mm (the length of a single cell), *j* = 0, 1,...,*m* and *m* = *L*/Δ*x*.  We next define the spatially discretized variables $a_{n}^{j} = a_{n}(x_{j})$, and similarly, $b_{n}^{j} = b_{n}(x_{j})$, $c_{n}^{j} = c_{n}(x_{j})$, $l_{n}^{j} = l_{n}(x_{j})$, and $t_{n}^{j} = t_{n}(x_{j})$.

For a cable stimulated with pacing period *T* at *x* = 0, we define the initial conditions at all spatial locations, ie, $a_{0}^{j}$, $b_{0}^{j}$, $c_{0}^{j}$, $l_{0}^{j}$, for all *j*, by first iterating the single cell map model ([Eq.  8](#fd8-cmc-suppl.1-2016-001){ref-type="disp-formula"}) for 1000 beats at stimulus period *T*.  We then define each variable as the average of the final four beats of their corresponding values in the single cell map model.  We also define $t_{- 1}^{j} = t_{0}^{j} = T$.

For the (*n* + 1)th beat, discretizing the spatial derivatives and imposing no-flux boundary conditions at *x* = 0 and *x* = *L*, we can write and solve [Eq.  9b](#fd13-cmc-suppl.1-2016-001){ref-type="disp-formula"} as a tridiagonal linear system of equations $$M\mathbf{a}_{n + 1} = \mathbf{A}_{*},$$where $$\begin{array}{l}
{\mathbf{a}_{n + 1} = \left( {a_{n + 1}^{- 1},a_{n + 1}^{0},a_{n + 1}^{1},\ldots,a_{n + 1}^{m},a_{n + 1}^{m + 1}} \right)^{T},} \\
{\mathbf{A}_{*} = \left( {0, - F_{a}\left( {a_{n}^{0},c_{n}^{0},l_{n}^{0},t_{n}^{0};\lambda} \right), - F_{a}\left( {a_{n}^{1},c_{n}^{1},l_{n}^{1},t_{n}^{1};\lambda} \right),\ldots, - F_{a}\left( {a_{n}^{m},c_{n}^{m},l_{n}^{m},t_{n}^{m};\lambda} \right),0} \right)^{T},} \\
{M = \left( \begin{array}{llllllll}
1 & 0 & {- 1} & 0 & \cdots & 0 & 0 & 0 \\
K_{1} & K_{2} & K_{3} & 0 & \cdots & 0 & 0 & 0 \\
0 & K_{1} & K_{2} & K_{3} & \cdots & 0 & 0 & 0 \\
 \vdots & ~ & \ddots & \ddots & \ddots & \vdots & \vdots & \vdots \\
0 & ~ & ~ & \ddots & \ddots & K_{3} & 0 & 0 \\
0 & ~ & ~ & ~ & \ddots & K_{2} & K_{3} & 0 \\
0 & \cdots & ~ & ~ & ~ & K_{1} & K_{2} & K_{3} \\
0 & \cdots & ~ & ~ & ~ & 1 & 0 & {- 1} \\
\end{array} \right),} \\
\end{array}$$*K*~1~ = (*ξ*/A*x*)^2^ + *w*/(2Δ*x*), *K*~2~ = −1--2(*ξ*/Δ*x*)^2^, and *K*~3~ = (*ξ*/Δ*x*)^2^−*w*/(2Δ*x*).  Note that *we* have included fictitious boundary grid points $a_{n + 1}^{- 1}$ and $a_{n + 1}^{m + 1}$.

We next use a forward Euler scheme to numerically integrate [Eq.  9a](#fd12-cmc-suppl.1-2016-001){ref-type="disp-formula"}: $$t_{n + 1}^{j + 1} = t_{n + 1}^{j} + \Delta x\left( {\frac{1}{v\left( {t_{n + 1}^{j} - a_{n + 1}^{j}} \right)} - \frac{1}{v\left( {t_{n}^{j} - a_{n}^{j}} \right)}} \right),$$with boundary condition $t_{n + 1}^{0} = T$ imposed at *x* = 0 and where we have used $d_{n + 1}^{j} = t_{n + 1}^{j} - a_{n + 1}^{j}$ and similarly for $d_{n}^{j}$. Finally, we iterate [Eqs.  9c](#fd14-cmc-suppl.1-2016-001){ref-type="disp-formula"}--[9e](#fd16-cmc-suppl.1-2016-001){ref-type="disp-formula"}: $$b_{n + 1}^{j} = F_{b}\left( {a_{n}^{j},b_{n}^{j},c_{n}^{j},l_{n}^{j},t_{n}^{j},t_{n + 1}^{j};\lambda} \right),$$ $$c_{n + 1}^{j} = F_{c}\left( {a_{n}^{j},b_{n}^{j},c_{n}^{j},l_{n}^{j},t_{n}^{j},t_{n + 1}^{j};\lambda,\nu} \right),$$ $$l_{n + 1}^{j} = F_{l}\left( {a_{n}^{j},c_{n}^{j},l_{n}^{j},t_{n}^{j},t_{n + 1}^{j};\lambda,\nu} \right),$$

Results for both single myocyte and cable simulations are presented for 1000 beats.  Loss of capture in the single myocyte model is defined as occurring when DI is less than 2 ms.  Conduction block in the cable simulations is defined when there is a loss of capture at any spatial location.  The incidence of alternans is defined as an APD alternans amplitude greater than 10 ms. 

Results
=======

Single cell excitation--contraction coupling map model
------------------------------------------------------

The single cell map model developed by Qu et al has a wide range of parameters.  Here, we focus on two critical SR Ca properties, specifically SR Ca uptake and release, and then extend the single cell model to a spatially extended cable to explore how these properties alter the formation of alternans and conduction block. 

We iterate the single cell map for 1000 beats and plot APD, peak intracellular Ca, and SR Ca as a function of the stimulus period *T* ([Fig.  2](#f2-cmc-suppl.1-2016-001){ref-type="fig"}).  For SR uptake rate *ν* = 0.7, APD and Ca alternans are present at faster stimulus rates, between *T* of 265 and 320 ms, below which there is a loss of capture ([Fig.  2A](#f2-cmc-suppl.1-2016-001){ref-type="fig"}).  Further, alternans is EMC, ie, a long APD corresponds with a large peak Ca, and vice versa.  SR Ca is out of phase with the peak intracellular Ca, consistent with a high SR Ca load leading to a larger SR Ca release on the subsequent beat.  For a smaller uptake rate, *ν* = 0.4, both the peak intracellular Ca and SR Ca load are reduced ([Fig.  2B](#f2-cmc-suppl.1-2016-001){ref-type="fig"}).  EMC alternans is present at a longer stimulus period, while loss of capture occurs at a shorter stimulus period. 

For a slightly smaller uptake rate, *ν* = 0.28, as *T* decreases, EMC alternans is present at *T* = 295 ms.  However, as *T* decreases further, APD and peak intracellular Ca become out of phase at *T* = 275 ms, such that the long APD corresponds with the small peak Ca, and vice versa, ie, EMD alternans.  EMD alternans is present until the loss of capture at *T* = 235 ms.  Reducing uptake rate further, *ν* = 0.1, alternans is present for *T* between 235 and 285 ms, EMD over the entire range of stimulus rates. 

For all uptake rates, SR Ca load decreases as *ν* decreases and is out of phase with APD.  APD and Ca alternans are present over the widest range of values for *T* for an SR uptake rate *ν* near 0.4.  For *ν* = 0.4, SR Ca load is between 70 and 95 µM, which corresponds to the steepest region of the SR Ca release vs. SR Ca load relationship *g*(*l~n~*) ([Fig.  A1](#fA1-cmc-suppl.1-2016-001){ref-type="fig"}), which has been previously shown to promote Ca alternans via the following mechanism.[@b9-cmc-suppl.1-2016-001] At the end of *n*th beat, if SR Ca load is reduced slightly, the subsequent SR Ca release on the (*n* + 1) th beat is also reduced.  As such, at the end of the (*n* + 1)th beat, the SR Ca load is larger (due to the prior smaller Ca release), and the subsequent SR Ca release on the (*n* + 2)th beat is larger, leading to a further reduction of the SR Ca load.  This alternating pattern is promoted by a steep SR Ca release vs. SR Ca load relationship.  As *ν* increases above or decreases below 0.4, the SR Ca load increases or decreases, respectively, such that the corresponding SR Ca load operates in a shallower portion of the SR Ca release relationship, and as such, alternans is promoted to a lesser extent. 

Importantly, note that both the Ca--APD coupling parameter γ and the APD--Ca coupling parameter *η* are positive, which both promote EMC alternans.  Thus, the incidence of EMD alternans in the setting of reduced SR Ca uptake is an interesting dynamical behavior ([Fig.  3](#f3-cmc-suppl.1-2016-001){ref-type="fig"}).  If APD on the *n*th beat is short, the DI on this beat will be long, and the diastolic Ca at the end of the beat will be low.  Thus, on the (*n* + 1)th beat, APD will be long, due to APD restitution ([Eq.  A1](#fd23-cmc-suppl.1-2016-001){ref-type="disp-formula"}), and SR Ca release will be large due to SR Ca release restitution ([Eq.  A2](#fd24-cmc-suppl.1-2016-001){ref-type="disp-formula"}).  The large SR Ca release can lead to a large peak intracellular Ca ([Eq.  2](#fd2-cmc-suppl.1-2016-001){ref-type="disp-formula"}), such that alternans is EMC ([Fig.  3A](#f3-cmc-suppl.1-2016-001){ref-type="fig"}).  However, when SR uptake is reduced, the SR Ca load is reduced, such that SR Ca release is also reduced, due to the steep dependence on SR Ca load ([Eq.  A3](#fd25-cmc-suppl.1-2016-001){ref-type="disp-formula"}).  Therefore, for reduced SR Ca uptake, on the (*n* + 1)th beat, even though the SR Ca release is relatively large (compared with the prior beat), peak intracellular Ca can still be small, due to the low diastolic Ca at the end of the prior beat, such that alternans is EMD ([Fig.  3B](#f3-cmc-suppl.1-2016-001){ref-type="fig"}).  This situation is further promoted at shorter stimulus periods due to reduced SR Ca restitution. 

We next investigate to what extent SR Ca release parameter *λ* alters the incidence of EMC and discordant alternans.  For large *ν* = 0.7, we find that, for all values of *λ*, there is a transition from regimes of no alternans, EMC alternans, and loss of capture as *T* decreases ([Fig.  4A](#f4-cmc-suppl.1-2016-001){ref-type="fig"}).  Decreasing *λ* leads to alternans and loss of capture at slightly shorter stimulus periods *T*.  For a smaller *ν* = 0.4, decreasing *λ* leads to a more pronounced decrease in the onset of alternans, with a smaller decrease in the period for loss of capture ([Fig.  4B](#f4-cmc-suppl.1-2016-001){ref-type="fig"}).  Critically, there is a transition from EMC to EMD alternans as *λ* decreases.  For a slightly smaller *ν* = 0.28, the transition from EMC to EMD alternans occurs as *λ* or *T* decreases ([Fig.  4C](#f4-cmc-suppl.1-2016-001){ref-type="fig"}), while for *ν* = 0.1, only EMD alternans is present ([Fig.  4D](#f4-cmc-suppl.1-2016-001){ref-type="fig"}). 

In summary, we find in the single cell map model that reduced SR Ca uptake and release can lead to EMD alternans.  In the next section, we extend our analysis to consider the spatially extended one-dimensional cable map model.  We investigate whether the results of the single cell simulations extend to the tissue model and further to what extent the SR Ca properties alter electromechanical and spatial discordance and conduction block. 

One-dimensional myocyte cable map model
---------------------------------------

In [Figure 5](#f5-cmc-suppl.1-2016-001){ref-type="fig"}, we plot APD and peak intracellular Ca as a function of spatial location on even (black) and odd (red) beats, for *T* = 288 ms, varying SR Ca uptake rate *ν*.  For *ν* = 0.7, APD and peak Ca are SD, with a node approximately at *x* = 5 cm, while at essentially all spatial locations, APD and peak Ca are EMC ([Fig.  5A](#f5-cmc-suppl.1-2016-001){ref-type="fig"}).  For smaller *ν* = 0.4, alternans become SC ([Fig.  5B](#f5-cmc-suppl.1-2016-001){ref-type="fig"}), and then transition to SD again for *ν* = 0.35 ([Fig.  5C](#f5-cmc-suppl.1-2016-001){ref-type="fig"}).  For *ν* = 0.3, alternans remain SD but become EMD as well ([Fig.  5D](#f5-cmc-suppl.1-2016-001){ref-type="fig"}). 

This example demonstrates the complex relationship between intracellular and SR Ca handling and APD dynamics at the cellular level and cell coupling and CV restitution at the tissue-level properties, in particular how variation in SR Ca uptake rate can alter spatial and electromechanical concordance and discordance, which we explain as follows.  Decreasing *ν* from 0.7 to 0.4 results in a transition from SD to SC alternans (via a transition at *ν* ≈ 0.672).  Decreasing *ν* decreases SR Ca uptake, and thus the SR Ca load.  Reduced SR Ca load decreases both the large and small peak intracellular Ca level (via the SR Ca release vs. SR Ca load relationship), which in turn shortens the long APD (due to positive Ca--APD coupling).  However, via APD restitution, the short APD is lengthened.  Collectively, this results in a longer *short* DI, attenuating CV slowing via CV restitution, and thus promoting SC alternans via the mechanism described in the "Introduction" section.  Interestingly, decreasing *ν* from 0.4 to 0.35 results in a transition from SC back to SD (via a transition at *ν* ≈ 0.358).  This decrease in *ν* similarly results in reduced SR Ca uptake and load, peak intracellular Ca, and long APD and an increase in the short APD.  However, in this example, the collective dynamics result in a shorter *long* DI, which promotes CV slowing via CV restitution, and thus SD alternans. 

We demonstrate this complexity in another example, in which we plot APD and peak intracellular Ca as a function of spatial location, for *T* = 280 ms, varying SR Ca release parameter *λ*.  For *λ* = 1, alternans is SC and EMC ([Fig.  6A](#f6-cmc-suppl.1-2016-001){ref-type="fig"}).  For smaller *λ* = 0.5, the APD spatial pattern changes only slightly, while peak Ca decreases at all spatial locations ([Fig.  6B](#f6-cmc-suppl.1-2016-001){ref-type="fig"}), such that for *λ* = 0.36, alternans is EMD, while SC ([Fig.  6C](#f6-cmc-suppl.1-2016-001){ref-type="fig"}).  For smaller *λ* = 0.3, alternans transitions to SD ([Fig.  6D](#f6-cmc-suppl.1-2016-001){ref-type="fig"}). 

In [Figure 7](#f7-cmc-suppl.1-2016-001){ref-type="fig"}, we plot the parameter regimes for spatial and electromechanical concordance and discordance for different stimulus rates *T*, and SR Ca uptake and release parameters *ν* and *λ*.  For *ν* = 0.7, as in the single myocyte model, alternans is EMC, for all values of *λ* ([Fig.  7A](#f7-cmc-suppl.1-2016-001){ref-type="fig"}).  Further, as *T* decreases, there is a transition from SC to SD alternans.  There is a small difference between the onsets of alternans and conduction block stimulus period in the cable simulations (black lines), in comparison with the single cell simulations (gray lines). 

For smaller *ν* = 0.5, the parameter regimes are more complex ([Fig.  7B](#f7-cmc-suppl.1-2016-001){ref-type="fig"}).  As in the single cell model, in general, as *λ* decreases, EMD alternans is promoted, while SD alternans is promoted at shorter values of *T*, such that all four possible combinations of spatial and electromechanical concordance and discordance occur for some values of *λ* and *T*.  As *ν* decreases further to 0.4 ([Fig.  7C](#f7-cmc-suppl.1-2016-001){ref-type="fig"}) and 0.35 ([Fig.  7D](#f7-cmc-suppl.1-2016-001){ref-type="fig"}), the parameter regimes for SC and EMC alternans become smaller, such that for *ν* = 0.3 ([Fig.  7E](#f7-cmc-suppl.1-2016-001){ref-type="fig"}) and *ν* = 0.2 ([Fig.  7F](#f7-cmc-suppl.1-2016-001){ref-type="fig"}), alternans is SD and EMD for all values of *T* at which alternans is present.  Interestingly, we note that for *ν* = 0.3, the single cell model predicts a parameter regime for EMC alternans, while in the cable model, alternans is EMD for all parameter values considered, demonstrating that spatial coupling promoted EMD. 

We next measure how both SR Ca parameter *ν* and *λ* alter the stimulus period for alternans onset (*T*~Alt~, defined as APD alternans greater than 10 ms) and conduction block (*T*~CB~).  We find that both *T*~Alt~ and *T*~CB~ decrease as *ν* decreases, and in general, *T*~Alt~ and *T*~CB~ decrease as *λ* decreases, although this dependence is weaker ([Fig.  8A and B](#f8-cmc-suppl.1-2016-001){ref-type="fig"}).  Thus, reduced SR Ca uptake and release decrease both the onset of alternans and conduction block, which would appear to be antiarrhythmogenic.  Since both the onset of alternans and conduction block are functions of *ν* and *λ*, we can define an *alternans window*, ie, the size of the stimulus period range for which alternans is induced, given by *T~A~*~lt~ -- *T*~CB~.[@b40-cmc-suppl.1-2016-001] We find that the alternans window is largest for larger *λ* and *ν* ≈ 0.55 ([Fig.  8C](#f8-cmc-suppl.1-2016-001){ref-type="fig"}), while the alternans window is smaller for reduced SR Ca uptake and release. 

To further measure propensity for arrhythmias, we can divide the alternans window into windows for spatial and electromechanical concordance and discordance alternans.  We measure the alternans windows for the four possible cases for spatial and electromechanical concordance and discordance ([Fig.  9](#f9-cmc-suppl.1-2016-001){ref-type="fig"}).  The window for SC-EMC alternans is largest for larger *λ* and *ν* between 0.4 and 0.65 ([Fig.  9A](#f9-cmc-suppl.1-2016-001){ref-type="fig"}), while the window for SC-EMD alternans is largest for smaller *λ* less than 0.5 and larger *ν* between 0.55 and 0.65 ([Fig.  9B](#f9-cmc-suppl.1-2016-001){ref-type="fig"}).  The parameter regime for the SD-EMC alternans window ([Fig.  9C](#f9-cmc-suppl.1-2016-001){ref-type="fig"}) is similar to that in panel A (ie, larger *λ* and *ν* between 0.4 and 0.65), however with the window largest for intermediate values of *λ* and *ν*.  Critically, the SD-EMD alternans window is largest for smaller *ν* and *λ* ([Fig.  9D](#f9-cmc-suppl.1-2016-001){ref-type="fig"}). 

Next, we investigate to what extent SR Ca handling properties influence the size of spatial and electromechanical concordant and discordant alternans windows.  We first consider the influence of a steeper SR Ca release vs. SR Ca load relationship, by decreasing parameter *β* ([Eq.  A3](#fd25-cmc-suppl.1-2016-001){ref-type="disp-formula"}) by a factor of 2 (dashed line in [Figure A1](#fA1-cmc-suppl.1-2016-001){ref-type="fig"}, *g*(*l~n~*)).  In [Figure 10A](#f10-cmc-suppl.1-2016-001){ref-type="fig"}, we measure the alternans windows for different values of *λ* and *ν*, for reduced *β*, and we find that, in general, the parameter spaces are very similar, in comparison with the baseline value for *β*.  The SC-EMC alternans window is largest for larger *λ* and *ν*, and SD-EMD alternans windows is largest for reduced *λ* and *ν*.  Consistent with prior studies,[@b9-cmc-suppl.1-2016-001],[@b41-cmc-suppl.1-2016-001] a steeper SR Ca release relationship promoted alternans and led to a much larger SC-EMC alternans window.  Further, via the mechanism described in the "Single cell excitation--contraction coupling map model" section, steeper SR Ca release also promoted larger EMD alternans windows, for some values of *λ* and *ν*. 

We next consider the influence of a steeper SR Ca release vs. peak intracellular Ca relationship, by decreasing parameter *θ* ([Eq.  A4](#fd26-cmc-suppl.1-2016-001){ref-type="disp-formula"}) by a factor of 2 (dashed line in [Figure A1](#fA1-cmc-suppl.1-2016-001){ref-type="fig"}, $u\left( c_{n + 1}^{p} \right)$).  We find that a steeper SR Ca uptake relationship has a similar effect as that as a steeper SR Ca release relationship, in that alternans are promoted, most prominently increasing the size of the SC-EMC alternans window and to a lesser extent, promoting larger EMD alternans windows ([Fig.  10B](#f10-cmc-suppl.1-2016-001){ref-type="fig"}).  In general, in the setting of steeper SR Ca uptake, the dependence on SR Ca release and uptake parameters, *λ* and *ν*, respectively, is similar, compared with the baseline value for *θ*. 

Finally, we investigate to what extent key tissue properties influence the size of spatial and electromechanical concordant and discordant alter nans windows.  We first consider the influence of reduced cell coupling on alternans formation.  Echebarria and Karma previously derived and showed that APD coupling parameters, *w* and *ξ*, are proportional to *D* and $\sqrt{D}$, respectively, where *D* is the coefficient parameter governing the diffusion of the membrane potential in the standard cable equation.[@b35-cmc-suppl.1-2016-001] In [Figure 11A](#f11-cmc-suppl.1-2016-001){ref-type="fig"}, we measure the alternans windows, as above, after reducing both *w* and *ξ* by a factor of 2 and $\sqrt{2}$, respectively.  Interestingly, reducing the cell coupling parameters has a minimal influence on the alternans windows parameter space, with approximately 11 ms being the largest difference between the alternans windows for the baseline cell coupling values in [Figure 9](#f9-cmc-suppl.1-2016-001){ref-type="fig"} and reduced cell coupling in [Figure 11A](#f11-cmc-suppl.1-2016-001){ref-type="fig"}.  In general, the SC alternans windows are slightly smaller and the SD alternans windows are slightly larger, for reducing cell coupling.  However, we note that reducing *w* and *ξ* only accounts for the effects of reduced cell coupling on APD.  The effects of reducing cell coupling on CV and CV restitution may further alter the size of alternans windows and will be investigated in further studies. 

We finally investigate the influence of reducing the size of the cardiac tissue, by reducing the cable length *L* by a factor of 2 ([Fig.  11B](#f11-cmc-suppl.1-2016-001){ref-type="fig"}).  Previous studies have shown that tissue of a sufficiently large size is needed for the formation of SD alternans, and indeed, longer cable lengths promote the formation of more nodes or phase reversals.[@b5-cmc-suppl.1-2016-001] We similarly find that reducing the cable length reduces the size of SD alternans windows and increases the size of SC alternans windows, with alternans window differences between the (*ν, λ*)--parameter spaces as large as 40 ms for some specific parameter sets.  Interestingly, EMC and EMD alternans windows are not greatly influenced by tissue size, as the general shape of the parameter space maps is similar, with the exception of larger SC-EMC alternans windows for small values of *ν*.  Thus, in general, reducing tissue size reduces the formation of SD alternans windows, consistent with previous work,[@b5-cmc-suppl.1-2016-001] with additional modulation by SR Ca uptake and release. 

Discussion
==========

Summary of main findings
------------------------

In this study, we incorporate a detailed single cell map model of Ca cycling and bidirectional APD--Ca coupling into a spatially extended tissue model to investigate the influence of SR Ca uptake and release properties on the incidence and characteristics of alternans and conduction block.  We find that an intermediate SR Ca uptake rate and larger SR Ca release resulted in the widest range of stimulus periods that promoted alternans, measured by the alternans window ([Fig.  8C](#f8-cmc-suppl.1-2016-001){ref-type="fig"}).  However, reduced SR Ca uptake and release lead to the largest window for arrhythmogenic SD and EMD alternans ([Fig.  9D](#f9-cmc-suppl.1-2016-001){ref-type="fig"}). 

Our findings illustrate the complex interaction between SR Ca handling, APD and CV restitution, and cell coupling, influencing the formation of alternans and electromechanical and spatial discordance at the cellular and tissue level, respectively.  As illustrated in [Figure 5](#f5-cmc-suppl.1-2016-001){ref-type="fig"}, subtle changes in SR Ca handling can alter APD via Ca--APD coupling, which in turn alters propagation via CV restitution, and can lead to spatial discordance.  In particular, our results suggest the interesting prediction that an intervention that acts to inhibit alternans onset and conduction block is, in fact, not less arrhythmogenic.  By facilitating conduction at faster rhythms, ie, shorter *T*, reduced SR Ca uptake and release engage dynamical instabilities at the cellular level, promoting electromechanical discordance, and tissue level, promoting spatial discordance, that are not present for normal SR Ca handling.  Critically, we find that normal SR Ca uptake and release play pivotal roles in preventing arrhythmias. 

Further, we performed additional simulations in both the setting of a steeper SR Ca release vs. SR Ca load relationship and a steeper SR Ca uptake vs. peak intracellular Ca relationship, in order to highlight the generality of our findings ([Fig.  10](#f10-cmc-suppl.1-2016-001){ref-type="fig"}).  Steepening both of these relationships led to larger alternans windows, in particular the window for SC-EMC alternans, and to a lesser extent, larger EMD alternans.  However, importantly, we show that the general conclusions of the study did not depend on the steepness of these relationships: reduced SR Ca uptake and release promote SD and EMD alternans. 

Our findings demonstrate that reduced cellular coupling during repolarization has minimal influence on EMD ([Fig.  11A](#f11-cmc-suppl.1-2016-001){ref-type="fig"}), suggesting that EMD is primarily governed by SR Ca handling properties.  However, differences between single cell and cable map model results emphasize that tissue-level properties do play a role in governing EMD, suggesting that CV and CV restitution are also critical factors.  Future work will investigate the role of CV restitution on electromechanical discordance, in conjunction with altered SR Ca handling.  Finally, in agreement with prior work,[@b5-cmc-suppl.1-2016-001] we find that tissue size is also a critical factor governing the formation of SD alternans, while demonstrating the novel result that SR Ca handling properties modulate this relationship ([Fig.  11B](#f11-cmc-suppl.1-2016-001){ref-type="fig"}). 

Comparison with prior results
-----------------------------

Many prior studies have utilized one-dimensional cable-coupled map models to investigate the arrhythmogenicity of different cellular- and tissue-level properties.  Fox et al previously showed that spatiotemporal heterogeneity in electrical properties promoted large APD alternans, potentially leading to conduction block.[@b26-cmc-suppl.1-2016-001] Conduction block was less likely in the setting of a shallow APD restitution curve, slow conduction at short DIs, and when short-term memory provided a larger contribution to APD dynamics.[@b36-cmc-suppl.1-2016-001] Cell--cell coupling, or electrotonic interactions, during repolarization have been shown to suppress alternans in the setting of steep APD restitution[@b42-cmc-suppl.1-2016-001] and the onset of SD alternans[@b39-cmc-suppl.1-2016-001] and, similarly, suppress incidences of conduction block.[@b36-cmc-suppl.1-2016-001] Henry and Rappel, confirming behavior observed in an ionic model, demonstrated a relationship between tissue size and the dynamics of conduction block.[@b37-cmc-suppl.1-2016-001] The aforementioned studies investigated alternans and conduction block properties in models accounting for APD and CV restitution, cell--cell coupling, and in some studies, short-term memory, but did not explicitly account for cellular Ca dynamics. 

In an excellent 2005 study, Shiferaw et al investigated the coupled dynamics of the transmembrane potential and Ca cycling in both a detailed ionic model and reduced single cell coupled map, demonstrating instabilities leading to EMC and EMD alternans and quasiperiodic oscillation.[@b11-cmc-suppl.1-2016-001] A further reduction investigating Ca-driven alternans identified the slope of SR Ca release vs. SR Ca load and the efficiency of SR Ca uptake as two key parameters governing instabilities.[@b9-cmc-suppl.1-2016-001],[@b41-cmc-suppl.1-2016-001] Shiferaw and Karma showed that in a single cell model coupling APD and Ca, with subcellular diffusion between neighboring sarcomeres, a Turing-type symmetry breaking instability driven by Ca diffusion can lead to out-of-phase Ca alternans within a single cell.[@b43-cmc-suppl.1-2016-001] In the single cell map model that is used in this study,[@b25-cmc-suppl.1-2016-001] Qu et al showed that SR Ca uptake rate *ν* influenced the parameter regime for alternans and complex dynamics, particularly for *ν* in the range of 0.1--0.6, consistent with our findings, although this study was only performed in the setting of the Ca cycling subsystem for fixed total cellular Ca.  In the map model with fully coupled APD and Ca dynamics, the authors demonstrated that negative Ca--APD coupling could lead to EMC or EMD alternans, as well as quasiperiodicity and chaos. 

Recently, Skardal et al formulated a spatially extended coupled map model, incorporating APD and intracellular Ca coupling, demonstrating that for Ca-driven instabilities, large discontinuous jumps in the Ca alternans amplitude can be observed, while APD alternans vary smoothly, a phenomenon not observed for voltage-driven instabilities.[@b44-cmc-suppl.1-2016-001],[@b45-cmc-suppl.1-2016-001] While providing significant insight into alternans incidence in the setting of Ca-driven instabilities in tissue, the authors used phenomenological forms of the maps governing APD and Ca dynamics.  To our knowledge, our study is the first to incorporate a detailed representation of intracellular and SR Ca cycling[@b25-cmc-suppl.1-2016-001] into a spatially extended coupled map model. 

Physiological implications
--------------------------

In this study, we investigated the influence of impaired SR Ca uptake and release on cellular- and tissue-level alternans.  Impaired SR Ca uptake prevents complete refilling of the SR stores, in particular at fast heart rates, promoting alternans via the mechanism described in the "Single cell excitation--contraction coupling map model" section.[@b46-cmc-suppl.1-2016-001],[@b47-cmc-suppl.1-2016-001] Many studies have demonstrated that SR Ca uptake is impaired in heart failure, a consequence of reduced expression of the SR Ca-ATPase (SERCA) pump.[@b28-cmc-suppl.1-2016-001],[@b31-cmc-suppl.1-2016-001],[@b32-cmc-suppl.1-2016-001],[@b48-cmc-suppl.1-2016-001],[@b49-cmc-suppl.1-2016-001] Reduced SR Ca uptake in heart failure sub sequently results in reduced SR Ca load, consistent with our findings shown in [Figure 2](#f2-cmc-suppl.1-2016-001){ref-type="fig"} (bottom panels), a reduction that is also promoted by enhanced diastolic SR Ca leak,[@b30-cmc-suppl.1-2016-001],[@b50-cmc-suppl.1-2016-001] and subsequently, results in alternans onset at slower heart rates and higher arrhythmia incidence.[@b47-cmc-suppl.1-2016-001] Further, upregulation of the SR Ca pump, SERCA2a, via gene transfer into failing hearts has been shown to suppress the incidence of alternans and critically reduce incidence of cardiac arrhythmias,[@b51-cmc-suppl.1-2016-001],[@b52-cmc-suppl.1-2016-001] consistent with our simulation predictions as well. 

Experimental evidence is less clear as to whether SR Ca release is impaired during heart failure.  Protein kinase A hyperphosphorylation of ryanodine receptors during heart failure has been associated with increased open probability of the ryanodine receptors, leading to the aforementioned diastolic SR Ca leak.[@b32-cmc-suppl.1-2016-001] However, there is also evidence that ryanodine receptor expression is reduced in human ventricular myocytes from heart failure patients.[@b29-cmc-suppl.1-2016-001],[@b33-cmc-suppl.1-2016-001],[@b53-cmc-suppl.1-2016-001] Further, there is significant evidence that T-tubule structure is disorganized in heart failure, leading to so-called *orphaned* ryanodine receptors and slower and less synchronous Ca release throughout the cell.[@b33-cmc-suppl.1-2016-001],[@b54-cmc-suppl.1-2016-001] Thus, in the context of our present study, in which SR Ca release parameter *λ* governs systolic SR Ca release, we investigated the influence of impaired systolic SR Ca release, ie, *λ* \< 1. 

Overall, our simulations of impaired SR Ca release and uptake are consistent with the clinical presentation of T-wave alternans in heart failure patients.  T-wave alternans, a precursor for lethal cardiac arrhythmias, is the clinical manifestation of cellular-level APD alternans and is associated with many cardiac diseases, including heart failure.[@b55-cmc-suppl.1-2016-001]--[@b57-cmc-suppl.1-2016-001] Consistent with these findings, we previously demonstrated in an *in vivo* tissue model that a larger SD alternans window was associated with a greater incidence of spiral wave formation.[@b40-cmc-suppl.1-2016-001] Due to the simplicity of the coupled map model, it is not possible to directly correlate the heart failure disease state with the SR Ca uptake and release parameters that are the focus of this study.  However, our results highlight the complex and highly nonlinear relationship between impaired SR Ca handling and the manifestation of arrhythmogenic EMD and SD alternans, consistent with the highly unpredictable clinical presentation of arrhythmias, in particular in heart failure patients.  The unpredictable nature of arrhythmias in heart failure patients is further complicated by the fact that, as noted by Bers et al, heart failure as a disease is heterogeneous, complex, and the relative contributions of impaired SR Ca uptake and leak likely vary in different origins and stages of the disease.[@b58-cmc-suppl.1-2016-001] Thus, the overall qualitative conclusions of our results, ie, impaired SR Ca uptake and release promote arrhythmogenic alternans, agree with previous findings that these properties are critical factors in promoting arrhythmias in heart failure. 

Limitations
-----------

The coupled map model formulation is inherently a reduction of complex subcellular-, cellular-, and tissue-level dynamics.  Our model is an attempt to incorporate detailed Ca cycling into a tissue model, while still maintaining the computational efficiency of a reduced coupled map.  However, the model is not without limitations.  In particular, our model does not account for the detailed kinetics of the transmembrane potential and ion channel gating[@b59-cmc-suppl.1-2016-001],[@b60-cmc-suppl.1-2016-001] and neglects stochastic aspects of subcellular Ca release, such as Ca sparks and microdomain Ca fluctuations.[@b61-cmc-suppl.1-2016-001]--[@b64-cmc-suppl.1-2016-001] Stochasticity may be particularly important to account for in future work, as Sato et al showed that Ca fluctuations can govern the phase of Ca alternans, which in turn may lead to the formation of SD alternans.[@b8-cmc-suppl.1-2016-001] Ma and Xiao showed that Ca fluctuations can alter alternans formation and lead to chaotic or disordered dynamics in a single cell-coupled map.[@b65-cmc-suppl.1-2016-001] It is worth noting that the use of the discrete time-coupled map model provides a significant computational savings over these more detailed continuous-time ordinary and stochastic differential equation formulations.  The present study performed considered over 825,000 parameter combinations (varying *λ*, *v, β, θ, L, w, ξ*, and *T*), summarized in [Figures 8](#f8-cmc-suppl.1-2016-001){ref-type="fig"}[](#f9-cmc-suppl.1-2016-001){ref-type="fig"}[](#f10-cmc-suppl.1-2016-001){ref-type="fig"}--[11](#f11-cmc-suppl.1-2016-001){ref-type="fig"}, a parameter investigation that would be essentially computationally infeasible using the detailed models mentioned above.  Importantly, the coupled map model simulations highlight potentially significant parameter regimes that can be investigated further in more detailed model formulations. 

Further, our model is limited to a one-dimensional tissue and does not account for the complex three-dimensional architecture of cardiac tissue, in particular, heterogeneities in anatomical and functional properties.[@b60-cmc-suppl.1-2016-001] Finally, our study focuses on two critical SR Ca handling properties, and thus, the conclusions of our study may depend on other model parameters, such as APD and CV restitution and cell coupling properties, that may vary in settings of disease.  However, these limitations do not detract from the overall conclusion that impaired SR Ca release and uptake can lead to an arrhythmogenic state, promoting instabilities at the cellular and tissue level. 

**ACADEMIC EDITOR:** Thomas E. Vanhecke, Editor in Chief

**PEER REVIEW:** Five peer reviewers contributed to the peer review report.  Reviewers' reports totaled 2,389 words, excluding any confidential comments to the academic editor. 

**FUNDING:** Author discloses no external funding sources. 

**COMPETING INTERESTS:** Author discloses no potential conflicts of interest. 

Paper subject to independent expert blind peer review.  All editorial decisions made by independent academic editor.  Upon submission manuscript was subject to anti-plagiarism scanning.  Prior to publication all authors have given signed confirmation of agreement to article publication and compliance with all applicable ethical and legal requirements, including the accuracy of author and contributor information, disclosure of competing interests and funding sources, compliance with ethical requirements relating to human and animal study participants, and compliance with any copyright requirements of third parties.  This journal is a member of the Committee on Publication Ethics (COPE). 

**Author Contributions**

Conceived and designed the experiments: SHW.  Analyzed the data: SHW.  Wrote the first draft of the manuscript: SHW.  Developed the structure and arguments for the paper: SHW.  Made critical revisions: SHW.  The author reviewed and approved of the final manuscript. 

The author is thankful for the support of the Turing High Performance Computing Cluster at Old Dominion University. 

Model equations
===============

In this appendix, the full equations for APD and CV restitution and Ca cycling are given, as described in the "Methods" section, originally from the study by Qu et al.[@b25-cmc-suppl.1-2016-001] (unless otherwise noted), and are plotted in [Figure A1](#fA1-cmc-suppl.1-2016-001){ref-type="fig"}.  All model parameters are given in [Table A1](#tA1-cmc-suppl.1-2016-001){ref-type="table"}.  The APD restitution function is given by $$f(d_{n}) = \begin{cases}
{a_{0}\left\lbrack {1 - \left\{ {1 + \exp((d_{n} - d_{0})/\tau_{0})} \right\}^{- 1}} \right\rbrack} & {\text{if}\ d_{n} \geq d_{\min}} \\
{0,} & \text{otherwise} \\
\end{cases}$$

The SR Ca release function that describes the restitution properties is given by $$q(d_{n}) = 1 - \sigma exp( - d_{n}/\tau_{q}),$$and the function that describes the SR Ca load dependence is given by $$g(l_{n}) = l_{n}\left\lbrack {1 - \frac{1 - \alpha}{1 + \exp((l_{n} - l_{c})/\beta)}} \right\rbrack.$$

The SR Ca uptake function that describes the dependence on the stimulus period is given by $$u(t_{n + 1}) = 1 - \rho\exp( - t_{n + 1}/\tau_{u}),$$and the function that describes the dependence on peak intracellular Ca is given by $$h(c_{n + 1}^{p}) = c_{n + 1}^{p}\left\lbrack {1 - \frac{1}{1 + \exp((c_{n + 1}^{p} - c_{0})/\theta)}} \right\rbrack.$$

The steady-state total Ca concentration is given by $$c(t_{n + 1}) = c_{0}\lbrack 1 + \varepsilon\exp( - t_{n + 1}/\tau_{c})\rbrack.$$

The CV restitution curve, from the study by Fox et al.[@b36-cmc-suppl.1-2016-001], is given by $$v(d_{n}) = v_{\max}\lbrack 1 - \exp( - (d_{n} + d_{\upsilon})/\delta)\rbrack.$$

[Eqs.  8a](#fd8-cmc-suppl.1-2016-001){ref-type="disp-formula"}--[8d](#fd11-cmc-suppl.1-2016-001){ref-type="disp-formula"} are given by $$\begin{array}{l}
{F_{a}(a_{n},c_{n},l_{n},t_{n};\lambda) = \frac{f(t_{n} - a_{n})}{1 - \gamma\lbrack c_{n} + \lambda_{q}(t_{n} - a_{n})g(l_{n})\rbrack},} \\
{F_{b}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda) = b_{n} - \kappa\lbrack c_{n} - c(t_{n + 1})\rbrack + \eta\lbrack F_{a}(a_{n},c_{n},l_{n},t_{n};\lambda) - a_{n}\rbrack,} \\
{F_{c}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda,\nu) = F_{b}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda) - F_{l}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda,\nu),} \\
{F_{l}(a_{n},c_{n},l_{n},t_{n},t_{n + 1};\lambda,\nu) = l_{n} - \lambda q(t_{n} - a_{n})g(l_{n}) + \nu u(t_{n + 1})h\lbrack c_{n} + \lambda q(t_{n} - a_{n})g(l_{n})\rbrack.} \\ 
\end{array}$$

###### 

Model parameters. 

  PARAMETER   DEFINITION                   UNITS      VALUE
  ----------- ---------------------------- ---------- -----------------------------------------------------------
  *a*~0~      APD restitution parameter    ms         220
  *d*~0~      APD restitution parameter    ms         40
  **τ**~0~    APD restitution parameter    ms         30
  *d~min~*    APD restitution parameter    ms         2
  *v~max~*    CV restitution parameter     mm/ms      0.55
  *d~v~*      CV restitution parameter     ms         2
  δ           CV restitution parameter     ms         14
  *σ*         SR Ca release parameter      --         0.5
  *t~q~*      SR Ca release parameter      ms         80
  *a*         SR Ca release parameter      --         0.036
  *l~c~*      SR Ca release parameter      µM         93.5
  *β*         SR Ca release parameter      µM         5[\*](#tfn2-cmc-suppl.1-2016-001){ref-type="table-fn"}
  λ           SR Ca release parameter      --         Varied
  *ν*         SR Ca uptake parameter       --         Varied
  *ρ*         SR Ca uptake parameter       --         0.15
  τ*~u~*      SR Ca uptake parameter       ms         200
  c~0~        SR Ca uptake parameter       *µ*M       28
  θ           SR Ca uptake parameter       *µ*M       20[\*](#tfn2-cmc-suppl.1-2016-001){ref-type="table-fn"}
  γ           Ca--APD coupling parameter   *µ*M^−1^   10^−3^
  *e*         Ca accumulation parameter    --         2
  **τ**~c~    Ca accumulation parameter    ms         300
  *k*         Ca accumulation parameter    --         0.1
  η           APD--Ca coupling parameter   *µ*M/ms    0.1
  ξ           Cell coupling parameter      mm         1[\*](#tfn2-cmc-suppl.1-2016-001){ref-type="table-fn"}
  *w*         Cell coupling parameter      mm         0.35[\*](#tfn2-cmc-suppl.1-2016-001){ref-type="table-fn"}
  *L*         Cable length                 mm         100[\*](#tfn2-cmc-suppl.1-2016-001){ref-type="table-fn"}

**Notes:** All parameters pertaining to APD and Ca dynamics are from the study by Qu et al.[@b25-cmc-suppl.1-2016-001], and CV restitution and cell coupling parameters are from the study by Fox et al.[@b36-cmc-suppl.1-2016-001]

Value used in simulations, unless otherwise noted. 

![ Nonlinear map model functions.  The model equations, given by [Eqs.  A1](#fd23-cmc-suppl.1-2016-001){ref-type="disp-formula"}--[A7](#fd29-cmc-suppl.1-2016-001){ref-type="disp-formula"}, are shown, in addition to the APD restitution slope.  The dashed line in the panels showing *g*(*l~n~*) and $h(c_{n + 1}^{p})$ correspond with decreasing the scaling parameters, *β* and *θ*, respectively, by a factor of 2.](cmc-suppl.1-2016-001fA1){#fA1-cmc-suppl.1-2016-001}

![ Single cell and one-dimensional cable simulations. (** A**, top) Illustration of the APD ***a****~n~*, DI *d~n~*, and IBI *t~n~*; (middle) peak intracellular Ca $(c_{n}^{p})$ at the end of each beat (*c~n~*), SR Ca release (*r~n~*); (bottom) and SR Ca load at the end of each beat (*l~n~*). (** B**) Illustration of action potential wavefronts (solid) and wavebacks (dashed), propagating along a cable, denoting APD, DI, and IBI as ***a****~n~*(*x*), *d~n~*(*x*), and *t~n~*(*x*), respectively, as functions of space *x*.](cmc-suppl.1-2016-001f1){#f1-cmc-suppl.1-2016-001}

![ Reduced SR Ca uptake can promote EMD alternans.  APD, peak intracellular Ca, and SR Ca load are shown as a function of stimulus period *T* on even (black) and odd (red) beats, for SR Ca uptake rate *v* = (**A**) 0.7, (**B**) 0.4, (**C**) 0.28, and (**D**) 0.1.  Other parameters: *λ* = 1.](cmc-suppl.1-2016-001f2){#f2-cmc-suppl.1-2016-001}

![ Mechanism underlying reduced SR Ca uptake-driven EMD alternans.  Map model variables, *a~n~*, *c^p^~n~*, *r~n~*, *c~n~*, *u~n~*, and *l~n~*, are shown as a function of beat number, for SR Ca uptake rate *v* = (**A**) 0.7 and (**B**) 0.1.  Other parameters: *T* = 270 ms, *λ* = 1.](cmc-suppl.1-2016-001f3){#f3-cmc-suppl.1-2016-001}

![ Reduced SR Ca uptake and release can promote EMD alternans.  Parameter regimes for no alternans, EMC alternans, EMD alternans, and loss of capture are shown as a function of stimulus period *T* and SR Ca release parameter *λ*, for SR Ca uptake rate *v* = (**A**) 0.7, (**B**) 0.4, (**C**) 0.28, and (**D**) 0.1.](cmc-suppl.1-2016-001f4){#f4-cmc-suppl.1-2016-001}

![ Reduced SR Ca uptake promotes complex spatial and electromechanical alternans patterns.  APD and peak intracellular Ca are shown as a function of spatial location on even (black) and odd (red) beats, for SR Ca uptake rate *v* = (**A**) 0.7, (**B**) 0.4, (**C**) 0.35, and (**D**) 0.3.  Other parameters: *λ* = 1, *T* = 288 ms.](cmc-suppl.1-2016-001f5){#f5-cmc-suppl.1-2016-001}

![ Reduced SR Ca release promotes complex spatial and electromechanical alternans patterns.  APD and peak intracellular Ca are shown as a function of spatial location on even (black) and odd (red) beats, for SR Ca release parameter *λ* = (**A**) 1, (**B**) 0.5, (**C**) 0.36, and (**D**) 0.3.  Other parameters: *v* = 0.4, *T* = 280 ms.](cmc-suppl.1-2016-001f6){#f6-cmc-suppl.1-2016-001}

![ Reduced SR Ca uptake and release promote complex alternans patterns.  Parameter regimes for no alternans, SC, SD, EMC, and EMD alternans, and conduction block are shown as a function of stimulus period *T* and SR Ca release parameter *λ*, for SR Ca uptake rate *v* = (**A**) 0.7, (**B**) 0.5, (**C**) 0.4, (**D**) 0.35, (**E**) 0.3, and (**F**) 0.2.  The regime boundary between SC and SD alternans regimes is indicated by a thick black line, and between EMC and EMD alternans, regime is indicated by a thin black line.  The parameter regimes for cable simulations are bounded by black lines.  Parameter regimes are identified as SC-EMC (blue), SD-EMC (red), SC-EMD (yellow), and SD-EMD (green).  The boundaries for the single cell simulation regimes are shown in gray, for comparison.](cmc-suppl.1-2016-001f7){#f7-cmc-suppl.1-2016-001}

![ Reduced SR Ca uptake and release alter alternans onset and conduction block.  The stimulus period for alternans onset (**A**, *T*~Alt~) and conduction block (**B**, *T*~CB~) is shown as a function of SR Ca uptake rate *v* and release parameter *λ*.  The alternans window (*T*~Alt~ − *T*~CB~) is shown in **C**.](cmc-suppl.1-2016-001f8){#f8-cmc-suppl.1-2016-001}

![ SC, EMC, SD, and EMD alternans windows.  The alternans windows for (**A**) SC-EMC, (**B**) SC-EMD, (**C**) SD-EMC, and (**D**) SD-EMD are shown as a function of SR Ca uptake rate *v* and SR Ca release parameter *λ*.](cmc-suppl.1-2016-001f9){#f9-cmc-suppl.1-2016-001}

![ SC, EMC, SD, and EMD alternans windows for steeper SR Ca release and uptake.  The alternans windows for SC-EMC, SC-EMD, SD-EMC, and SD-EMD are shown as a function of SR Ca uptake rate *v* and SR Ca release parameter *λ*, for (**A**) a steeper SR Ca release vs. SR Ca load relationship and (**B**) a steeper Ca uptake vs. peak intracellular Ca relationship.](cmc-suppl.1-2016-001f10){#f10-cmc-suppl.1-2016-001}

![ SC, EMC, SD, and EMD alternans windows for reduced cell coupling and tissue size.  The alternans windows for SC-EMC, SC-EMD, SD-EMC, and SD-EMD are shown as a function of SR Ca uptake rate *v* and SR Ca release parameter *λ*, for (**A**) reduced cell coupling and (**B**) reduced tissue size.](cmc-suppl.1-2016-001f11){#f11-cmc-suppl.1-2016-001}
