Figures
Abstract
Early afterdepolarizations (EADs) are abnormal depolarizations during the plateau phase of the action potential, which are known to be associated with lethal arrhythmias in the heart. There are two major hypotheses for EAD genesis based on experimental observations, i.e., the voltage (Vm)-driven and intracellular calcium (Ca)-driven mechanisms. In ventricular myocytes, Ca and Vm are bidirectionally coupled, which can affect each other’s dynamics and result in new dynamics, however, the roles of Ca cycling and its coupling with Vm in the genesis of EADs have not been well understood. In this study, we use an action potential model that is capable of independent Vm and Ca oscillations to investigate the roles of Vm and Ca coupling in EAD genesis. Four different mechanisms of EADs are identified, which are either driven by Vm oscillations or Ca oscillations alone, or oscillations caused by their interactions. We also use 5 other ventricular action potential models to assess these EAD mechanisms and show that EADs in these models are mainly Vm-driven. These mechanistic insights from our simulations provide a theoretical base for understanding experimentally observed EADs and EAD-related arrhythmogenesis.
Author summary
Early afterdepolarizations (EADs) are dangerous abnormal electrical activities in the heart, which may cause lethal arrhythmias. Although EADs have been widely investigated in both experimental and computer simulation studies, their mechanisms remain incompletely understood. In the present work, we carry out computer simulations using action potential models with detailed formulations of ionic currents and calcium cycling to investigate the genesis of EADs. Different mechanisms and causes are identified in our simulations, which agree with experimental observations. The mechanistic insights from our work provide a theoretical base for understanding the mechanisms of EADs and EAD-related arrhythmogenesis.
Citation: Wang R, Qu Z, Huang X (2024) Dissecting the roles of calcium cycling and its coupling with voltage in the genesis of early afterdepolarizations in cardiac myocyte models. PLoS Comput Biol 20(2): e1011930. https://doi.org/10.1371/journal.pcbi.1011930
Editor: Jeffrey J. Saucerman, University of Virginia, UNITED STATES
Received: November 4, 2023; Accepted: February 19, 2024; Published: February 28, 2024
Copyright: © 2024 Wang et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All relevant data are within the manuscript and its Supporting information files.
Funding: This work is supported by the Guangdong Basic and Applied Basic Research Foundation under Grant No. 2021A1515010500 (XH); National Institutes of Health grants R01 HL134709, R01 HL139829, R01 HL157116, and P01 HL164311 (ZQ). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Early afterdepolarizations (EADs) are secondary transmembrane voltage (Vm) depolarizations of a cardiac myocyte during the plateau and repolarizing phases [1, 2]. EADs are associated with lethal arrhythmias under certain diseased conditions, such as Torsade de Pointes in long QT syndrome (LQTS) [3–7]. Due to its importance in cardiac arrhythmogenesis, the causes and mechanisms of EADs have been comprehensively investigated in experiments [8–13], computer simulations [14–21], and bifurcation theories [22–27].
Hitherto there are two dominant hypotheses of EAD generation mechanisms [2, 28, 29]. The first one is the voltage-driven (Vm-driven) mechanism in which EADs are secondary depolarizations driven by instabilities in the voltage system. In this mechanism, the induction of EADs is attributed to two major factors, i.e., reduced repolarization reserve (RRR) and reactivation of L-type calcium (Ca) current (ICa,L) [5, 9]. RRR is a condition in which the outward currents are reduced and/or the inward currents are increased. For example, in LQTS type 1 and 2, potassium currents are reduced, while in LQTS type 3, the late component of the sodium current is increased. RRR decelerates repolarization and thereby provides enough time for reactivation of ICa,L to result in depolarizations in the repolarization phase. Based on nonlinear dynamics analyses, RRR combining proper ICa,L kinetics can give rise to a Hopf bifurcation to lead to limit cycle oscillations in the plateau or repolarizing phase of the action potential (AP), manifesting as EADs [22, 24].
The second one is the Ca-driven mechanism in which instabilities in the intracellular Ca cycling system give rise to Ca oscillations or spontaneous Ca releases to generate EADs. Evidence of this mechanism has been demonstrated in several experimental studies [29–32]. The physiological hypothesis for this mechanism is as follows: Spontaneous Ca releases elevate the intracellular Ca level which increases the sodium-calcium-exchange (NCX) current (INCX) to depolarize Vm, manifesting as EADs. In this case, Ca oscillations are the primary cause and Vm oscillations are a passive result of Ca oscillations. Simulation studies (e.g., Ref. [33–35]) have also been carried out to investigate this mechanism.
Although rigorous bifurcation analyses [20, 22–27] have revealed the Vm-driven instabilities induced bifurcations leading to EADs, the dynamical mechanisms and bifurcations for Ca-driven EADs have not been well understood. Moreover, Vm and Ca are coupled via Ca-dependent ionic currents, and the coupling can either potentiate or suppress each other’s instabilities, or bring in completely new dynamics [36, 37]. Therefore, a more rigorous analysis of Ca-driven instabilities and the roles of the interaction of Ca and Vm in the genesis of EADs are necessary. To investigate the roles of Ca-driven instabilities in EAD genesis, one needs an AP model capable of exhibiting spontaneous Ca oscillations. Most of the AP models do not exhibit this behavior, but a recent model developed by Wilson et al [34] (called WG model in the present paper) incorporates a Ca cycling model that exhibits spontaneous oscillations and was used to investigate the Ca-driven mechanism of EADs. In this study, we use the WG model to perform systematic analyses on the roles of Ca and Vm coupling in the genesis of EADs. Our analysis reveals 4 mechanisms of EADs: 1) EADs driven by Vm oscillations alone; 2) EADs driven by Ca oscillations alone; 3) EADs driven by a Ca-ICa,L-Vm-Ca feedback loop; and 4) EADs driven by a large Ca transient. We also carry out simulations of 5 other ventricular AP models and observe that the majority of EADs in these models are via the Vm-driven mechanism, only a very small portion are related to the Ca-driven mechanism.
Models and methods
The WG model is a rabbit ventricular myocyte AP model modified from the Shannon and Bers (SB) model [38], which can be presented in the general formulation as follows: (1) where Vm is the transmembrane voltage, and y is a vector representing the gating variables (i.e., y = {m, h, d, f, …}, and each one of them controls the activation or inactivation of the corresponding ion channel). C represents the set of concentrations of different ions in different compartments of the cell (i.e., C = {[Ca2+]JSR, [K+]JSR, [Na+]JSR, [Ca2+]JXN, [K+]JXN, [Na+]JXN, …}, and the compartments are schematically shown in Fig 1A). Iz is the set of transmembrane ionic currents which modulates the corresponding intracellular ion cyclings (i.e., Iz = {ICaL,JXN, IKs,JXN, INCX,JXN, ICaL,SL, IKs,SL, INCX,SL, …}). P is the set of opening probabilities of the ryanodine receptor (RyR) gates controlling Ca release from junctional sarcoplasmic reticulum (JSR). Iion is the sum of the transmembrane ionic currents, and F and H stand for the functions of the corresponding variables. A schematic diagram of the WG model is shown in Fig 1A. For detailed mathematical formulations and definitions of the variables, one is referred to the original paper by Wilson et al [34].
A. A schematic plot of the WG model. The outer solid line represents the membrane. The space between the dashed line and the solid line is the submembrane space which is divided into two compartments: the junctional (JXN) and subsarcolemmal (SL) space. The JXN space (grey) is the cleft adjacent to the junctional sarcoplasmic reticulum (JSR, green), and the SL space (light grey) is the remaining submembrane space. NSR is the network sarcoplasmic reticulum. B. The AP of the control case.
In order to investigate the roles of Ca cycling and its coupling with Vm in the genesis of EADs, we modulate several major parameters associated with Ca cycling and Vm. Note that the model divides the submembrane space into two compartments, i.e., the junctional (JXN) space and subsarcolemmal (SL) space (as colored by grey and light grey in Fig 1A). The formulation and maximal conductance of each current in both compartments are identical, i.e., PCa,JXN = PCa,SL = PCa, GKs,JXN = GKs,SL = GKs and . The parameters we vary are: kmax (the maximal SR release rate constant controlling releasing flux via RyR, control kmax = 0.2 ms−1), the diffusive strength of Ca2+ between SL and cytoplasm, control JCaslmyo = 7.4485*10−13 l/ms), PCa (permeability to Ca which controls the intensity of ICa,L, control PCa=0.00027 cm/s), GKs (maximal conductance of slowly inward rectified potassium channel, control GKs = 1.23 mS/μF), and (NCX current density, control ). These parameters are denoted in Fig 1A by red, and their associated formulations are given in the original publications [34, 38]. The AP of the control model is shown in Fig 1B. An adjustable factor α(z) (z stands for any one of the adjustable parameters) is multiplied to the control value. The other parameters are set as default unless specified otherwise.
In this model, due to fast diffusion between the two submembrane compartments, [Ca2+]JXN and [Ca2+]SL keep nearly identical, and thus we just treat them as one variable as [Ca]sub (i.e., [Ca2+]JXN = [Ca2+]SL = [Ca]sub). Clamp of [Ca]sub means simultaneous clamp of [Ca2+]JXN and [Ca2+]SL. The explicit Euler method is used to integrate Eq 1 with a time step 0.01 ms. The stimulus current density Isti is a pulse of 1 ms duration and 40 pA/pF magnitude.
For simplicity, we use a single set of initial condition and apply a single stimulus at t = 0 for most of the simulations to assess the EAD behaviors. The initial condition is presented in Supplemental Information (S1 Text). We have also used other pacing protocols (i.e., S1S1 and S1S2 pacing) to assess the EAD behaviors. To dissect out the EAD mechanisms, we use voltage clamp, Ca clamp, and ionic current clamp in a case-by-case manner, as will be shown in each EAD mechanism.
We also perform simulations using 5 other AP models. The control parameters and initial conditions for each model are presented in S1 Text.
Results
EADs driven by Vm oscillations (type I)
We first set kmax to be relatively small and JCaslmyo to be relatively large, which give rise to a low [Ca]sub due to an attenuated release and a rapid diffusion between the submembrane space and cytoplasm. Fig 2A is a phase diagram showing the behaviors in the α(PCa) − α(GKs) plane. The upper left grey region represents repolarization failure, and the lower right white region is normal repolarization. Three mechanisms of EADs are identified in this phase diagram, which are colored by dark yellow, green, and magenta. In the dark yellow region, the EADs are caused by the traditional Vm-driven oscillations. Fig 2B shows a representative AP. During the EAD phase, the Vm turning point from the decreasing phase to the rising phase occurs earlier than that of [Ca]sub, an indication that Vm leads the oscillations, i.e., the oscillations are originated from instabilities in the Vm-system. Moreover, the rising of Vm is with ICa,L but not INCX, indicating that EADs are caused by reactivation of ICa,L. To further show that the oscillations are originated from Vm, we clamp [Ca]sub at different values within the oscillation range of [Ca]sub (grey region in Fig 2B) during the EAD phase to dissect out the role of Ca cycling (Fig 2C). After clamping [Ca]sub, Vm oscillations still occur although Vm may either repolarize much later or fail to repolarize. These results indicate that under this parameter setting, the EADs are driven by oscillations from Vm, which is the well-known mechanism of EADs that has been investigated widely [5, 9, 20–22, 24].
A. Phase diagram showing different EAD behaviors in the α(PCa) − α(GKs) plane for α(kmax) = 2, α(JCaslmyo) = 2.64, and . The grey and white regions represent repolarization failure and normal repolarization, respectively. Insets show two representative APs. Three types of EADs are observed, as colored dark yellow, green, and magenta. B. A representative AP and associated currents in the dark yellow region (i.e., EADs driven by Vm oscillations). The parameter set is α(GKs) = 1.16 and α(PCa) = 5.4, as marked by the black circle in A. The dashed vertical lines indicate the takeoff moments of the two EADs. The grey patch indicates the range of [Ca]sub clamped used in the [Ca]sub-clamp simulations. C. An example showing Vm and [Ca]sub versus time with [Ca]sub being clamped at the lowest (red) and highest (blue) levels.
In Fig 2A, the EAD behaviors are obtained from the same initial condition. It is known that the cardiac AP dynamics depends on initial conditions as well as pacing protocols [19, 21, 39, 40]. To show the results in Fig 2A are representative, we also use two other pacing protocols to assess the EAD behaviors. In the first protocol, the cell is paced with a fixed pacing period, called the S1S1 pacing protocol. Fig 3A is a phase diagram obtained using this pacing protocol with the S1S1 interval being 3 s. Fig 3C plots Vm and [Ca]sub versus time showing an example of Vm-driven EADs. In the second protocol, the cell is first paced with a fixed pacing period for many beats and then the pacing interval is lengthened for the last pacing beat, called the S1S2 pacing protocol. Fig 3B is a phase diagram obtained using this pacing protocol with the S1S1 interval being 1 s and the S1S2 interval being 2 s. Fig 3D plots Vm and [Ca]sub versus time showing an example of Vm-driven EADs. Note that the grey regions in both cases are defined when the AP duration (APD) is longer than the S1S1 interval. This is because when APD is longer than the pacing period, it can be true repolarization failure, but it can also be 2:1 block and other complex behaviors (such as chaos) which complicate the identification of the EAD mechanisms. Nevertheless, although the parameters for the EAD behaviors are different, the structures of the phase diagrams obtained using these two pacing protocols are similar to the one in Fig 2A. This indicates that using the single initial condition with a single stimulus is an appropriate stimulation protocol for dissecting the EAD mechanisms. We will use this protocol for all simulations shown later in this study.
A. Phase diagram obtained with the S1S1 protocol. The parameter set is the same as Fig 2A. S1S1 interval = 3 s. The grey region is determined using the criterion of APD being greater than the S1S1 interval. The colored regions are different EAD regions classified by the same way as in Fig 2A. B. Phase diagram obtained with the S1S2 protocol. The EAD classification is done using the S2 beat. The S1S1 interval is 1 s and the S1S2 interval is 2 s. C. An example showing Vm and [Ca]sub versus time from the regime of Vm-driven EADs with the S1S1 protocol. α(GKs)=0.8 and α(PCa)=4.4. The traces in the right panels are those when [Ca]sub is clamped in the last pacing beat. D. An example showing Vm and [Ca]sub versus time from the regime of Vm-driven EADs with the S1S2 protocol. α(GKs)=0.8 and α(PCa)=4.6. The traces in the right panels are those when [Ca]sub is clamped in the S2 beat.
EADs driven by Ca oscillations (type II)
In the magenta region in Fig 2A where PCa is low, EADs are caused by Ca oscillations alone. Since the EAD amplitude in the parameter regime of Fig 2A is relatively small, we select another parameter set for α(kmax), α(JCaslmyo), (all are increased) to better illustrate the properties of this type of EADs. Fig 4A is a phase diagram showing the EAD behaviors in the α(PCa) − α(GKs) plane. A representative case (star in Fig 4A) is shown in Fig 4B, where Vm, [Ca]sub, and INCX oscillate completely in-phase (aligned by the vertical dashed lines), whereas ICa,L oscillates out-of-phase with them. This implies that the EADs are originated from the Ca oscillations not by reactivation of ICa,L as in the EADs driven by Vm oscillations shown in Fig 2. To further verify this, we carried out simulations under different clamp conditions. Fig 4C shows Vm and [Ca]sub versus time for [Ca]sub is clamped during the EAD phase, and no Vm oscillations could occur. On the other hand, if Vm is clamped during the EAD phase (Fig 4D), [Ca]sub still oscillates. Moreover, if ICa,L is clamped at a constant level (horizontal dashed line in Fig 4B), [Ca]sub and Vm oscillations still maintain (upper panel in Fig 4E). If INCX is clamped at a constant level (horizontal dashed line in Fig 4B), [Ca]sub still oscillates but not Vm (lower panel in Fig 4E). Therefore, these results demonstrate that Ca oscillations are the primary driver of the EADs via INCX, agreeing with the experimental observations [32].
α(kmax) = 7, α(JCaslmyo) = 3.6, and α(INCX)=2.2. A. Phase diagram in the α(PCa) − α(GKs) plane. EADs occur in the magenta region. B. An example showing Vm, [Ca]sub, ICa,L, and INCX versus time from the EAD region [α(GKs) = 0.56, α(PCa)=2.8, marked as star in A]. The first two vertical dashed lines indicate the takeoff of the first two EADs and the third vertical dashed line marks the peak of the EAD prior to the last one. C. The same case as in B but [Ca]sub is clamped at two different levels in the plateau (grey patch in the [Ca]sub panel in B). D. The same case as in B but Vm is clamped at two different levels in the plateau (grey patch in the Vm panel in B). E. Upper panel: Vm and [Ca]sub versus time when ICa,L is clamped at a constant in the plateau (dashed line in the ICa,L panel in B). Lower panel: Vm and [Ca]sub versus time when INCX is clamped at a constant in the plateau (dashed line in the INCX panel in B).
EADs driven by the Ca-ICa,L-Vm-Ca feedback loop (type III)
When PCa is large (the green region in Fig 2A), a new type of EAD mechanism is observed. Fig 5A is a representative case of this type of EADs (the square in Fig 2A). In this case, unlike the one driven by Ca oscillations, Vm and [Ca]sub oscillations are out-of-phase. When either [Ca]sub (Fig 5B) or Vm (Fig 5C) is clamped, neither Vm nor Ca exhibits oscillations. This implies that neither the Vm system alone nor the Ca system alone can oscillate, but coupling of Vm and Ca promotes the instability for oscillations.
A. Vm, [Ca]sub, ICa,L, and INCX versus time for a representative case [α(GKs) = 1.6 and α(PCa) = 8.2, square in Fig 2A]. The left vertical dashed line indicates the peak and the right one indicates the valley of [Ca]sub. B. Vm and [Ca]sub versus time when [Ca]sub is clamped at the lowest (red), middle (green), and highest (blue) levels during the plateau phase. C. Vm and [Ca]sub versus time when Vm is clamped at lowest (red) and highest (blue) levels in the plateau phase. D. [Ca]sub versus Vm from simulations using the Vm clamp protocol in which Vm is switched from -80 mV to different clamped levels. Black dots: INCX = 0. Red dots: ICa,L = 0. The inset is a schematic plot of the feedback loop.
We hypothesize that the mechanism of this type of EADs is driven by a feedback loop between Ca and Vm mediated by ICa,L. Specifically, elevation of [Ca]sub suppresses ICa,L via Ca-dependent inactivation (via fCa gate in the model). In the supplemental Fig A in S1 Text, we plot peak ICa,L versus [Ca]sub under different holding Vm for the WG model and other models, which shows that the ICa,L amplitude decreases as [Ca]sub increases for holding Vm < 40 mV. Therefore, depending on the relative contributions of ICa,L and INCX (INCX is overall smaller than ICa,L, compare Fig A with Fig D in S1 Text), elevation of [Ca]sub may affect Vm negatively. Conversely, Vm also affects [Ca]sub by modulating ICa,L and INCX. In the WG model, increasing holding Vm from -40 mV to +40 mV first decreases [Ca]sub (Vm from -40 mV to -25 mV), then increases [Ca]sub (Vm from -25 mV to +10 mV) and decreases [Ca]sub after Vm is above +10 mV (see the black dots in Fig 5D). This is caused by the effect of Vm on ICa,L in which the ICa,L magnitude exhibits the same relationship with Vm [see Fig B(f) in S1 Text]. This characteristic stems from the steady-state inactivation curve (fss) of ICa,L. These results show that in the Vm range from -40 mV to +10 mV, [Ca]sub affects ICa,L negatively (Fig A in S1 Text) and thus the [Ca]sub effects on Vm are negative, but the Vm effects on [Ca]sub can be either positive or negative (Fig 5D). Therefore, the feedback between Vm and Ca in the Vm range from -40 mV to +10 mV can be either positive or negative. Without a rigorous nonlinear dynamics analysis (the AP model is too complex to be useful for performing such an analysis), we cannot pinpoint the exact underlying mechanism. Nevertheless, one can postulate that these feedback loops may play important roles in the genesis of this type of EADs and understand this loosely as follows. In the AP plateau, a high ICa,L maintains Vm in the plateau. Elevation of [Ca]sub suppresses ICa,L so that Vm repolarizes partially (see the left vertical line in Fig 5A which aligns with the peak of [Ca]sub). Then [Ca]sub falls back in response to Vm repolarization. Lower Vm can cause ICa,L reactivation, and thereby depolarizes Vm back and elevates [Ca]sub again (see the right vertical line in Fig 5A which aligns with the valley of [Ca]sub). This chain process repeats until [Ca]sub accumulates to a high level which overly suppresses ICa,L and repolarizes the AP. The effects of [Ca]sub on ICa,L and thus Vm can be seen in the [Ca]sub clamp simulation shown in Fig 5B, in which when [Ca]sub is clamped at a low level, Vm fails to repolarize, but clamped at a high level, Vm repolarizes normally.
EADs driven by a large Ca transient (type IV)
If we increase kmax (increase SR release) and/or reduce JCaslmyo (slow Ca diffusion from the submembrane space to cytoplasm) to increase Ca concentration in the submembrane space, another type of EADs occurs. This type of EADs occurs in the cyan strip in the phase diagram shown in Fig 6A. The characteristic feature of this type of EADs is that there is no oscillation in [Ca]sub but an EAD occurs in the AP (Fig 6B). The occurrence of the EAD is sensitive to [Ca]sub and if we reduce [Ca]sub (which is a clamped trace from the trace in Fig 6B with a 50% amplitude), no EAD occurs (Fig 6C). If we reduce either ICa,L or INCX, the EAD amplitude is attenuated. Therefore, the mechanism of this type of EADs is as follows: Enhanced [Ca]sub increases INCX, which helps the reactivation of ICa,L to reverse the Vm decay, generating an EAD. This mechanism of EADs may be responsible for the late phase 3 EADs observed in experiments [41, 42] in which a large contraction or prolonged Ca transient causes phase 3 triggered activity or EADs.
α(kmax) = 7, α(JCaslmyo) = 1, and α(INCX) = 1.2. A. Phase diagram showing different mechanisms of EAD on the α(PCa) − α(GKs) plane. The Ca transient driven EADs occur in the cyan region. Other colors are the same as indicated in the phase diagram in Fig 2A. B. Vm, [Ca]sub, ICa,L, and INCX versus time for a representative case [α(GKs) = 1 and α(PCa) = 1, diamond in A]. The vertical dashed line indicates the EAD takeoff moment. C. Same as B but [Ca]sub is clamped to 0.5 times of the trace in B (highlighted by red). The blue dashed Vm trace in the Vm panel is the original one for comparison. D. Same as B but ICa,L is reduced linearly from the EAD takeoff moment. E. Same as B but INCX is reduced linearly from the EAD takeoff moment.
EAD properties of other ventricular AP models
We also investigate EAD properties using other ventricular AP models, including the phase II Luo and Rudy guinea pig model (LRd [43], later modified by Faber and Rudy [44]), Mahajan rabbit model [45] (later modified by Huang et al [46], HUCLA), Ten Tuscher and Panfilov human model (TP04) [47], O’Hara and Rudy human model (ORd) [48], and Grandi and Bers human model (GB) [49]. These models belong to the second generation of cardiac AP models which incorporate detailed intracellular Ca cycling. Instead of plotting/scanning phase diagrams as for the WG model investigated above, we use a Monte Carlo approach by randomly assigning parameters between 0.1–10 times of their control values (see S1 Text). The types of EADs are also identified by the analyses similar to that used for the WG model. The results are summarized in Table 1, and the typical APs are shown in supplemental Figs F-M in S1 Text. The WG model can exhibit 4 types of EADs, as already shown above. As expected, all models exhibit Vm-driven (type I) EADs. The HUCLA model can exhibit Ca transient driven (type IV) EADs with a very low probability, and the EAD amplitudes are extremely small (see Fig H in S1 Text). The GB model can exhibit both type II and type III EADs (see Figs L and M in S1 Text).
For each model, a total of 10000 randomly selected parameter sets is simulated, and the number of parameter sets for which the APs exhibit EADs is indicated in the parentheses.
According to the results in Table 1, one can divide the models into two groups: I) the LRd model, the HUCLA model, the TP04 model and the ORd model; and II) the GB model and the WG model. In group I, the Ca system is stable, which cannot exhibit spontaneous Ca oscillations (it can passively follow Vm oscillations as in type I EADs), whereas in group II, the Ca system can oscillate spontaneously. A major difference distinguishing the two groups is the steady-state inactivation curve (fss) for ICa,L. In group I, fss is a typical decreasing sigmoidal function, but in group II, fss first decreases but then increases with Vm, which is a non-monotonic function (see the different fss curves in Fig B in S1 Text). This difference in fss causes different responses of ICa,L to Vm and thus [Ca]sub to Vm in the GB model and WG model from those in the group I (see Figs B and C in S1 Text), which allow spontaneous Ca oscillations to occur in group II. Moreover, the two models in group II also have similar Ca cycling models which are based on the SB model [38]. If we flatten the late elevated portion of fss in the GB model and the WG model, type II and type III EADs cannot be observed (see Table 1).
Discussion and conclusions
In this study, we carry out computer simulations of AP models with detailed intracellular Ca cycling to investigate the mechanisms of EADs. More specifically, we investigate the roles of Ca cycling and its coupling with Vm in the genesis of EADs in cardiac myocytes. We first use the WG model and identify 4 mechanisms of EADs:
- The first type of EADs is driven by oscillations originating from the Vm system. The condition for this to occur is called RRR in which outward currents are reduced and/or inward currents are increased. RRR combined with proper ICa,L kinetics provides a condition for a Hopf bifurcation to occur, leading to Vm oscillations that manifest as EADs. The rise of Vm in the EADs is mediated by the reactivation of ICa,L and this type of EADs needs a high ICa,L conductance (see Fig 2A). This is the type of EAD mechanism widely investigated in theoretical and simulation studies [20–22, 24], which may account for the mechanism for the majority of EADs observed in experimental studies [8, 9].
- The second type of EADs is driven by oscillations originating from the Ca cycling system. This type of EADs tends to occur when SR release and INCX are strong, but a larger than control ICa,L is still needed (Figs 2A and 4A). Ca oscillations are promoted by strong SR release and the elevated Ca combined with a strong INCX gives rise to Vm depolarizations to manifest as EADs. Reactivation of ICa,L is not required, however, a strong ICa,L is still needed to hold Vm at the plateau phase for Ca and INCX to cause Vm depolarizations. This type of EADs has been investigated in computer simulation studies [33–35] and experimental evidence has been shown in several experimental studies [29, 32]. However, a rigorous nonlinear dynamics analysis is still lacking to pinpoint the exact nonlinear dynamics/bifurcations for this mechanism of EADs.
- The third type of EADs is driven by a feedback loop between Vm and Ca mediated mainly by ICa,L. This type of EADs tends to occur when ICa,L conductance is very high (Fig 2A). In this case, oscillations cannot be generated by either Vm or Ca alone but are caused by their coupling. The exact mechanism of oscillations is unknown, and we postulate that the feedback between Vm and Ca combined with a special steady-state inactivation curve of ICa,L is the underlying driver for this type of EADs. To the best of our knowledge, there is no experimental evidence on this mechanism of EADs. Moreover, even if EADs of this mechanism occur in an experimental system, it will be difficult to distinguish them from the ones driven by the Vm oscillations or the Ca oscillations.
- The fourth type of EADs is driven by an enhanced Ca transient. This type of EADs occurs when the SR Ca release is very strong to result in an elevated and prolonged Ca transient. However, a stronger than normal ICa,L or INCX is not a requirement (Fig 6A) but both play an important role. Only a single EAD can be generated in an AP, and it is not an oscillation but a transient phenomenon. Experimental evidence for this EAD mechanism may be supported by the late phase 3 EADs observed in experimental studies [41, 42] in which the occurrence of a phase 3 EAD or triggered activity is companied by a large and prolonged contraction or Ca transient.
We also carry out simulations of 5 other AP models to assess the 4 mechanisms of EADs. We use a Monte Carlo approach to randomly select parameters in pre-assigned wide ranges. All models can exhibit the Vm oscillation driven (type I) EADs but only the GB and WG models can exhibit the Ca oscillation driven (type II) and Vm-Ca feedback loop driven (type III) EADs. We identify that the difference in the steady-state inactivation curve of ICa,L may be responsible for the different behaviors of the models since if we change the steady-state inactivation curve in WG (or GB) model to the same as in the other models, we are not able to observe both type II and type III EADs. Another difference is that the Ca cycling models in the WG model and the GB model are based from the one in the SB model which allows spontaneous Ca oscillations. The results of these simulations also raise a question on proper modeling of Ca cycling in the AP models, which needs to be addressed properly using advanced modeling approaches [50, 51].
Although we use computer simulations to scan parameters for phase diagrams and use a Monte Carlo approach to select parameters in wide ranges, rigorous nonlinear dynamics analyses are still needed to pinpoint the exact dynamics for EAD genesis. So far, such analyses have been mainly done for Vm oscillation driven EADs. However, such analyses have been done mostly in simplified models with no Ca cycling. In the presence of Ca cycling, the models become too complex to be useful for such analyses, and simplified models are needed in future bifurcation analyses. In most of the simulations, we use a single set of initial condition and a single stimulus to elicit an AP. We then scan the parameters to assess the EAD behaviors. However, the heart is under constant pacing, and a change in parameters will alter its equilibrium state, particularly those of the intracellular ions. Since the equilibrium state is heartrate dependent, and thus the EAD behaviors can also be heartrate dependent. Nevertheless, using different pacing protocols (Fig 3), we obtained similar phase diagrams to that using the single initial condition, indicating that although the EAD behaviors may be heartrate dependent, our simple pacing protocol can still reveal the underlying mechanisms of EADs.
Besides the insights gained on dynamical mechanisms of EAD genesis, our study may also provide insights into EAD genesis in pathological conditions. EADs are a hallmark of LQTS in which either inward currents are increased or outward currents are decreased, or both occur. Due to the main changes are membrane currents to prolong APD, the EAD mechanisms are mainly via Vm-driven oscillations (type I) [52]. However, experimental evidence of Ca-driven oscillations might also be responsible for EADs in LQTS [32]. EADs are widely observed in heart failure [53, 54]. In heart failure, remodeling occurs in both membrane ionic currents and Ca cycling [55], it is possible that types I-III EAD mechanisms can occur depending on specific heart failure conditions. EADs are also observed in catecholaminergic polymorphic ventricular tachycardia in which RyR mutations cause leaking RyRs to alter the Ca cycling behaviors and promote spontaneous Ca release [56, 57]. EADs under this diseased condition are likely caused by the Ca-driven oscillations (type II mechanism), however, this mechanism still requires an increase in ICa,L (Fig 2A). In another diseased condition in which RyR mutations causes a reduction of RyR open probability and EADs [58]. Reduced RyR open probability increases the SR load which then results in a large Ca transient, and EADs can be caused by the large Ca transient as in the type IV EAD mechanism.
In conclusion, using computer simulations of AP models with detailed Ca cycling, we identify 4 mechanisms of EAD genesis arising from Vm oscillations or Ca oscillations, or their coupling and feedback loops. The insights from this study provide theoretical bases for the EADs observed in experiments and deepen our understanding of the mechanisms of EADs and cardiac arrhythmogenesis.
Supporting information
S1 Text. Supplemental information.
The supplemental information contains the following information: i) Simulation results of bidirectional regulations of Vm and submembrane Ca in different models; ii) the types of EADs; and iii) Control parameters and initial values of each model.
https://doi.org/10.1371/journal.pcbi.1011930.s001
(PDF)
References
- 1. Cranefield PF. Action potentials, afterdepolarizations, and arrhythmias. Circ Res. 1977; 41(4):415–423. pmid:409566
- 2. Qu Z, Xie L-H, Olcese R, Karagueuzian HS, Chen P-S, Garfinkel A, Weiss JN. Early afterdepolarizations in cardiac myocytes: beyond reduced repolarization reserve. Cardiovasc Res. 2013; 99:6–15. pmid:23619423
- 3. Rosen MR, Moak JP, Damiano B. The clinical relevance of afterdepolarizations. Ann N Y Acad Sci. 1984; 427:84–93. pmid:6378020
- 4. Cranefiled PF, Aronson RS. Torsade de pointes and other pause-induced ventricular tachycardias: The short-long-short sequence and early afterdepolarizations. Pacing Clin Electrophysiol. 1988; 11(6):670–678.
- 5. Roden DM. Taking the “Idio” out of “Idiosyncratic”: Predicting Torsadesde Pointes. Pacing Clin Electrophysiol. 1998; 21(5):1029–1034. pmid:9604234
- 6. Weiss JN, Garfinkel A, Karagueuzian HS, Chen P-S, Qu Z. Early afterdepolarizations and cardiac arrhythmias. Heart Rhythm. 2010; 7(12):1891–1899. pmid:20868774
- 7. Qu Z, Weiss JN. Mechanisms of ventricular arrhythmias: From molecular fluctuation to electrical turbulence. Ann Rev Physiol. 2015; 77(1):29–55. pmid:25340965
- 8. January CT, Riddle JM, Salata JJ. A model for early afterdepolarizations: induction with the Ca2+ channel agonist Bay K 8644. Circ Res. 1988; 62(3):563–571. pmid:2449297
- 9. January CT, Riddle JM. Early afterdepolarizations: mechanism of induction and block. A role for L-type Ca2+ current. Circ Res. 1989; 64(5):977–990. pmid:2468430
- 10. Song Y, Thedford S, Lerman BB, Belardinelli L. Adenosine-sensitive afterdepolarizations and triggered activity in guinea pig ventricular myocytes. Circ Res. 1992;70:743–753. pmid:1551200
- 11. Gilmour RF, Moise NS. Triggered activity as a mechanism for inherited ventricular arrhythmias in german shepherd dogs. J Am Coll Cardiol. 1996; 27:1526–1533. pmid:8626969
- 12. Li GR, Lau CP, Ducharme A, Tardif JC, Nattel S. Transmural action potential and ionic current remodeling in ventricles of failing canine hearts. Am J Physiol Heart Circ Physiol. 2002; 283(3):H1031–1041. pmid:12181133
- 13. Sanguinetti MC, Tristani-Firouzi M. M. hERG potassium channels and cardiac arrhythmia. Nature. 2006; 440:463–469. pmid:16554806
- 14. January CT, Chau V, Makielski JC. Triggered activity in the heart: cellular mechanisms of early afterdepolarizations. Eur Heart J. 1991; 12 Suppl F:4–9. pmid:1725155
- 15. Zeng J, Rudy Y. Early afterdepolarizations in cardiac myocytes: mechanism and rate dependence. Biophys J. 1995; 68(3):949–964. pmid:7538806
- 16. Clancy CE and Rudy Y. Linking a genetic defect to its cellular phenotype in a cardiac arrhythmia. Nature. 1999; 400:566–569. pmid:10448858
- 17. Saucerman JJ, Healy SN, Belik ME, Puglisi JL, McCulloch AD. Proarrhythmic Consequences of a KCNQ1 AKAP-Binding Domain Mutation. Circ Res. 2004; 95(12):1216–1224. pmid:15528464
- 18. Tanskanen AJ, Greenstein JL, O’Rourke B, Winslow RL. The Role of Stochastic and Modal Gating of Cardiac L-Type Ca2+ Channels on Early After-Depolarizations. Biophys J. 2005; 88(1):85–95. pmid:15501946
- 19. Sato D, Xie L-H, Sovari AA, Tran DX, Morita N, Xie F, Karagueuzian HS, Garfinkel A, Weiss JN, Qu Z. Synchronization of chaotic early afterdepolarizations in the genesis of cardiac arrhythmias. Proc Natl Acad Sci U S A. 2009; 106(9):2983–2988. pmid:19218447
- 20. Huang X, Song Z, Qu Z. Determinants of early afterdepolarization properties in ventricular myocyte models. PLoS Comput Biol. 2018; 14(11):e1006382. pmid:30475801
- 21. Barrio R, Martinez MA, Serrano S, Pueyo E. Dynamical mechanism for generation of arrhythmogenic early afterdepolarizations in cardiac myocytes: Insights from in silico electrophysiological models. Phys Rev E. 2022; 106:024402. pmid:36109976
- 22. Tran DX, Sato D, Yochelis A, Weiss JN, Garfinkel A, Qu Z. Bifurcation and Chaos in a Model of Cardiac Early Afterdepolarizations. Phys Rev Lett. 2009; 102(25):258103. pmid:19659123
- 23. Xie Y, Izu LT, Bers DM, Sato D. Arrhythmogenic transient dynamics in cardiac myocytes. Biophys J. 2014; 106(6):1391–1397. pmid:24655514
- 24. Kugler P. Early afterdepolarizations with growing amplitudes via delayed subcritical Hopf bifurcations and unstable manifolds of saddle foci in cardiac action potential dynamics. PLoS ONE. 2016; 11(3):e0151178. pmid:26977805
- 25. Kurata Y, Tsumoto K, Hayashi K, Hisatome I, Tanida M, Kuda Y. Dynamical mechanisms of phase-2 early afterdepolarizations in human ventricular myocytes: insights from bifurcation analyses of two mathematical models. Am J Physiol Heart Circ Physiol. 2017; 312(1):H106–H127. pmid:27836893
- 26. Kimrey J, Vo T, Bertram R. Canard analysis reveals why a large Ca2+ window current promotes early afterdepolarizations in cardiac myocytes. PLoS Comput Biol. 2020; 16(11):e1008341. pmid:33147207
- 27. Barrio R, Martinez MA, Pueyo E, Serrano S. Dynamical analysis of early afterdepolarization patterns in a biophysically detailed cardiac model. Chaos. 2021; 31(7):073137. pmid:34340346
- 28. Volders PG, Kulcsar A, Vos MA, Sipido KR, Wellens HJ, Lazzara R, Szabo B. Similarities between early and delayed afterdepolarizations induced by isoproterenol in canine ventricular myocytes. Cardiovasc Res. 1997; 34:348–359. pmid:9205549
- 29. Zhao Z, Wen H, Fefelova N, Allen C, Baba A, Matsuda T, Xie LH. Revisiting the ionic mechanisms of early afterdepolarizations in cardiomyocytes: predominant by Ca waves or Ca currents? Am J Physiol Heart Circ Physiol. 2012; 302:H1636–1644. pmid:22307670
- 30. Priori SG, Corr PB. Mechanisms underlying early and delayed afterdepolarizations induced by catecholamines. Am J Physiol Heart Circ Physiol. 1990; 258:H1796–1805. pmid:2163219
- 31. Volders PGA, Vos MA, Szabo B, Sipido KR, Marieke de Groot SH, Gorgels APM, Wellens HJJ, Lazzara R. Progress in the understanding of cardiac early afterdepolarizations and torsades de pointes: time to revise current concepts. Cardiovasc Res. 2000; 46(3):376–392. pmid:10912449
- 32. Choi BR, Burton F, Salama G. Cytosolic Ca2+ triggers early afterdepolarizations and Torsade de pointes in rabbit hearts with type 2 long QT syndrome. J Physiol. 2002; 543:615–631. pmid:12205194
- 33. Song Z, Ko CY, Nivala M, Weiss JN, Qu Z. Calcium-voltage coupling in the genesis of early and delayed afterdepolarizations in cardiac myocytes. Biophys J. 2015; 108:1908–1921. pmid:25902431
- 34. Wilson D, Ermentrout B, Nemec J, Salama G. A model of cardiac ryanodine receptor gating predicts experimental Ca2+-dynamics and Ca2+-triggered arrhythmia in the long QT syndrome. Chaos. 2017; 27(9):093940. pmid:28964110
- 35. Kurata Y, Tsumoto K, Hayashi K, Hisatome I, Kuda Y, Tanida M. Multiple dynamical mechanisms of phase-2 early afterdepolarizations in a human ventricular myocyte model: Involvement of spontaneous SR Ca2+ release. Front Physiol. 2020; 10:1545. pmid:31998140
- 36. Shiferaw Y, Sato D, Karma A. Coupled dynamics of voltage and calcium in paced cardiac cells. Phys Rev E. 2005; 71:021903. pmid:15783348
- 37. Qu Z, Shiferaw Y, Weiss JN. Nonlinear dynamics of cardiac excitation-contraction coupling: an iterated map study. Phys Rev E. 2007; 75:011927. pmid:17358204
- 38. Shannon TR, Wang F, Puglisi J, Weber C, Bers DM. A mathematical treatment of integrated Ca dynamics within the ventricular myocyte. Biophys J. 2004; 87:3351–3371. pmid:15347581
- 39. Huang X, Qian Y, Zhang X, Hu G. Hysteresis and bistability in periodically paced cardiac tissue. Phys Rev E. 2010; 81:051903. pmid:20866257
- 40. Tsumoto K, Kurata Y, Furutani K, Kurachi Y. Hysteretic Dynamics of MultiStable Early Afterdepolarisations with Repolarisation Reserve Attenuation: A Potential Dynamical Mechanism for Cardiac Arrhythmias. Sci Reps. 2017; 7:10711. pmid:28883639
- 41. Burashnikov A, Antzelevitch C. Reinduction of atrial fibrillation immediately after termination of the arrhythmia is mediated by late phase 3 early afterdepolarization-induced triggered activity. Circulation. 2003; 107:2355–2360. pmid:12695296
- 42. Maruyama M, Ai T, Chua S-K, Park H-W, Lee Y-S, Shen MJ, Chang P-C, Lin S-F, Chen P-S. Hypokalemia promotes late phase 3 early afterdepolarization and recurrent ventricular fibrillation during isoproterenol infusion in Langendorff perfused rabbit ventricles. Heart Rhythm. 2014; 11:697–706. pmid:24378768
- 43. Luo CH, Rudy Y. A modelof the ventricular cardiac actionpotential: depolarization, repolarization, and their interaction. Circ Res. 1991; 68(6):1501–1526. pmid:1709839
- 44. Faber GM, Rudy Y. Action Potential and Contractility Changes in [Na+]i Overloaded Cardiac Myocytes: A Simulation Study. Biophys J. 2000; 78(5):2392–2404. pmid:10777735
- 45. Mahajan A, Shiferaw Y, Sato D, Baher A, Olcese R, Xie L-H, Yang M-J, Chen P-S, Restrepo JG, Karma A, Garfinkel A, Qu Z, Weiss JN. A Rabbit Ventricular Action Potential Model Replicating Cardiac Dynamics at Rapid Heart Rates. Biophys J. 2008; 94:392–410. pmid:18160660
- 46. Huang X, Kim TY, Koren G, Choi BR, Qu Z. Spontaneous initiation of premature ventricular complexes and arrhythmias in type2 long QT syndrome. Am J Physiol Heart Circ Physiol. 2016; 311(6):H1470–H1484. pmid:27765749
- 47. ten Tusscher KH, Noble D, Noble PJ, Panfilov AV. A model for human ventricular tissue. Am J Physiol Heart Circ Physiol. 2004; 286(4):H1573–H1589. pmid:14656705
- 48. O’Hara T, Virag L, Varro A, Rudy Y. Simulation of the undiseased human cardiac ventricular action potential: model formulation and experimental validation. PLoS Comput Biol. 2011; 7(5):e1002061. pmid:21637795
- 49. Grandi E, Pasqualini FS, Bers DM. A novel computational modelof the human ventricular action potential and Ca transient. J Mol Cell Cardiol. 2010; 48(1):112–121. pmid:19835882
- 50. Qu Z, Yan D, Song Z. Modeling Calcium Cycling in the Heart: Progress, Pitfalls, and Challenges. Biomolecules. 2022; 12:1686. pmid:36421700
- 51. Colman MA, Alvarez-Lacalle E, Echebarria B, Sato D, Sutanto H, Heijman J. Multi-scale computational modeling of spatial calcium handling from nanodomain to whole-heart: Overview and perspectives. Front Physiol. 2022; 13:836622. pmid:35370783
- 52. Liu GX, Choi BR, Ziv O, Li W, de Lange E, Qu Z, Koren G. Differential conditions for early after-depolarizations and triggered activity in cardiomyocytes derived from transgenic LQT1 and LQT2 rabbits. J Physiol. 2012; 590:1171–1180. pmid:22183728
- 53. Nuss HB, Kaab S, Kass DA, Tomaselli GF, Marban E. Cellular basis of ventricular arrhythmias and abnormal automaticity in heart failure. Am J Physiol Heart Circ Physiol. 1999; 277:H80–H91. pmid:10409185
- 54. Maltsev VA, Silverman N, Sabbah HN, Undrovinas AI. Chronic heart failure slows late sodium current in human and canine ventricular myocytes: Implications for repolarization variability. European Journal of Heart Failure. 2007; 9:219–227. pmid:17067855
- 55. Tomaselli GF, Marban E. Electrophysiological remodeling in hypertrophy and heart failure. Cardiovasc Res. 1999; 42:270–283. pmid:10533566
- 56. Liu N, Denegri M, Dun W, Boncompagni S, Lodola F, Protasi F, Napolitano C, Boyden PA, Priori SG. Abnormal Propagation of Calcium Waves and Ultrastructural Remodeling in Recessive Catecholaminergic Polymorphic Ventricular Tachycardia. Circ Res. 2013; 113:142–152. pmid:23674379
- 57. Zhao YT, Valdivia CR, Gurrola GB, Powers PP, Willis BC, Moss RL, Jalife J, Valdivia HH. Arrhythmogenesis in a catecholaminergic polymorphic ventricular tachycardia mutation that depresses ryanodine receptor function. Proc Natl Acad Sci U S A. 2015; 112:E1669–1677. pmid:25775566
- 58. Sun B, Yao J, Ni M, Wei J, Zhong X, Guo W, Zhang L, Wang R, Belke D, Chen Y-X, Lieve KVV, Broendberg AK, Roston TM, Blankoff I, Kammeraad JA, von Alvensleben JC, Lazarte J, Vallmitjana A, Bohne LJ, Rose RA, Benitez R, Hove-Madsen L, Napolitano C, Hegele RA, Fill M, Sanatani S, Wilde AAM, Roberts JD, Priori SG, Jensen HK, Chen SRW. Cardiac ryanodine receptor calcium release deficiency syndrome. Science Translational Medicine. 2021; 13:eaba7287. pmid:33536282