- 1Department of Physics and Astronomy, Biological Physics Group, School of Physics & Astronomy, The University of Manchester, Manchester, United Kingdom
- 2Peng Cheng Laboratory, Shenzhen, China
- 3School of Computer Science and Technology, Harbin Institute of Technology, Harbin, China
- 4Shenzhen Key Laboratory of Visual Object Detection and Recognition, Harbin Institute of Technology, Shenzhen, China
Biophysically detailed mathematical models of cardiac electrophysiology provide an alternative to experimental approaches for investigating possible ionic mechanisms underlying the genesis of electrical action potentials and their propagation through the heart. The aim of this study was to develop a biophysically detailed mathematical model of the action potentials of mouse atrial myocytes, a popular experimental model for elucidating molecular and cellular mechanisms of arrhythmogenesis. Based on experimental data from isolated mouse atrial cardiomyocytes, a set of mathematical equations for describing the biophysical properties of membrane ion channel currents, intracellular Ca2+ handling, and Ca2+-calmodulin activated protein kinase II and β-adrenergic signaling pathways were developed. Wherever possible, membrane ion channel currents were modeled using Markov chain formalisms, allowing detailed representation of channel kinetics. The model also considered heterogeneous electrophysiological properties between the left and the right atrial cardiomyocytes. The developed model was validated by its ability to reproduce the characteristics of action potentials and Ca2+ transients, matching quantitatively to experimental data. Using the model, the functional roles of four K+ channel currents in atrial action potential were evaluated by channel block simulations, results of which were quantitatively in agreement with existent experimental data. To conclude, this newly developed model of mouse atrial cardiomyocytes provides a powerful tool for investigating possible ion channel mechanisms of atrial electrical activity at the cellular level and can be further used to investigate mechanisms underlying atrial arrhythmogenesis.
Introduction
Murine hearts are commonly used as animal models in cardiac research for investigating possible molecular and cellular bases of cardiac arrhythmias (Martin et al., 2012; Huang, 2016). Over years, plethora experimental data on cardiac electrophysiology of the mouse heart have been obtained in a wide range of physiological and pathological conditions. Now it is a challenge to integrate these experimental data obtained at different scales into an integrated mathematical model for systematically elucidating better functional roles of ion channels in atrial electrical excitations and arrhythmogenesis. Up to date, a number of mathematical models have been developed for cellular action potentials (APs) of different types of mouse cardiomyocytes. Examples include the first mathematical model of mouse sinoatrial node by Mangoni et al. (2006), which was further updated by incorporating biophysical properties of membrane ionic currents and intracellular Ca2+ handling mechanisms by Kharche et al. (2011), the first mouse ventricular cell models developed by Bondarenko et al. (2004), which were updated by incorporating newer experimental data by Li et al. (2010). Recently, mouse ventricular cell models with β-adrenergic signaling regulation and with calmodulin (CaM) mediating Ca2+-dependent regulation were developed by Yang and Saucerman (2012) and Morotti et al. (2014), respectively. These models provide useful tools for underpinning insights into the regulation of Ca2+ handling in physiological and pathological conditions.
However, experimental data of mouse atrial electrophysiology is spare due to its limited use for the study of atrial fibrillation, as the mouse heart is too small to initiate and maintain atrial arrhythmias (Huang, 2016) in normal conditions. Recently, however, sustained atrial tachycardia and fibrillation in mouse heart has been reported (Wakimoto et al., 2001), and intensely investigated with transgenic or wild type mice (Kovoor et al., 2001; Verheule et al., 2004; Chelu et al., 2009; Choi et al., 2012). These experimental data provide sufficient foundation to develop a mouse atrial cell model. In addition, after decades of experimental studies on the intracellular signaling network (Yang and Saucerman, 2011), sufficient data are also available to construct a biophysically detailed model of mouse atrial myocytes, with consideration of intracellular signaling pathways to simulate the regulation from nervous and endocrine systems (Ramos-Franco et al., 1998; Doi et al., 2005).
Furthermore, electrophysiological heterogeneity within or between the mouse left and right atria has been implicated as a factor initiating and/or maintaining atrial rhythm disturbances (Nygren et al., 2004). Such electrical heterogeneity between left and right mouse atria has also been found in many studies (Lomax et al., 2003a; Trépanier-Boulay et al., 2004; Hu et al., 2005; Nakamura et al., 2010), in which left-to-right differences in several K+ currents, as well as in the action potential duration (APD) of isolated atrial myocytes, were reported. A mouse atrial cell model considering the heterogeneity between left and right atria will be helpful on the study of the genesis of AF, such as to investigate how the fast rotors in the left atrium (LA) propagate toward the right atrium (RA) (Sarmast et al., 2003).
The main goal of this study is to develop a mathematical model of mouse atrial myocytes to simulate mouse atrial APs and intracellular Ca2+ handling mechanisms. In the developed model, Markov chain formalizations were used to describe membrane ion channels for three main potassium currents (the transient outward K+ current (Ito), the ultra-rapidly activating delayed rectifier K+ current (IKur) and the rapid delayed rectifier K+ current (IKr) and the L-type Ca2 + channels (LTCCs). Heterogeneity between left and right atria was considered to derive a left atrial cell model based on the baseline model of right atria. The characteristics of AP and Ca2+ transient (CaT) generated by our model were compared with mouse atrial experimental findings, which validated the model development. The Ca2+-calmodulin activated protein kinase II (CaMKII) and β-adrenergic signaling pathways were also included.
Model Development
Details of model implementation, pacing protocols, and performance of the model can be found in Supplementary Document 1. The source code in C++ of the model is available under the request to the authors (henggui.zhang@manchester.ac.uk).
General Structure
The mouse atrial myocyte model was developed based on the basal model of the mouse ventricular cell model of Morotti et al. (2014) (body temperature), which was extensively modified to match experimental data of mouse atrial cells. The main changes to the original model are described below, and more details are summarized in Supplementary Document 1. The model consists of 159 coupled ordinary differential equations, each of which describes kinetic properties of membrane ion channels, pump and exchanger currents, intracellular ionic homeostasis, and the CaMKII and β-adrenergic signaling pathway. To better describe gating kinetics of ionic channels of Ito, IKur, and IKr, Markov chain models were used. Equations for CaMKII and β-adrenergic signaling pathways are inherited from the Morotti et al. model, with some parameters of the activation of CaMKII and β-adrenergic activation of protein kinase A (PKA), and their phosphorylation on excitation-contraction coupling (ECC) targets, being modified based on atrial experimental findings. A schematic diagram of the mouse atrial myocyte model is shown in Figure 1.
Figure 1. A schematic diagram of the computational model of a mouse atrial myocyte. (A) Ionic channel currents, fluxes, signaling pathways and physical compartments of the mouse atrial model. (B) Diffusion of Ca2+ and CaM in three compartments. CaMKII is activated by Ca2+-CaM binding in the dyadic cleft, sub-sarcolemma area (SL), and bulk cytosol. Phosphatases 1 and phosphatases 2A oppose phosphorylation by either kinase. Ca2+ and CaM can diffuse across these three regions at different diffusion coefficients. The blue and red circle with “P” in it indicates the phosphorylation of the target by PKA or CaMKII, respectively. INa, fast Na+ current; INaL, late Na+ current; ICaL (or LTCC), L-type Ca2+ current; Ito, transient outward K+ current; IKur, ultra-rapidly activating delayed rectifier K+ current; IKr, rapid delayed rectifier K+ current; Iss, non-inactivating steady-state voltage-activated K+ current; IK1, time-independent inwardly rectifying K+ current; IKACh, acetylcholine-activated K+ current; IKCa, small conductance Ca2+-activated K+ current; INab, background Na+ current; ICab, background Ca2+ current; IKb, background K+ current; IClb, background Cl– current; IClCa, Ca2+-activated Cl– current; PMCA, plasma membrane Ca2+ pump; NKA, Na+/K+-ATPase; NCX, Na+/Ca2+ exchanger.
Modeling of Membrane Currents
LTCC
LTCC provides main Ca2+ influx for myocytes and acts as a trigger for Ca2+ release from the sarcoplasmic reticulum (SR). There are four types of LTCCs: Cav1.1 to Cav1.4, among which Cav1.2 plays a major role during cardiac ECC (Hofmann et al., 1999). ICaL has the kinetics of a fast activation and a slower inactivation time course.
The formulation of LTCC in our model was based on a previous Markov chain model developed by Mahajan et al. (2008) (Figure 2A), with CaMKII and PKA phosphorylation modules being updated for mouse ventricles by Morotti et al. (2014). The LTCC model has a minimum scheme of seven states, which incorporates both calcium-dependent inactivation and voltage-dependent inactivation, where the transitions between C2 and C1 are strongly voltage-dependent, and the transitions from C1 to O are voltage-independent and determine the steady-state open probability.
Figure 2. ICaL in the atrial myocyte model and comparison with experimental data. (A) A schematic diagram of the Markov chain model of the L-type Ca2+ channel (Mahajan et al., 2008). (B) Current traces from simulation and experimental records (inset) from Mancarella et al. (2008). (C) Normalized I–V relationship of ICaL. (D) Peak simulated current density compared with experimental data (Rose et al., 1992; Lomax et al., 2003a; Xiao et al., 2006; Mancarella et al., 2008; Xie et al., 2012). All experiments were conducted at room temperature except the data from Xie et al. at body temperature (marked with a red star).
In order to fit the experimental measures, the voltage dependence of the steady-state transition between C2 and C1 was manually shifted left by 4.4 mV (4.4 mV was added on the voltage in the voltage-dependent transition rates α and β shown in Figure 2A). The current traces shown in Figure 2B were simulated by using the same voltage-clamp protocol as described in Mancarella et al. (2008): the ICaL was activated by a series of 250-ms depolarization pulses from a holding potential of −90 mV, to test potentials ranging from −40 to +40 mV (10-mV steps). And the resultant normalized I–V relationship was in good agreements with experimental observations (Mancarella et al., 2008; Hua et al., 2012) (Figure 2C). The total conductance of ICaL was reduced to 1/4 of the ventricular cell model and the simulated peak current density is 4.76 pA/pF, which is in the range of various mouse atrial experimental data shown in Figure 2D.
Ito
Ito is a rapidly activating and inactivating K+ current, encoded by Kv4.2 and Kv4.3 (Nerbonne, 2000). Xu et al. (1999c) reported that the elimination of Ito gave rise to substantial prolonged atrial APD, suggesting that Ito plays a prominent role in mouse atrial repolarization.
The equations for Ito in the mouse atrial myocyte were formulated based on the Markov chain model of Campbell et al. (1993) that comprised seven states as schematically shown in Figure 3A. Transition rates were refitted with a non-manual method (details in Supplementary Document 1). Simulated current traces in Figure 3B were similar to those seen in experiments (the inset in Figure 3B). The current traces were obtained by using a protocol of 500-ms voltage steps in 10-mV increments between −70 and +50 mV from a holding potential of −75 mV. The normalized I–V relationship of our model was then calculated and compared with experimental data (Figure 3C). Simulation of Ito shows a peak current density of 11.6 pA/pF at +30 mV (Figure 3D).
Figure 3. Ito in the atrial myocyte model and comparison with experimental data. (A) Schematic diagram of seven-state Markov chain model of Ito (Campbell et al., 1993). (B) Current traces for the corresponding voltage clamp protocol. The inset shows the experimental record. (C) I–V relationship, current normalized to its amplitude at +30 mV. (D) Current density at +30 mV. (E) Steady-state inactivation curve. (F) Recovery curve from inactivation. (Experimental data: Lomax et al., 2003a; Trépanier-Boulay et al., 2004; Qin et al., 2012).
The curve of steady-state inactivation and recovery from inactivation of Ito were shown in Figures 3E,F. To elucidate the steady-state inactivation curve for Ito, a two-pulse protocol was applied (panel A of Supplementary Figure 1). The experimental current traces were shown in panel B of Supplementary Figure 1. The simulated current traces from our model shown in panel C of Supplementary Figure 1 were more similar to the experiment traces compared with the current traces of the Hodgkin and Huxley (HH) formalization from the parent model (panel D of Supplementary Figure 1), implying a better description of our Markov chain model on the fast inactivation of Ito. The steady-state inactivation curve and the recovery curve from inactivation produced from our model were all similar with experimental findings (Figures 3E,F). A temperature adjustment factor (Q10) used to rescale the model to body temperature was determined separately for channel conductance and kinetics. The Q10 of the channel kinetics were first determined as 1.8 from literature (Courtemanche et al., 1998; Brouillette et al., 2004), and then the Q10 factor for channel conductance was adjusted and set to 1.1 to match the current density at 32°C, which is ∼1.34 times of that at 22°C (Brouillette et al., 2004).
IKur
IKur encoded by Kv1.5, was identified in atrial myocytes from several species (Lomax et al., 2003a; Wettwer et al., 2004; Pandit et al., 2011). The kinetics of rapid activation and slow inactivation were reflected in the early current nomenclature, such as rapidly activating slowly inactivating current IK,s (Bou-Abboud et al., 2000), or ultra-rapid K+ current IKur (Trépanier-Boulay et al., 2004), which is the currently accepted name for this current.
The HH formulations of IKur for the mouse ventricular myocyte were replaced by a six-state Markov chain model from Zhou et al. (2012) (Figure 4A), in which the transition rates α and β were decreased to 1/8 of their original values in order to slow down the activation process of the IKur model to match to experimental data as shown in the Figure 4. The simulated current traces, normalized I–V curve, and peak current density of 5.1 pA/pF at +30 mV were similar to those seen in experiments (Figures 4B–D) (the current traces was simulated by using a protocol of 500-ms voltage steps in 10-mV increments between −70 and +50 mV from a holding potential of −80 mV). Since all experimental data were obtained under room temperature, the transition rates in the model were fitted to these data then adjust to the physiological temperature. Using the same method as for Ito, the Q10 factor of activation and inactivation transition rates of IKur was set to 1.9 according to literature (Brouillette et al., 2004). Then the Q10 factor of channel conductance was chosen as 1.37 to fit the increase in current density at higher temperature (Brouillette et al., 2004; Nakamura et al., 2010), which is about 35 pA/pF at 22°C and 49 pA/pF at 32°C found in mouse ventricles.
Figure 4. IKur in the atrial myocyte model and comparison with experimental data. (A) Schematic diagram of six-state Markov chain model of IKur (Zhou et al., 2012). (B) Current traces for the corresponding voltage clamp protocol. The inset shows the experimental record from Trepanier-Boulay et al. (C) I–V relationship, current densities were normalized to the current at +30 mV. (D) Current density when clamped to +30 mV (Experimental data: Xu et al., 1999c; Brouillette et al., 2004; Trépanier-Boulay et al., 2004; Nakamura et al., 2010; Qin et al., 2012).
IKr
IKr encoded by ERG (KCNH2, Kv11.1), has been proved to have a key role in the late repolarization of AP in several mammalian species (Sanguinetti and Jurkiewicz, 1990; Wang et al., 1994; Li et al., 1996). The presence and function of IKr in the mouse atrium has been demonstrated by Nakamura et al. (2010). Another study also shows that the relative mRNA level of ERG channel was higher and the amplitude of IKr was greater in mouse atria than in ventricles (Mewe et al., 2010).
The HH formulations of IKr in the mouse ventricular myocyte were replaced with the five-state Markov chain model of Clancy and Rudy (2001) (Figure 5A) with modifications on transition rates to fit experimental measures. The model was stimulated with a voltage-clamp protocol shown in Figure 5B and resultant current traces are in Figure 5C. The time constants for inactivation and activation at several test potentials (Figure 5D) and the resultant tail current I–V curve (Figure 5E) were in agreement with experimental data (Nakamura et al., 2010) (the tail current were measured when the clamp voltage returning to −40 mV). The simulated steady-state I–V curve were also shown in Figure 5F. The model was fitted to experiment data at 37°C, so no temperature scaling factor was needed.
Figure 5. IKr in the atrial myocyte model and comparison with experimental data. (A) Five-state Markov chain model of IKr. (B) The protocol used in the simulation. (C) Current traces for the corresponding voltage clamp protocol. (D) Time constants for inactivation and activation. (E) Tail I–V relationship. (F) Steady-state I–V relationship (Experimental data: Nakamura et al., 2010).
Iss
Iss was found in the mouse ventricles and atria by Xu et al. (1999b,c). It is a 4-aminopyridine (4-AP) resistant current, which is distinctively different from Ito and IKur (Brouillette et al., 2004). The activation of Iss is slow and it undergoes little or no inactivation (Xu et al., 1999a).
Iss in our model was reformulated in the form of an HH formula to match its kinetics observed from atrial experimental measures. Since Iss has no inactivating process, one activation gate xss is sufficient to reflect its relatively slow activation. The simulated current traces and normalized I–V relationship shown in Figures 6A,B were determined by using a protocol of 500-ms voltage steps in 10-mV increments between −80 and +60 mV from a holding potential of −80 mV. Considering the wide range of current density from experiments, conductance of Iss was set to 0.049 nS/pF. The resultant peak current density 3.32 pA/pF at 30 mV and the comparison with experimental measurements are shown in Figure 6C. The simulated peak current density, normalized I–V curve and related current traces are similar to those seen in experiments. Experimental data were mostly acquired at 22°C, so a Q10 factor of 2.4 (Brouillette et al., 2004) was applied to the time constant of xss and channel conductance.
Figure 6. Iss, IKb and IK1 in the atrial myocyte model and comparison with experimental data. (A) Simulated current traces of Iss. (B) I–V relationship of Iss, values normalized to the current at +30 mV. (C) Current density of Iss at +40 mV (experimental data: Xu et al., 1999c; Brouillette et al., 2004; Hu et al., 2005). (D) I–V relationship of IKb. (E) I–V relationship for normalized IK1. (F) Current density of IK1 at –110 mV (experimental data: Lomax et al., 2003a; Trépanier-Boulay et al., 2004; Qin et al., 2012).
IKb
In cardiac myocytes, the background potassium current IKb plays an important role in maintaining the polarized resting membrane potential (RMP) and intracellular K+ homeostasis. The recently discovered family of two pore domain K+ channels (K2Ps) have properties well suited to a role in mediating background K+ conductance (Gurney and Manoury, 2008). The model of IKb was described by a single Boltzmann equation fitted to experimental data measured by Syeda et al. (2016) at 37°C (Figure 6D). They applied 10 mM Ba2+ in isolated mouse myocytes and attributed the change in background current to twik-related acid-sensitive K+ (TASK) channels, which are in the K2P channel family. Since TASK channels display very fast activation and inactivation processes (Gurney and Manoury, 2008), no time constant was incorporated in our model.
IK1
The channel of IK1 encoded by Kir2.x stabilizes the resting potential and is responsible for shaping the initial depolarization and final repolarization of the action potential.
The IK1 equation in our model was based on that from DiFrancesco and Noble (1985) with model parameters adapted to mouse experimental data. The 500-ms pulse in the voltage-clamp protocol was from the holding potential of −80 mV to voltages from −150 to −40 mV. Currents were normalized to the magnitude at −100 mV. Simulated results were similar to the experimental data (Figures 6E,F).
Other Currents
INa
The INa accounts for the fast upstroke of APs. Due to the absence of experimental data, in mouse atria, the channel kinetics and conductance were unaltered from the mouse ventricular model of Morotti et al. (2014). The resultant maximal upstroke velocity (dV/dtmax) and AP amplitude were in good agreement with experiments.
INaL
The kinetics of INaL were identical to the formulation described in Morotti et al. (2014), due to the lack of experimental data for the mouse atrium. The conductance of INaL was increased by 68% from the mouse ventricular model to fit the peak current density observed in the experiment (Lemoine et al., 2011).
IKACh
Acetylcholine activates IKACh, an inward rectifying K+ current. The molecular composition of cardiac IKACh channels was identified as a heterotetramer consisting of two Kir3.1 (GIRK1) and two Kir3.4 (GIRK4) channel subunits. In the heart, IKACh channels are primarily expressed in the atrioventricular and sinoatrial node and in the atria (Voigt et al., 2014). The IKACh in our model was reparameterized from a HH formulation in a rabbit sinoatrial node model (Zhang et al., 2002). The simulated I–V relationship was compared with mouse atrial experimental data (Supplementary Figure 4).
IKCa
The small conductance Ca2+-activated K+ current, IKCa, is shown to contribute to repolarization in mouse atrial myocytes (Hancock et al., 2015). Three subtypes of the small conductance Ca2+-activated channels (SK1, SK2, SK3) are the carrier of IKCa, which are all sensitive to a selective SK channel blocker apamin (Adelman et al., 2012). We adopted the model of IKCa from another mouse atrial cell model (Asfaw et al., 2020). The channel trafficking effect stimulated by β-adrenergic signaling system was removed.
INCX
NCX (Na+-Ca2+ exchanger) encoded by the NCX family (NCX1, NCX2 and NCX3), is responsible for maintaining steady intracellular Ca2+ balance. It extrudes Ca2+ from the cytosol along with the plasma membrane Ca2+ pump (PMCA). As a reversible transporter, it also mediates Ca2+ entry in parallel with LTCCs (Blaustein and Lederer, 1999). Due to the lack of experimental results on NCX electrophysiological data for the mouse atria, the formulations for the NCX were kept unchanged.
INaK
The INaK is responsible for maintaining the Na+ and K+ gradients between the cytosol and extracellular medium, and plays an important role in the restoration of the resting potential, and the generation and propagation of APs (Glitsch, 2001).
The conductance of INaK was adjusted according to the different regional expression of INaK in ventricles and atria. Schmidt et al. (1990) found that there were 50% fewer INaK in the atria than in the ventricles or septa by measuring the number of ouabain binding sites in porcine and canine hearts. In experiments conducted by Wang et al. (1993), the INaK activity in guinea pig atria was only one third of that in ventricles. In human atrial myocytes, the INaK activity was nearly 50% lower than that in ventricular myocytes as well (Wang et al., 1996). To adapt this change and also keep the intracellular Na+ concentration in physiological range, the maximum rate for the INaK was decreased by 40% in our model.
Modeling of Intracellular Ca2+ Handling
The intracellular Ca2+ handling in this model is based on the framework developed by Shannon et al. (2004) and adopted by Morotti et al. into their ventricular model. This framework was utilized with some modifications.
Similar to the model of Shannon et al., the cell volume in our model was divided into four compartments (Figure 1B): the SR, the dyadic cleft, the sub-sarcolemma area (SL) and the bulk cytosol. Compared with ventricular cells, atrial cells usually have smaller size, higher percentage of cytosolic compartment, and lower percentage volume of other compartments. Thus, dimensions and parameters for the related structure in our model were modified accordingly (see Supplementary Document 1).
RyR Ca2+ Release (Jrel)
The formulations for ryanodine receptor (RyR) were based on the model of mouse ventricular myocyte (Morotti et al., 2014), where the half maximal effective concentration (EC50) of the Ca2+ concentration in SR ([Ca]SR) on RyR was reduced by 10% to match the fractional Ca2+ release data from experiments.
The RyR Markov chain model has four states (panel A of Supplementary Figure 2), with strong dependence of RyR gating upon dyadic Ca2+ concentration ([Ca]cleft). During a normal calcium induced calcium release (CICR), some RyR channels promptly open (O) and rapidly transit into the inactivated states (I) as Ca2+ diffuses out of the dyadic cleft (panel B of Supplementary Figure 2). These channels then become resting inactivated states (RI) and eventually return to the resting states (R) (panel C of Supplementary Figure 2). The rate of Ca2+ release from the SR was 5.2 mM/s with a Ca2+ leakage rate of 0.36 μM/s (panel D of Supplementary Figure 2).
SR Ca2+ Reuptake (Jup)
The main function of sarco/endoplasmic reticulum Ca2+-ATPase (SERCA) is the uptake of Ca2+ into SR. A proper model of SERCA is important for sustained normal CaTs and is highly related to the frequency dependent acceleration of relaxation (FDAR) of the cell model. The SERCA in our mouse atrial model was based on the model of Morotti et al. (2014) with a few modifications.
There is study showing that the decay of intracellular CaT in mouse atrial cells is far slower than in ventricular cells (Shanmugam et al., 2011) at 0.5-Hz pacing. Furthermore, Picht et al. (2007) suggested in their experiment in mouse ventricles that the Ca2+ removal rate via SERCA was greatly elevated at 4 Hz (2.1 times of that at 0.5 Hz) due to increased maximum turnover rate of SERCA (Vmax) without change in Ca2+ affinity of SERCA (Km). In our model, the Vmax was reduced by 50% and a coefficient KSERCA describing the effect of CaMKII at high pacing frequencies was adopted which is calculated by:
where PhPLB denotes the level of phosphorylation of phospholamban (PLB) in percentage. And the forward mode Ca2+ affinity (Kmf) was increased by 13%.
Heterogeneity of Left and Right Atria
Electrophysiological heterogeneity between mouse left and right atria has been investigated in many experimental works (Lomax et al., 2003a, b; Nygren et al., 2004; Verheule et al., 2004; Nakamura et al., 2010). Specifically, the current density differences between left and right atria of ICaL, Ito, IKr, IKur, and IK1 were measured by Lomax et al. (2003a) and Nakamura et al. (2010), which are listed in Figure 7. The most prominent dissimilarity between left and right atrial cells lies in the current densities of IKur and IK1. There are no significant differences in the recovery from inactivation of Ito (Lomax et al., 2003a) and the time constant of activation of IKr (Nakamura et al., 2010). Considering these experimental findings, a model describing the mouse left atrial AP was derived from the constructed right atrial cell model by increasing the channel conductances of IKur and IK1 by 110 and 70%, respectively. Validation of the model is shown in the “Result” section.
Figure 7. Heterogeneity of current densities between left and right atria. The current densities in left atria were normalized to their counterparts in the right atrium. The values of ICaL, Ito, IKr, IKur, and IK1 were extracted at clamp voltages of 0, 35, 25, 35, and –105 mV, respectively.
Unless stated otherwise, simulations in this paper were based on the right atrial model.
Validation of the Model
The developed model was validated by comparing the characteristics of simulated action potentials with those of experimental data obtained from mouse atrial cells in several perspectives, including the AP characteristics, CaT characteristics (with their rate dependences). Simulation results on the role of K+ currents in the AP were also compared to experimental data for validation or calibration purposes. The rate dependence of intracellular ion concentrations are also shown although cannot be treated as true validation, as the INaK was calibrated to reproduce physiological [K+]i and [Na+]i. Because a large number of experimental data are referred in the “Result” section for validation, in order to avoid confusion, only experimental data explicitly from the mouse atrial cells will be presented in simulation figures for comparison purposes, while the experimental data from other species or other parts of the heart will only be referred in the text.
Results
Simulated AP and Currents
The simulated time course of the mouse atrial AP and major underlying currents are shown in Figure 8A. The AP of left atrial cell model was slightly shorter than that of right atrial cell model due to the significantly greater IKur and IK1, while other potassium currents were smaller during the repolarization, producing a partial compensation to the large IKur and IK1. Figure 8B shows characteristics of simulated AP under 1-Hz pacing frequency, compared with experimental data. The dV/dtmax is about 210 V/s in the atrial model (not shown), which is similar to experimental results (Bao et al., 2016; Syeda et al., 2016). Experimentally measured APDs at 50 and 90% repolarization (APD50 and APD90, respectively) distributed in a wide range, which may be due to different experimental environments, such as temperature, pacing frequency, genetic backgrounds or strains. The computed APD50 by our model, though was marginally short, matched to the low bound of the experimental data range. The normalized APD based on the ratio of APD50 over APD90 were measured and presented (middle panel of Figure 8B), showing that the normalized APD50 from experiments and our model are in good agreement. Experimental data shown in Figure 8B are listed in Supplementary Table 1.
Figure 8. Simulated mouse atrial APs, currents and AP characteristics. (A) Mouse left and right atrial APs (shown as the membrane potential (Vm) against time) and major underlying currents simulated at steady state with 1-Hz pacing. (B) Comparison between simulated (bars) and experimental records (circles) of AP characteristics.
Intracellular Ca2+ Handling
Our model well reproduced intracellular CaTs in mouse atrial cells. The CaTs under steady stimuli of 0.5 Hz were depicted in Figure 9A. When paced at 0.5 Hz, the diastolic Ca2+ concentration in the bulk cytosol ([Ca]i) was about 0.11 μM in our model, which was consistent with ventricular experimental findings (Glukhov et al., 2015). As suggested by Walden et al. (2009), there was no apparent difference in the diastolic [Ca]i between murine ventricular and atrial myocytes. The caffeine-induced Ca2+ transient (Figure 9B) was similar with experiments by Xie et al. (2012). The effect of caffeine was simulated by increasing the opening probability of RyRs by 7.5 times and complete block of SERCA and ICab [the same as configurations in the parent mouse ventricular model (Morotti et al., 2014)]. Characteristics of CaTs were calculated from Figures 9A,B and shown in Figure 9C. The peak of amplitude of CaTs were ∼2.6 times higher than the basal level, which is comparable with various experimental data (Li et al., 2005; Shanmugam et al., 2011; Xie et al., 2012). The time constant for the decay phase and time to peak of CaT were also similar with the values reported in Li et al. (2005); Mancarella et al. (2008), and Shanmugam et al. (2011). The fractional SR Ca2+ release was calculated as the fraction of the maximum [Ca]i before and after application of caffeine.
Figure 9. Atrial intracellular Ca2+ transient from simulation and experiments. (A) CaT under the pacing rate of 0.5 Hz. (B) Simulated caffeine-induced CaTs. (C) [Ca]i amplitude, time to the peak of [Ca]i, time constant for the decay phase of [Ca]i, and fractional Ca2+ release, calculated at 0.5 Hz. (D) Time course of CaTs under 1-, 2-, and 4-, and 8-Hz pacing rate. (E) Amplitude of CaTs at various pacing rates. The diastolic Ca2+ concentration (Dia [Ca]i) was fixed as the value when the cell was paced at 0.5 Hz [experimental data: Xie et al. (2012) at body temperature, Li et al. (2005); Mancarella et al. (2008), Escobar et al. (2004); Walden et al. (2009), and Shanmugam et al. (2011)].
Figure 9D shows the simulation results of CaTs at various pacing frequencies. Although no direct atrial data available to compare, the amplitude of CaT (Figure 9E) decreased along with the increasing pacing frequency as seen in ventricular experiments (Picht et al., 2007). When the pacing frequency was between 1 and 3 Hz, there was no obvious change in the maximum [Ca]i, the reduced amplitude of CaT was due to the elevated diastolic [Ca]i. From 3 to 8 Hz, the total SR Ca2+ release was larger than the Ca2+ uptake by SERCA, leading to lower [Ca]SR which in turn reduced the amplitude of every Ca2+ release until new steady states were reached. The new states had a markedly decreased maximum [Ca]i and increased diastolic [Ca]i, therefore, a frequency-dependent attenuation of CaTs.
Ca2+ fluxes play a key role in the regulation of ECC (Bers, 2000; Eisner et al., 2000). The mouse atrial model can closely resemble the experimental data of Ca2+ fluxes. Supplementary Figure 2 shows the behavior of RyRs during a cardiac cycle. The peak value of Jrel was 95 mmol/(L SR)/s or 2.12 mmol/(L cytosol)/s, which was similar to experimental observation (Lukyanenko et al., 1998). The Ca2+ influx through LTCC peaked at 0.21 mmol/(L cytosol)/s, which was close to the experimental estimates of 0.30 mM/s (Bers, 2001). The Ca2+ fluxes of SERCA, NCX, and sarcolemmal Ca2+ pump are shown in Supplementary Figure 3. In our atrial model, the predominant removal mediator, SERCA, contributed 91.4% of the overall Ca2+ removal, with NCX significantly dropping to 7.2%. These results were closely similar to the experimental observations in the rat atrium (SERCA and NCX contributed 92.6 and 6.13% of the total Ca2+ removal) (Walden et al., 2009).
Rate Dependence of Intracellular Ion Concentrations
The rate dependence in concentrations of intracellular Na+, K+ is shown Figures 10A,B. For the Na+ concentration in bulk cytosol ([Na]i), although no experimental data for mouse atria, our simulated results well fitted to the values found in murine ventricular myocytes. Without pacing, [Na]i of the model stabilized at about 10.8 mM, which is quite similar to that of 10 mM in the rat ventricle (Despa et al., 2002). It has been found in experiments that [Na]i increases along with the increasing pacing rates. In mouse ventricle, [Na]i was 8.24 ± 4.9, 12.3 ± 4, and 15.1 ± 5.5 mM at 0.1, 0.5, and 3 Hz pacing rates, respectively (Wagner et al., 2006). In rat ventricle, [Na]i was 11.2 ± 2.3 and 15 ± 1 mM at 0.5 and 2 Hz pacing rates, respectively. The [Na]i of our model at various pacing frequencies were well distributed in the physiological range of 10–17 mM. The K+ concentration in the bulk cytosol ([K]i) of our model was about 140 mM at resting, and dropped to around 110 mM when paced at 10 Hz. This change in [K]i significantly influence the equilibrium potential of K+ and then the RMP (details in the following section). These results implied that INaK magnitude was properly implemented, since INaK is responsible for pumping in K+ and Na+ out thus keeping the transmembrane gradient of Na+ and K+ concentration.
Figure 10. Steady-state rate dependence of intracellular ion concentrations. Rate dependence of (A) [Na]i, (B) [K]i, and (C) [Cl]i. (D) Rate dependence of the time constant of CaT decay.
The Cl– concentration in the bulk cytosol ([Cl]i) is also dynamic in our model. It changed from 8 mM to 14 mM when the pacing rate was increased from 1 Hz to 10 Hz (Figure 10C). The change in [Cl]i generated a higher equilibrium potential of Cl– and played a role in modifying the behavior of Cl– currents (i.e., Ca2+-activated Cl– current and background Cl– current in our model). At 1 Hz, there was only an influx of the Cl– forming a transient outward current when the cell was activated. At higher pacing rates, equilibrium potential of Cl– was elevated to about −50 mV, leading to a Cl– efflux in the diastolic phase of AP which slightly depolarized the RMP.
Figure 10D shows the rate dependence of the time constant of CaT decay. The decay became more rapid with increasing frequency, which reflected FDAR. At high pacing rates, there was enough CaMKII activated to phosphorylate PLBs (see the section of CaMKII-mediated phosphorylation in Supplementary Document 1 for more details) and enhance SERCA (by a larger coefficient KSERCA described in section “SR Ca2+ Reuptake (Jup)”). This boosted Ca2+ uptake contributed to the faster drop of CaT.
Rate Dependence of Mouse Atrial APs
We further validated our mouse atrial model against experimental data of AP characteristics at different pacing frequencies (Figure 11). According to experimental data from Bao et al. (2016), the APD adaptation to pacing frequency of mouse atrial cells (shown in Figure 11A) is not distinct, which is quite similar with that in mouse ventricular cells (Wu and Anderson, 2002; Wagner et al., 2009). For restitution, we compared APD30 after S1 pacing at CL = 130 ms, followed by a single S2 extra-systolic stimulus delivered at various S1S2 intervals (Figure 11B). The APD30 restitution curve well reproduces the trend of increase when S1S2 intervals were being shortened. The steady-state rate dependence of mouse atrial dV/dtmax and RMP are shown in Figures 11C,D. At pacing frequencies lower than 4 Hz (CL > 250 ms), dV/dtmax and RMP remained almost stable; at higher pacing frequencies, dV/dtmax markedly decreased along with the increase of RMP, which is consistent with experimental data.
Figure 11. Rate dependence of mouse atrial APs. (A) Steady-state rate dependence of simulated APDs (lines) compared with experimental data (dots) from Bao et al. (2016). Red, blue, and green represent APD90, APD50, and APD30, respectively. (B) Simulated APD restitution curve under an S1S2 protocol compared with experimental data from Knollmann et al. (2007). The BCL of S1 was 130 ms. (C) Rate dependence of the maximum upstroke velocity. Data from Syeda et al. are all listed without using means and error bars. (D) Rate dependence of resting membrane potential (RMP) (All data from Bao et al. in panels (A) and (C) have symmetric errors, but some of them are shown with one-side error bars for clarity).
There are two likely causes for the elevated RMP at high pacing frequencies: slower late repolarization (defined as the phase of repolarization when membrane potential is below −60 mV) or depletion of [K]i. To investigate which one contributed most to the elevated RMP, we first clamped the [K]i to 140 mM, which is the steady [K]i under 1-Hz pacing, then stimulated the cell with 100-ms CL until the cell reached steady state. In this condition, the final RMP of the cell stayed below −78 mV, which means the rate dependence of the RMP disappeared. Moreover, this relatively low RMP also promoted the recovery of INa to such an extent that it eliminated the phenomenon of reduced INa at high pacing frequencies found in experiments by Syeda et al. (2016), and also led to more Na+ influx and [Na]i overload. On the other hand, slower late repolarization may also raise the RMP due to not fully repolarized membrane potential when next stimulus comes. The late repolarization was slowed in our model because of reduced IKr and Iss at high pacing rates. Below −60 mV, IKr and Iss were reduced by about 50 and 40% at CL = 100 ms compared with those at CL = 500 ms, respectively. According to these changes in current amplitude, we tried to augment the IKr and Iss to the level of slow pacing (CL = 500 ms) but the resultant RMP (−70 mV) was still elevated. Since it was the clamped [K]i but not faster late repolarization that eliminated the rate dependence of RMP, our simulation suggested the [K]i depletion was the predominant mechanism underlying the increase of RMP when the pacing rate was high.
Rate Dependence of Ion Channel Currents
For better explaining the mechanism underlying the rate dependence of AP, we investigated the role of each main channel current in the model. In Figures 12A,B, ion channel currents under 500, 250, and 100-ms CL were plotted against time or Vm. These three CLs were chosen since the rate dependence of currents was unobvious when the CL was larger than 500 ms. There is still no obvious variation of the morphology of the AP with the CL from 500 to 250 ms, whereas distinct changes in IK1, INCX, and INaK could be observed. Due to more APs evoked, much more Na+ entered the cell by INa, resulting in the increased [Na]i and then INaK. On the other hand, the elevated [Na]i suppressed the INCX and led to a rise of [Ca]i (as in Figure 9D). IK1 was decreased because its role in the late repolarization was partly taken over by the larger INaK.
Figure 12. Steady-state rate dependence of currents and corresponding I–V relationship. The unit of all currents is pA/pF. (A) Currents against time at various CLs. The traces of INa were slightly moved along time to avoid overlapping. (B) I–V relationship of currents at various CLs.
With a further decrease in CL (from 250 to 100 ms), distinct changes on AP morphology could be observed. The RMP was elevated by about 5 mV, resulting from the decrease of [K]i. The overshoot (OS) decreased from 30 to 16 mV, which was mainly due to the decrease of INa. INa decrease was also found by O’Hara et al. (2011) in their human ventricular cell model, but the reason of the decrease in our model was the elevated RMP leading to insufficient recovery of INa and limited Na channel availability, which was different to what was believed to be insufficient time for the INa to recover. For other currents, because of the smaller OS, ICaL and most currents for repolarization (Ito, IKur, IKr) attenuated in varying degrees. It is also worth noting that the decrease in [K+]i not only elevated the RMP, but also led to a reverse of IK1 due to the altered equilibrium potential of the potassium channel.
Role of K+ Currents in the Genesis of Mouse Atrial APs
To further validate the model, the functional roles of main K+ currents were investigated by ion channel blocking simulations. Results are compared with experimental data.
Role of Ito and IKur
Ito and IKur are both 4-AP-sensitive currents. Studies have demonstrated that relatively low concentrations (≤100 μM) of 4-AP predominantly blocks IKur while slightly affecting other types of K+ currents including Ito (Xu et al., 1999b). Since we cannot find experimental data for only Ito block, we consider the block of Ito and IKur by 4-AP in a combination manner.
Studies have shown that 4-AP prolongs APs of mouse atrial cell (Lomax et al., 2003a; Nygren et al., 2004; Trépanier-Boulay et al., 2004; Nakamura et al., 2010). Experimental data from 3 different studies are shown in Figures 13Ai,Bi,Ci. Since different amount of 4-AP and pacing frequencies were used in these experiments, to make meaningful comparison with our model, it is necessary to first determine how much potassium currents were affected by different concentration of 4-AP. In experiments, Trépanier-Boulay et al. (2004) found that Ito was only blocked 5% by 100-μM 4-AP, whereas Nakamura et al. (2010) found that 100-μM 4-AP may affect Ito to a greater extent, so they chose to add 50-μM 4-AP to avoid its effect on Ito. We also referred to the dose-dependent data of the effect of 4-AP on Ito and IKur in mouse ventricular cells (Xu et al., 1999b) and assumed that 4-AP should have similar effect in atrial cells. As a result, the percentages of Ito and IKur block were determined to: 16% Ito block and 70% IKur block under 50-μM 4-AP, or 30% Ito block and 100% IKur block under 100-μM 4-AP. The simulated APs and corresponding normalized APDs are shown in Figure 13. The block of Ito and IKur obviously prolongs APD20, APD50, and APD90. The simulation results were in good agreement with experimental data.
Figure 13. Role of Ito and IKur in the mouse atrial cell model. (Ai,Aii) Comparison of normalized APDs with experimental data (Ai) and simulated APs (Aii) before and after the application of 100-μM 4-AP. The cell was paced at 4 Hz. Experimental data was extracted from Trépanier-Boulay et al. (2004). (Bi,Bii) Comparison of normalized APDs with experimental data (Bi) and simulated APs (Bii) before and after the application of 50-μM 4-AP. The cell was paced at 2 Hz. Experimental data was extracted from Nakamura et al. (2010). (Ci,Cii) Comparison of normalized APDs with experimental data (Ci) and simulated APs (Cii) before and after the application of 100-μM 4-AP. The cell was paced at 1 Hz. Experimental data was extracted from Lomax et al. (2003a).
Role of IKr
The role of IKr in the present model is to modulate the late phase of APs in mouse atrial cells. Normalized APDs and APs with or without application of E-4031 are shown in Figure 14A. IKr is extremely sensitive to E-4031 and can be fully blocked by only 5 μM of E-4031 (Nakamura et al., 2010). In our simulation, when paced at 2 Hz, the APD20, APD50, and APD90 of the cell model were 4.9, 8.1, and 24.4 ms under control condition, respectively. Full block of IKr didn’t change the APD20 but prolonged APD50 by 0.3 ms and APD90 by 2.6 ms. Exposure to E-4031 only affected the repolarizing of APs at levels lower than −25 mV, and it prolonged APD90 much more than APD20 and APD50. These results were in agreement with experimental findings (Nakamura et al., 2010).
Figure 14. Role of IKr and Iss in the mouse atrial cell model. (Ai,Aii) Comparison of normalized APDs with experimental data (Ai) and simulated APs (Aii) before and after the application of 5-μM E-4031. The cell was paced at 2 Hz. Experimental data was extracted from Nakamura et al. (2010). (Bi,Bii) Comparison of normalized APDs with experimental data (Bi) and simulated APs (Bii) before and after the application of 5-mM TEA. The CL of the pacing protocol was 130 ms. Experimental data was extracted from Nygren et al. (2004).
Role of Iss
Iss is not affected by 4-AP lower than 1 mM, but it is sensitive to tetraethylammonium (TEA). Nygren et al. (2004) examined the effect of 5-mM TEA on APDs in isolated mouse atrial preparations. Their results are shown in Figure 14Bi. Since we cannot find dose-dependent measurements in mouse atrial cells, we assume that the sensitivity of Iss on TEA should be similar in the atrium and ventricle. It is reported that 25-mM TEA blocks 58% both Iss and IKur in mouse ventricular cells (Xu et al., 1999a), but their dose dependences on TEA are quite different. Increasing the concentration of TEA to 125 mM will block IKur completely but only block 61% Iss in mouse ventricular cells (Xu et al., 1999b), implying that the block of IKur by TEA has a steeper dose-dependent curve and lowering the dose of TEA should significantly decrease the extent of block on IKur. We did a linear regression on the log scaled dose-dependent curve of the amplitude of Iss and IKur, and estimated that 5-mM TEA should block about 55% of Iss and nearly no IKur. Our simulation results were in good agreement with experimental data (Figures 14Bi,Bii).
Discussion
Novel Achievements
A new biophysically detailed model for the mouse atrial myocyte with considerations of electrical heterogeneity between the left and right atrial cells has been developed, which is the first model representing the spatially heterogeneous mouse atrial APs in the atria. Based on the modification on the mouse ventricular cell model developed by Morotti et al. (2014), the mouse atrial cell model presented here describes cellular Ca2+ and Na+ handling and their regulation by CaMKII and PKA. Three main potassium currents (Ito, IKur, and IKr) were developed using Markov chain formalization and fitted to atrial experimental data, which allowed a better representation of atrial specific channel kinetics. Other membrane currents were also updated (INa, INaL, ICaL, Iss, IK1, IKACh, IKCa, IPMCA, INCX, and INaK) for atrial myocytes. All membrane currents were adjusted to physiological temperature either by fitting to data at 37°C or using temperature adjustment factors. Intracellular Ca2+ handling (RyRs, SERCA) and the CaMKII phosphorylation module were modified. Considering the heterogeneity between left and right atria, the atrial model was further extended to describe APs of mouse left atrial myocytes. The model was validated by its ability to reproduce the morphology and characteristics of APs and CaTs under various pacing frequencies, including APA, RP, APDs, etc. Characteristics of Ca2+ transients were comparable with those recorded from mouse atrial cells (see Figure 9). The role of main potassium currents in the genesis of the APs was investigated by simulation and validated with experimental data. The present model provides a basis for further simulating the conduction of excitation waves in the tissue of mouse atria in the future.
In the mouse atrial cell model constructed, we implemented Markov chain models instead of traditional HH models to model some of ion channels. Although the Markov chain model is more complex in the process of parameter fitting, and usually has more differential equations than the HH model and requires more simulation time, it helps to avoid limitations of the HH model in terms of model structure (Rudy and Silva, 2006). On the other hand, the Markov chain model is closely related to the conformation of ion channel proteins, which will be helpful in studying the change of channel function caused by drug action or gene defect in the future.
Comparison With Other Models
The morphology of the simulated AP in our atrial myocyte was relatively short, similar to the report that mouse AP is different from that in most other species (e.g., dog, pig or human) (Kaese and Verheule, 2012). The main factor responsible for short mouse AP is the relative magnitude of the currents present. Mouse myocytes have large outward K+ currents, which ensure rapid repolarization. In our model, APD90 was shorter than that in the ventricular model (Morotti et al., 2014), which was in agreement with experimental observation (Knollmann et al., 2007). In most mammalian species, the atrial AP is characterized by a shorter early repolarization phase, resulting in a triangular wave shape, compared with the ventricular AP with its usually pronounced plateau phase (Grand et al., 1990). Our atrial model is able to reproduce this triangular atrial AP wave shape. Unlike other mammals, the atrial mouse AP was longer during the early repolarization phase (APD25), compared with that in the ventricle. Our model suggests that it is likely to be the result of reduced current density of Ito in atria.
Ca2+ handling in mouse myocardium has a very rapid mechanism (Bers, 2001). Unlike other mammals with longer APDs, the percentage of Ca2+ released from intracellular stores (mainly stored in SR) to cytosol during a Ca2+ transient is much higher in mouse myocytes, and the contribution of Ca2+ influx via LTCC smaller. Due to the relatively short APD of mouse myocytes, the Ca2+ removal during diastole is also very fast to ensure a sufficient refilling before the next contraction. Compared with the mouse ventricular cell model (Bondarenko et al., 2004; Morotti et al., 2014), the atrial model presented here has smaller ICaL and INCX, indicating that the normal ECC coupling requires less Ca2+ influx and thus less Ca2+ needs to be extruded from the cell. Therefore, CaTs in mouse atrial myocytes are more dependent on the intracellular Ca2+ store. Due to rapid SR Ca2+ release and reuptake, CaTs reach the peak and decay faster than in ventricular myocytes. This is in line with the observation in rat experiments that the time to peak of CaTs is earlier and the decay is faster in atria than in ventricles (Escobar et al., 2004). The difference between atrial and ventricular cells is also reflected in the contribution of three main proteins (SERCA, NCX, and PMCA) to Ca2+ removal. SERCA contributes more to calcium clearance in the atrial calcium transient, which is also consistent with experimental measurements (Walden et al., 2009).
Recently, another mouse atrial model was developed by Asfaw et al. (2020). Both our model and the Asfaw model provide a set of mathematical equations to simulate the action potential and CaT of mouse atrial myocytes. But there are some distinctive differences between the two as listed below. (1) Our model considers IKACh, and three main K+ currents (Ito, IKur, and IKr), which are based on Markov chain formulations, which may be used as a tool for simulating functional impacts of gene mutations and drug interaction in future studies; (2) Our model also considers the electrical heterogeneity of left and right atrial myocytes (Figure 7), and thus extend the model to the left and right atrial cell model; (3) From the perspective of signal pathway, our model includes the effects of CaMKII on a variety of membrane currents, RyR, and SERCA, which plays an important role in the frequency dependent characteristics of the cell; (4) Though both models consider β-adrenergic receptor pathway, while the Asfaw model constructs β1- and β2-adrenergic receptors separately; and (5) In the subcellular structure, the Asfaw model considers caveolae and separates the SR into junctional SR and network SR, while our model does not. The simplification of our model on subcellular structure is described in more details in the limitation. From the perspective of model validation, our model was compared with more experimental data. Specifically, we collected the APDs of the left and right atrial myocytes, the characteristics of the CaT, the rate dependence of APD, RMP, dV/dtmax, and also the data of the effects of K+ current block on the AP. In addition, the normal heart rate of mice is around 600 beats per minute (BPM) (Donner et al., 2011). Our model is stable under 12.5-Hz stimulation (or 750 BPM), while the Asfaw model was only tested with pacing rates up to 6 Hz.
Due to the limited kinetic data on the CaMKII phosphorylation in mouse atrial cells, the CaMKII module is mainly based on the previously published mouse ventricular model (Morotti et al., 2014). Experimental studies have shown that CaMKII can phosphorylate the α1C and the β2a subunits of the LTCC (Grueter et al., 2006). Our model is able to reproduce the quick CaMKII-dependent LTCC facilitation found in many experiments (Valverde et al., 2005; Huke and Bers, 2007). There are also experimental studies showing that CaMKII can phosphorylate the RyR and enhance its sensitivity to Ca2+ (Wehrens et al., 2004), and the SR Ca2+ leak and open probability of RyR increase after being phosphorylated by CaMKII (Chelu et al., 2009). Our model is able to mimic these effects by implementing a CaMKII-dependent factor of SR sensitivity and a factor on SR Ca2+ leak. PLB resides in the SR membrane and negatively regulates the activity of SERCA. PLB phosphorylation will release its inhibition on SERCA, allowing SERCA to more freely uptake Ca2+ from cytosol. Our simulation shows slow kinetics of PLB phosphorylation (see the section of CaMKII-mediated phosphorylation in Supplementary Document 1 for more details), which is observed experimentally (Huke and Bers, 2007). Compared with the ventricular model, the inhibitory effect of PLB on SERCA in our atrial model is relatively small, because the expression of PLB in murine atrial tissue is only half of that in ventricular tissue (Freestone et al., 2000).
Intracellular Ionic Homeostasis
A valid cell model should be able to not only reproduce the action potential waveform but also ensure the ionic homeostasis at physiological heart rates. Under the prerequisite of keeping 4 types of ions dynamically changing (Na+, K+, Ca2+, Cl–), our model achieves ionic homeostasis under various pacing rates while all ion concentrations were in the physiological range. Dynamically changing ion concentrations enabled us to investigate the relationship between pacing frequency and cell behaviors. For example, we found in simulation that the [K]i depletion at high pacing rates led to a series of influences which was reported by Syeda et al. (2016), including elevated RMP, decreased INa, and reduced potassium currents. In all the main currents in our model, INaK and INCX play a key role of maintaining ionic homeostasis. INaK is the only current in charge of pumping Na+ out and K+ in, and INCX contributes 91% of the efflux of Ca2+ (the other part is extruded by PMCA). Along with increasing pacing frequency, more Na+ influx and K+ efflux lead to a higher [Na]i hence larger INaK until new balance is constructed between the influx and efflux of Na+.
Because of the conservation law, when the steady state is achieved, the net fluxes of each type of ions should be 0 in a whole cycle of an AP. As adopted by other cell models, the stimulus current in our simulation is treated as a K+ current to avoid long-term K+ depletion. Without doing this, the final balanced [K]i will decrease from 140 to 100 mM, leading to more than 5 mV elevation of the RMP, which is quite prominent. Our atrial cell model was carefully tuned to reproduce reasonable steady-state ion concentrations.
In simulations, a decrease of [K]i was observed during fast pacing. To test whether the rundown of the intracellular K may occur in models as an artifact of the pacing stimulus, we further checked if any change in the stimulus waveform (i.e., different amplitudes of the stimulus current) would affect the intracellular K+ concentration. In the model, to evoke an action potential a pacing stimulus current as an analog of K+ current with a square wave form, similar to that used in the human ventricular cell model developed by O’Hara et al. (2011) was implemented. In the control case, the amplitude and duration of the stimulus current implemented were 10 pA/pF and 4 ms respectively, which was about 1.5 times of the excitation threshold of the model (at 7 pA/pF with 4-ms duration). In the test conditions, six different amplitudes of stimulus ranging from 80 to 200% of the original stimulus amplitude were implemented at the pacing rate of 10 Hz, and the steady state values of intracellular K concentration were measured. It was shown that up to 200% of the stimulus amplitude used in the study had no obvious effects on the computed [K]i (<5%), suggesting that the rundown of [K]i was not attributed to the artifact of stimulus.
SERCA and FDAR
Frequency dependent acceleration of relaxation (FDAR) in murine hearts has been reported by many researchers (Kassiri et al., 2000; DeSantiago et al., 2002). This phenomenon is related to changes in myofilament properties (Chung et al., 2016) and in CaTs under high pacing rates. The changes of myofilament properties are beyond the aim of this paper, hereby we mainly discuss other possible roles contributing to frequency dependent acceleration of CaT decay (FDAD).
CaMKII plays an important role in the FDAD. The level of activated CaMKII is higher when the pacing is faster, and most targets phosphorylated by CaMKII directly influence the CaT. Among these targets, SERCA is predominant as it contributes more than 90% in the removal of [Ca]i, i.e., CaT decay, in murine cardiac myocytes (Walden et al., 2009). It is natural to guess that CaMKII amplifies SERCA under high pacing rates through some kind of mechanism. PLB may contribute to FDAD as it is a target phosphorylated by CaMKII and can regulate SERCA. There are reports showing that CaMKII phosphorylates PLB in a frequency-dependent manner (Hagemann et al., 2000; Valverde et al., 2005). But there are also studies suggesting temporal mismatch of the phosphorylation of PLB and the occurrence of FDAD (Valverde et al., 2005; Huke and Bers, 2007), which weakens the view that the frequency-dependent phosphorylation of PLB is the cause of FDAD. On the other hand, Picht et al. (2007) found that FDAD was caused by an increase in the Vmax of SERCA with presence of CaMKII and unaltered Km (Ca2+ affinity of SERCA), and was not related to the phosphorylation of PLB. From the perspective of modeling, if we want to guarantee the model has a 200-ms time constant of CaT decay under 0.5-Hz pacing, which is the value found in experiments (Li et al., 2005), the amplitude of SERCA has to be reduced to a degree that is too small at fast pacing rates to pump enough Ca2+ back in to the SR. In this case, the Ca2+ concentration in the cytosol and junctional area will be larger, leading to extra RyR Ca2+ releases during diastolic periods. Considering the exact mechanism of FDAD is still controversy, we chose to amplify the SERCA in an intuitive way that correlates the level of PLB phosphorylation with the Vmax of SERCA as described in the “Model Development” section. It is proved by the simulation results that our model shows normal FDAD and CaTs.
Limitations
Lack of Experimental Data
Although the mouse is a widely used animal model in the study of atrial arrhythmias, the quantification of the biophysical properties of mouse electrophysiology is still limited. In some instances, we need to use data obtained either at low temperature (lower than 36–37°C), from other regions of heart, or from different species, such as mouse ventricles, rats, canines, and rabbits. Most of the ion channels have to be fitted to experimental data acquired at non-physiological temperatures. Those experiments might be conducted at a temperature within the range of 20–37°C and by different groups. For Iss and IK1, especially, only data at room temperature were available. The model of ICaL was reformulated to fit with the mouse atrial I–V relationship, whereas the voltage-dependence of inactivation, the time constants of inactivation and recovery from inactivation were not validated due to the lack of available data. Also due to the lack of experimental data, the CaMKII and β-adrenergic modules were adapted from the parent mouse ventricular model. Another issue was the uncertainty about the accuracy of experimental data, which can vary by an order of magnitude, further limited the functionality of our model.
Limitation in the Cellular Structure
In ventricular myocytes, the well-developed transverse tubules (T-tubules) uniformly spread depolarization into the cell, resulting in a nearly synchronous SR Ca2+ release throughout the entire cell (Cheng et al., 1993). Atrial myocytes in small species, however, either lack a well-developed T-tubule network or have an irregular internal transverse-axial tubular system (Kirk et al., 2003) (Supplementary Figure S1B in Supplementary Document 1). As a result, the initial rise of the Ca2+ transient starts from the periphery of the cell and then the inner Ca2+ are released sequentially toward the cell center with decreasing amplitude (Trafford et al., 2013). Since the Ca2+ handling in our atrial model is based on the framework of ventricular myocytes (Shannon et al., 2004), SR Ca2+ releases only happen in the junctional area, and Ca2+ fluxes from the junctional area to SL and SL to cytosol are proportional with the Ca2+ gradients. This structure cannot mimic the Ca2+ sparks propagating in the cytosol, therefore underestimates the diffusion speed of Ca2+. Increasing the diffusion coefficient cannot be a solution because the Ca2+ sparks exist in only a short period, which cannot be mimicked by a larger but fixed diffusion coefficient through the whole AP.
Computational Cost
The Markov chain model as a widely used modeling paradigm, its advantages have been discussed extensively (Rudy and Silva, 2006; Fink and Noble, 2009). For example, Markov chain model can better mimic the kinetics of ion channels or proteins, and can be more easily extended to simulate the effects of drugs or gene mutations. However, Markov chain model is generally more complex than HH model, which consists of more differential equations and needs more computation. In the absence of extensive experimental data there are also some limitations in determining the model parameters for the complex system as Markov chain models, and the predictive capabilities of models fit only to experimental I–V curves are limited (Clerx et al., 2019; Whittaker et al., 2020).
Our first intention in constructing the presented model with Markov chain formulations for some ion channels was to describe the electrophysiological properties of mouse atrial myocytes as well as possible, and modulations to them by signaling pathways as detailed as possible, at a cost of non-optimized computational performance of the model. In fact, it costs about 0.8 s for 1 s simulation (pacing the model at 1 Hz) on a 2.8 GHz Core i7 CPU with one core. In multi-dimensional simulation, the computation can be paralleled to significantly reduce the time consumption. The computational time is nearly linear to the number of myocytes, which means that, it may take about 10 s for simulating 1 s electrical activity of a network of 100 myocytes with an 8-core CPU. Since modern CPUs are more and more powerful, we think the performance of our model is acceptable for most applications. If higher performance is needed, one could choose to replace some ion channels with HH models or remove the CaMKII or β-adrenergic pathways from the model.
Conclusion
A new model accounting for the heterogeneous electrical action potentials between the left and the right atrial cells has been developed. The model was derived from, and validated against experimental data obtained from mouse atrial cells. The developed model provides a novel basis for further study of possible mechanisms underlying atrial fibrillations.
Data Availability Statement
All datasets presented in this study are included in the article/Supplementary Material.
Author Contributions
HZ conceived the study. WS, SZ, and HZ designed most of the study, performed the simulations and analyses, and wrote most of the manuscript. WW contributed to the design of simulations, figure design, and manuscript writing. WS aided SZ with reuse of codes from his previous study. HZ, WS, and KW supervised the project. All authors contributed to the article and approved the submitted version.
Funding
This work was supported by the China Scholarship Council (CSC) and grants from EPSRC (United Kingdom) (EP/J00958X/1 and EP/I029826/1) and NSFC (61179009).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2020.00972/full#supplementary-material
References
Adelman, J. P., Maylie, J., and Sah, P. (2012). Small-conductance Ca2+-activated K+ channels: form and function. Annu. Rev. Physiol. 74, 245–269. doi: 10.1146/annurev-physiol-020911-153336
Asfaw, T. N., Tyan, L., Glukhov, A. V., and Bondarenko, V. E. (2020). A compartmentalized mathematical model of mouse atrial myocytes. Am. J. Physiol. 318, H485–H507. doi: 10.1152/ajpheart.00460.2019
Bao, Y., Willis, B. C., Frasier, C. R., Lopez-Santiago, L. F., Lin, X., Ramos-Mondragón, R., et al. (2016). Scn2b deletion in mice results in ventricular and atrial arrhythmias. Circ. Arrhythm. Electrophysiol. 9:e003923. doi: 10.1161/CIRCEP.116.003923
Bers, D. (2001). Excitation-Contraction Coupling and Cardiac Contractile Force. New York, NY: Springer Science & Business Media.
Bers, D. M. (2000). Calcium fluxes involved in control of cardiac myocyte contraction. Circ. Res. 87, 275–281. doi: 10.1161/01.RES.87.4.275
Blaustein, M. P., and Lederer, W. J. (1999). Sodium/calcium exchange: its physiological implications. Physiol. Rev. 79, 763–854. doi: 10.1152/physrev.1999.79.3.763
Bondarenko, V. E., Szigeti, G. P., Bett, G. C. L., Kim, S.-J., and Rasmusson, R. L. (2004). Computer model of action potential of mouse ventricular myocytes. Am. J. Physiol. 287, H1378–H1403. doi: 10.1152/ajpheart.00185.2003
Bou-Abboud, E., Li, H., and Nerbonne, J. M. (2000). Molecular diversity of the repolarizing voltage-gated K + currents in mouse atrial cells. J. Physiol. 529, 345–358. doi: 10.1111/j.1469-7793.2000.00345.x
Brouillette, J., Clark, R. B., Giles, W. R., and Fiset, C. (2004). Functional properties of K+ currents in adult mouse ventricular myocytes. J. Physiol. 559, 777–798. doi: 10.1113/jphysiol.2004.063446
Campbell, D. L., Rasmusson, R. L., Qu, Y., and Strauss, H. C. (1993). The calcium-independent transient outward potassium current in isolated ferret right ventricular myocytes. I. Basic characterization and kinetic analysis. J. Gen. Physiol. 101, 571–601. doi: 10.1085/jgp.101.4.571
Chelu, M. G., Sarma, S., Sood, S., Wang, S., van Oort, R. J., Skapura, D. G., et al. (2009). Calmodulin kinase II–mediated sarcoplasmic reticulum Ca2+ leak promotes atrial fibrillation in mice. J. Clin. Invest. 119, 1940–1951. doi: 10.1172/JCI37059
Cheng, H., Lederer, W. J., and Cannell, M. B. (1993). Calcium sparks: elementary events underlying excitation-contraction coupling in heart muscle. Science 262, 740–744. doi: 10.1126/science.8235594
Choi, E.-K., Chang, P.-C., Lee, Y.-S., Lin, S.-F., Zhu, W., Maruyama, M., et al. (2012). Triggered firing and atrial fibrillation in transgenic mice with selective atrial fibrosis induced by overexpression of TGF-β1. Circ. J. Off. J. Jpn. Circ. Soc. 76, 1354–1362. doi: 10.1253/circj.cj-11-1301
Chung, J.-H., Biesiadecki, B. J., Ziolo, M. T., Davis, J. P., and Janssen, P. M. L. (2016). Myofilament calcium sensitivity: role in regulation of in vivo cardiac contraction and relaxation. Front. Physiol. 7:562. doi: 10.3389/fphys.2016.00562
Clancy, C. E., and Rudy, Y. (2001). Cellular consequences of HERG mutations in the long QT syndrome: precursors to sudden cardiac death. Cardiovasc. Res. 50, 301–313. doi: 10.1016/s0008-6363(00)00293-5
Clerx, M., Beattie, K. A., Gavaghan, D. J., and Mirams, G. R. (2019). Four ways to fit an ion channel model. Biophys. J. 117, 2420–2437. doi: 10.1016/j.bpj.2019.08.001
Courtemanche, M., Ramirez, R. J., and Nattel, S. (1998). Ionic mechanisms underlying human atrial action potential properties: insights from a mathematical model. Am. J. Physiol. 275, H301–H321.
DeSantiago, J., Maier, L. S., and Bers, D. M. (2002). Frequency-dependent acceleration of relaxation in the heart depends on CaMKII, but not Phospholamban. J. Mol. Cell. Cardiol. 34, 975–984. doi: 10.1006/jmcc.2002.2034
Despa, S., Islam, M. A., Pogwizd, S. M., and Bers, D. M. (2002). Intracellular [Na+] and Na+ pump rate in rat and rabbit ventricular myocytes. J. Physiol. 539, 133–143. doi: 10.1113/jphysiol.2001.012940
DiFrancesco, D., and Noble, D. (1985). A model of cardiac electrical activity incorporating ionic pumps and concentration changes. Philos. Trans. R. Soc. Lond. B Biol. Sci. 307, 353–398. doi: 10.1098/rstb.1985.0001
Doi, T., Kuroda, S., Michikawa, T., and Kawato, M. (2005). Inositol 1,4,5-trisphosphate-dependent Ca2+ threshold dynamics detect spike timing in cerebellar purkinje cells. J. Neurosci. 25, 950–961. doi: 10.1523/JNEUROSCI.2727-04.2005
Donner, B. C., Schullenberg, M., Geduldig, N., Hüning, A., Mersmann, J., Zacharowski, K., et al. (2011). Functional role of TASK-1 in the heart: studies in TASK-1-deficient mice show prolonged cardiac repolarization and reduced heart rate variability. Basic Res. Cardiol. 106, 75–87. doi: 10.1007/s00395-010-0128-x
Eisner, D. A., Choi, H. S., Díaz, M. E., O’Neill, S. C., and Trafford, A. W. (2000). Integrative analysis of calcium cycling in cardiac muscle. Circ. Res. 87, 1087–1094. doi: 10.1161/01.RES.87.12.1087
Escobar, A. L., Ribeiro-Costa, R., Villalba-Galea, C., Zoghbi, M. E., Pérez, C. G., and Mejía-Alvarez, R. (2004). Developmental changes of intracellular Ca2+ transients in beating rat hearts. Am. J. Physiol. - Heart Circ. Physiol. 286, H971–H978. doi: 10.1152/ajpheart.00308.2003
Fink, M., and Noble, D. (2009). Markov models for ion channels: versatility versus identifiability and speed. Philos. Trans. R. Soc. Math. Phys. Eng. Sci. 367, 2161–2179. doi: 10.1098/rsta.2008.0301
Freestone, N. S., Ribaric, S., Scheuermann, M., Mauser, U., Paul, M., and Vetter, R. (2000). Differential lusitropic responsiveness to β-adrenergic stimulation in rat atrial and ventricular cardiac myocytes. Pflüg. Arch. 441, 78–87. doi: 10.1007/s004240000397
Glitsch, H. G. (2001). Electrophysiology of the sodium-potassium-ATPase in cardiac cells. Physiol. Rev. 81, 1791–1826. doi: 10.1152/physrev.2001.81.4.1791
Glukhov, A. V., Balycheva, M., Sanchez-Alonso, J. L., Ilkan, Z., Alvarez-Laviada, A., Bhogal, N., et al. (2015). Direct evidence for microdomain-specific localization and remodeling of functional L-type calcium channels in rat and human atrial myocytes. Circulation 132, 2372–2384. doi: 10.1161/CIRCULATIONAHA.115.018131
Grand, B. L., Heuzey, J.-Y. L., Perier, P., Peronneau, P., Lavergne, T., Hatem, S., et al. (1990). Cellular electrophysiological effects of flecainide on human atrial fibres. Cardiovasc. Res. 24, 232–238. doi: 10.1093/cvr/24.3.232
Grueter, C. E., Abiria, S. A., Dzhura, I., Wu, Y., Ham, A.-J. L., Mohler, P. J., et al. (2006). RETRACTED: L-type Ca2+ channel facilitation mediated by phosphorylation of the β subunit by CaMKII. Mol. Cell 23, 641–650. doi: 10.1016/j.molcel.2006.07.006
Gurney, A., and Manoury, B. (2008). Two-pore potassium channels in the cardiovascular system. Eur. Biophys. J. 38:305. doi: 10.1007/s00249-008-0326-8
Hagemann, D., Kuschel, M., Kuramochi, T., Zhu, W., Cheng, H., and Xiao, R.-P. (2000). Frequency-encoding Thr17 phospholamban phosphorylation is independent of Ser16 phosphorylation in cardiac myocytes. J. Biol. Chem. 275, 22532–22536. doi: 10.1074/jbc.C000253200
Hancock, J. M., Weatherall, K. L., Choisy, S. C., James, A. F., Hancox, J. C., and Marrion, N. V. (2015). Selective activation of heteromeric SK channels contributes to action potential repolarization in mouse atrial myocytes. Heart Rhythm 12, 1003–1015. doi: 10.1016/j.hrthm.2015.01.027
Hofmann, F., Lacinová, L., and Klugbauer, N. (1999). “Voltage-dependent calcium channels: from structure to function,” in Reviews of Physiology, Biochemistry and Pharmacology, Volume 139 Reviews of Physiology, Biochemistry and Pharmacology, (Berlin: Springer), 33–87. doi: 10.1007/BFb0033648
Hu, Y., Jones, S. V. P., and Dillmann, W. H. (2005). Effects of hyperthyroidism on delayed rectifier K+ currents in left and right murine atria. Am. J. Physiol. 289, H1448–H1455. doi: 10.1152/ajpheart.00828.2004
Hua, R., Adamczyk, A., Robbins, C., Ray, G., and Rose, R. A. (2012). Distinct patterns of constitutive phosphodiesterase activity in mouse sinoatrial node and atrial myocardium. PLoS One 7:e47652. doi: 10.1371/journal.pone.0047652
Huang, C. L.-H. (2016). Murine electrophysiological models of cardiac arrhythmogenesis. Physiol. Rev. 97, 283–409. doi: 10.1152/physrev.00007.2016
Huke, S., and Bers, D. M. (2007). Temporal dissociation of frequency-dependent acceleration of relaxation and protein phosphorylation by CaMKII. J. Mol. Cell. Cardiol. 42, 590–599. doi: 10.1016/j.yjmcc.2006.12.007
Kaese, S., and Verheule, S. (2012). Cardiac electrophysiology in mice: a matter of size. Front. Physiol. 3:345. doi: 10.3389/fphys.2012.00345
Kassiri, Z., Myers, R., Kaprielian, R., Banijamali, H., and Backx, P. (2000). Rate-dependent changes of twitch force duration in rat cardiac trabeculae: a property of the contractile system. J. Physiol. 524, 221–231. doi: 10.1111/j.1469-7793.2000.t01-3-00221.x
Kharche, S., Yu, J., Lei, M., and Zhang, H. (2011). A mathematical model of action potentials of mouse sinoatrial node cells with molecular bases. Am. J. Physiol. 301, H945–H963. doi: 10.1152/ajpheart.00143.2010
Kirk, M. M., Izu, L. T., Chen-Izu, Y., McCulle, S. L., Wier, W. G., Balke, C. W., et al. (2003). Role of the transverse-axial tubule system in generating calcium sparks and calcium transients in rat atrial myocytes. J. Physiol. 547, 441–451. doi: 10.1113/jphysiol.2002.034355
Knollmann, B. C., Schober, T., Petersen, A. O., Sirenko, S. G., and Franz, M. R. (2007). Action potential characterization in intact mouse heart: steady-state cycle length dependence and electrical restitution. Am. J. Physiol. 292, H614–H621. doi: 10.1152/ajpheart.01085.2005
Kovoor, P., Wickman, K., Maguire, C. T., Pu, W., Gehrmann, J., Berul, C. I., et al. (2001). Evaluation of the role of IKAChin atrial fibrillation using a mouse knockout model. J. Am. Coll. Cardiol. 37, 2136–2143. doi: 10.1016/S0735-1097(01)01304-3
Lemoine, M. D., Duverger, J. E., Naud, P., Chartier, D., Qi, X. Y., Comtois, P., et al. (2011). Arrhythmogenic left atrial cellular electrophysiology in a murine genetic long QT syndrome model. Cardiovasc. Res. 92, 67–74. doi: 10.1093/cvr/cvr166
Li, G.-R., Feng, J., Yue, L., Carrier, M., and Nattel, S. (1996). Evidence for two components of delayed rectifier K+ current in human ventricular myocytes. Circ. Res. 78, 689–696. doi: 10.1161/01.RES.78.4.689
Li, L., Niederer, S. A., Idigo, W., Zhang, Y. H., Swietach, P., Casadei, B., et al. (2010). A mathematical model of the murine ventricular myocyte: a data-driven biophysically based approach applied to mice overexpressing the canine NCX isoform. Am. J. Physiol. 299, H1045–H1063. doi: 10.1152/ajpheart.00219.2010
Li, X., Zima, A. V., Sheikh, F., Blatter, L. A., and Chen, J. (2005). Endothelin-1–induced arrhythmogenic Ca2+ signaling is abolished in atrial myocytes of inositol-1,4,5-trisphosphate(IP3)–receptor type 2–deficient mice. Circ. Res. 96, 1274–1281. doi: 10.1161/01.RES.0000172556.05576.4c
Lomax, A. E., Kondo, C. S., and Giles, W. R. (2003a). Comparison of time- and voltage-dependent K+ currents in myocytes from left and right atria of adult mice. Am. J. Physiol. 285, H1837–H1848. doi: 10.1152/ajpheart.00386.2003
Lomax, A. E., Rose, R. A., and Giles, W. R. (2003b). Electrophysiological evidence for a gradient of G protein-gated K+ current in adult mouse atria. Br. J. Pharmacol. 140, 576–584. doi: 10.1038/sj.bjp.0705474
Lukyanenko, V., Wiesner, T. F., and Györke, S. (1998). Termination of Ca2+ release during Ca2+ sparks in rat ventricular myocytes. J. Physiol. 507, 667–677. doi: 10.1111/j.1469-7793.1998.667bs.x
Mahajan, A., Shiferaw, Y., Sato, D., Baher, A., Olcese, R., Xie, L.-H., et al. (2008). A rabbit ventricular action potential model replicating cardiac dynamics at rapid heart rates. Biophys. J. 94, 392–410. doi: 10.1529/biophysj.106.98160
Mancarella, S., Yue, Y., Karnabi, E., Qu, Y., El-Sherif, N., and Boutjdir, M. (2008). Impaired Ca2+ homeostasis is associated with atrial fibrillation in the α1D L-type Ca2+ channel KO mouse. Am. J. Physiol. 295, H2017–H2024. doi: 10.1152/ajpheart.00537.2008
Mangoni, M. E., Traboulsie, A., Leoni, A.-L., Couette, B., Marger, L., Quang, K. L., et al. (2006). Bradycardia and slowing of the atrioventricular conduction in mice lacking CaV3.1/α1G T-type calcium channels. Circ. Res. 98, 1422–1430. doi: 10.1161/01.RES.0000225862.14314.49
Martin, C. A., Matthews, G. D. K., and Huang, C. L.-H. (2012). Sudden cardiac death and inherited channelopathy: the basic electrophysiology of the myocyte and myocardium in ion channel disease. Heart 98, 536–543. doi: 10.1136/heartjnl-2011-300953
Mewe, M., Mauerhöfer, M., Wulfsen, I., Szlachta, K., Zhou, X.-B., Schwarz, J. R., et al. (2010). Modulation of cardiac ERG1 K(+) channels by cGMP signaling. J. Mol. Cell. Cardiol. 49, 48–57. doi: 10.1016/j.yjmcc.2010.02.015
Morotti, S., Edwards, A. G., McCulloch, A. D., Bers, D. M., and Grandi, E. (2014). A novel computational model of mouse myocyte electrophysiology to assess the synergy between Na+ loading and CaMKII. J. Physiol. 592, 1181–1197. doi: 10.1113/jphysiol.2013.266676
Nakamura, H., Ding, W.-G., Sanada, M., Maeda, K., Kawai, H., Maegawa, H., et al. (2010). Presence and functional role of the rapidly activating delayed rectifier K+ current in left and right atria of adult mice. Eur. J. Pharmacol. 649, 14–22. doi: 10.1016/j.ejphar.2010.08.025
Nerbonne, J. M. (2000). Molecular basis of functional voltage-gated K+ channel diversity in the mammalian myocardium. J. Physiol. 525, 285–298. doi: 10.1111/j.1469-7793.2000.t01-1-00285.x
Nygren, A., Lomax, A. E., and Giles, W. R. (2004). Heterogeneity of action potential durations in isolated mouse left and right atria recorded using voltage-sensitive dye mapping. Am. J. Physiol. 287, H2634–H2643. doi: 10.1152/ajpheart.00380.2004
O’Hara, T., Virág, L., Varró, A., and Rudy, Y. (2011). Simulation of the undiseased human cardiac ventricular action potential: model formulation and experimental validation. PLoS Comput. Biol. 7:e1002061. doi: 10.1371/journal.pcbi.1002061
Pandit, S. V., Zlochiver, S., Filgueiras-Rama, D., Mironov, S., Yamazaki, M., Ennis, S. R., et al. (2011). Targeting atrioventricular differences in ion channel properties for terminating acute atrial fibrillation in pigs. Cardiovasc. Res. 89, 843–851. doi: 10.1093/cvr/cvq359
Picht, E., DeSantiago, J., Huke, S., Kaetzel, M. A., Dedman, J. R., and Bers, D. M. (2007). CaMKII inhibition targeted to the sarcoplasmic reticulum inhibits frequency dependent acceleration of relaxation and Ca2+ current facilitation. J. Mol. Cell. Cardiol. 42, 196–205. doi: 10.1016/j.yjmcc.2006.09.007
Qin, M., Huang, H., Wang, T., Hu, H., Liu, Y., Gu, Y., et al. (2012). Atrial tachyarrhythmia in Rgs5-null mice. PLoS One 7:e46856. doi: 10.1371/journal.pone.0046856
Ramos-Franco, J., Fill, M., and Mignery, G. A. (1998). Isoform-specific function of single inositol 1,4,5-trisphosphate receptor channels. Biophys. J. 75, 834–839. doi: 10.1016/S0006-3495(98)77572-1
Rose, W. C., Balke, C. W., Wier, W. G., and Marban, E. (1992). Macroscopic and unitary properties of physiological ion flux through L-type Ca2+ channels in guinea-pig heart cells. J. Physiol. 456, 267–284. doi: 10.1113/jphysiol.1992.sp019336
Rudy, Y., and Silva, J. R. (2006). Computational biology in the study of cardiac ion channels and cell electrophysiology. Q. Rev. Biophys. 39, 57–116. doi: 10.1017/s0033583506004227
Sanguinetti, M. C., and Jurkiewicz, N. K. (1990). Two components of cardiac delayed rectifier K+ current. Differential sensitivity to block by class III antiarrhythmic agents. J. Gen. Physiol. 96, 195–215. doi: 10.1085/jgp.96.1.195
Sarmast, F., Kolli, A., Zaitsev, A., Parisian, K., Dhamoon, A. S., Guha, P. K., et al. (2003). Cholinergic atrial fibrillation: IK, ACh gradients determine unequal left/right atrial frequencies and rotor dynamics. Cardiovasc. Res. 59, 863–873. doi: 10.1016/S0008-6363(03)00540-6
Schmidt, T. A., Svendsen, J. H., Haunsø, S., and Kjeldsen, K. (1990). Quantification of the total Na, K-ATPase concentration in atria and ventricles from mammalian species by measuring 3 H-ouabain binding to intact myocardial samples. Stability to short term ischemia reperfusion. Basic Res. Cardiol. 85, 411–427. doi: 10.1007/bf01907133
Shanmugam, M., Gao, S., Hong, C., Fefelova, N., Nowycky, M. C., Xie, L.-H., et al. (2011). Ablation of phospholamban and sarcolipin results in cardiac hypertrophy and decreased cardiac contractility. Cardiovasc. Res. 89, 353–361. doi: 10.1093/cvr/cvq294
Shannon, T. R., Wang, F., Puglisi, J., Weber, C., and Bers, D. M. (2004). A mathematical treatment of integrated Ca dynamics within the ventricular myocyte. Biophys. J. 87, 3351–3371. doi: 10.1529/biophysj.104.047449
Syeda, F., Holmes, A. P., Yu, T. Y., Tull, S., Kuhlmann, S. M., Pavlovic, D., et al. (2016). PITX2 modulates atrial membrane potential and the antiarrhythmic effects of sodium-channel blockers. J. Am. Coll. Cardiol. 68, 1881–1894. doi: 10.1016/j.jacc.2016.07.766
Trafford, A. W., Clarke, J. D., Richards, M. A., Eisner, D. A., and Dibb, K. M. (2013). Calcium signalling microdomains and the t-tubular system in atrial mycoytes: potential roles in cardiac disease and arrhythmias. Cardiovasc. Res. 98, 192–203. doi: 10.1093/cvr/cvt018
Trépanier-Boulay, V., Lupien, M.-A., St-Michel, C., and Fiset, C. (2004). Postnatal development of atrial repolarization in the mouse. Cardiovasc. Res. 64, 84–93. doi: 10.1016/j.cardiores.2004.06.002
Valverde, C. A., Mundiña-Weilenmann, C., Said, M., Ferrero, P., Vittone, L., Salas, M., et al. (2005). Frequency-dependent acceleration of relaxation in mammalian heart: a property not relying on phospholamban and SERCA2a phosphorylation. J. Physiol. 562, 801–813. doi: 10.1113/jphysiol.2004.075432
Verheule, S., Sato, T., Everett, T., Engle, S. K., Otten, D., der Lohe, M. R., et al. (2004). Increased vulnerability to atrial fibrillation in transgenic mice with selective atrial fibrosis caused by overexpression of TGF-β1. Circ. Res. 94, 1458–1465. doi: 10.1161/01.RES.0000129579.59664.9d
Voigt, N., Abu-Taha, I., Heijman, J., and Dobrev, D. (2014). Constitutive activity of the acetylcholine-activated potassium current IK, ACh in cardiomyocytes. Adv. Pharmacol. 70, 393–409. doi: 10.1016/B978-0-12-417197-8.00013-4
Wagner, S., Dybkova, N., Rasenack, E. C. L., Jacobshagen, C., Fabritz, L., Kirchhof, P., et al. (2006). Ca2+/calmodulin-dependent protein kinase II regulates cardiac Na+ channels. J. Clin. Invest. 116, 3127–3138. doi: 10.1172/JCI26620
Wagner, S., Hacker, E., Grandi, E., Weber, S. L., Dybkova, N., Sossalla, S., et al. (2009). Ca/Calmodulin kinase II differentially modulates potassium currents. Circ. Arrhythm. Electrophysiol. 2, 285–294. doi: 10.1161/CIRCEP.108.842799
Wakimoto, H., Maguire, C. T., Kovoor, P., Hammer, P. E., Gehrmann, J., Triedman, J. K., et al. (2001). Induction of atrial tachycardia and fibrillation in the mouse heart. Cardiovasc. Res. 50, 463–473. doi: 10.1016/S0008-6363(01)00264-4
Walden, A., Dibb, K., and Trafford, A. (2009). Differences in intracellular calcium homeostasis between atrial and ventricular myocytes. J. Mol. Cell. Cardiol. 46, 463–473. doi: 10.1016/j.yjmcc.2008.11.003
Wang, G., Schmied, R., Ebner, F., and Korth, M. (1993). Intracellular sodium activity and its regulation in guinea-pig atrial myocardium. J. Physiol. 465, 73–84. doi: 10.1113/jphysiol.1993.sp019667
Wang, J., Schwinger, R., Frank, K., Müller-Ehmsen, J., Martin-Vasallo, P., Pressley, T. A., et al. (1996). Regional expression of sodium pump subunits isoforms and Na+-Ca++ exchanger in the human heart. J. Clin. Invest. 98:1650. doi: 10.1172/jci118960
Wang, Z., Fermini, B., and Nattel, S. (1994). Rapid and slow components of delayed rectifier current in human atrial myocytes. Cardiovasc. Res. 28, 1540–1546. doi: 10.1093/cvr/28.10.1540
Wehrens, X. H. T., Lehnart, S. E., Reiken, S. R., and Marks, A. R. (2004). Ca2+/Calmodulin-dependent protein kinase ii phosphorylation regulates the cardiac ryanodine receptor. Circ. Res. 94, e61–e70. doi: 10.1161/01.RES.0000125626.33738.E2
Wettwer, E., Hála, O., Christ, T., Heubach, J. F., Dobrev, D., Knaut, M., et al. (2004). Role of IKur in controlling action potential shape and contractility in the human atrium influence of chronic atrial fibrillation. Circulation 110, 2299–2306. doi: 10.1161/01.CIR.0000145155.60288.71
Whittaker, D. G., Clerx, M., Lei, C. L., Christini, D. J., and Mirams, G. R. (2020). Calibration of ionic and cellular cardiac electrophysiology models. Wiley Interdiscip. Rev. Syst. Biol. Med. 12:e1482.
Wu, Y., and Anderson, M. E. (2002). Reduced repolarization reserve in ventricular myocytes from female mice. Cardiovasc. Res. 53, 763–769. doi: 10.1016/S0008-6363(01)00387-X
Xiao, Y.-F., TenBroek, E. M., Wilhelm, J. J., Iaizzo, P. A., and Sigg, D. C. (2006). Electrophysiological characterization of murine HL-5 atrial cardiomyocytes. Am. J. Physiol. 291, C407–C416. doi: 10.1152/ajpcell.00020.2006
Xie, L.-H., Shanmugam, M., Park, J. Y., Zhao, Z., Wen, H., Tian, B., et al. (2012). Ablation of sarcolipin results in atrial remodeling. Am. J. Physiol. 302, C1762–C1771. doi: 10.1152/ajpcell.00425.2011
Xu, H., Barry, D. M., Li, H., Brunet, S., Guo, W., and Nerbonne, J. M. (1999a). Attenuation of the slow component of delayed rectification, action potential prolongation, and triggered activity in mice expressing a dominant-negative Kv2 α subunit. Circ. Res. 85, 623–633. doi: 10.1161/01.RES.85.7.623
Xu, H., Guo, W., and Nerbonne, J. M. (1999b). Four kinetically distinct depolarization-activated K+ currents in adult mouse ventricular myocytes. J. Gen. Physiol. 113, 661–678. doi: 10.1085/jgp.113.5.661
Xu, H., Li, H., and Nerbonne, J. M. (1999c). Elimination of the transient outward current and action potential prolongation in mouse atrial myocytes expressing a dominant negative Kv4 α subunit. J. Physiol. 519, 11–21. doi: 10.1111/j.1469-7793.1999.0011o.x
Yang, J. H., and Saucerman, J. J. (2011). Computational models reduce complexity and accelerate insight into cardiac signaling networks. Circ. Res. 108, 85–97. doi: 10.1161/circresaha.110.223602
Yang, J. H., and Saucerman, J. J. (2012). Phospholemman is a negative feed-forward regulator of Ca2+ in β-adrenergic signaling, accelerating β-adrenergic inotropy. J. Mol. Cell. Cardiol. 52, 1048–1055. doi: 10.1016/j.yjmcc.2011.12.015
Zhang, H., Holden, A. V., Noble, D., and Boyett, M. R. (2002). Analysis of the chronotropic effect of acetylcholine on sinoatrial node cells. J. Cardiovasc. Electrophysiol. 13, 465–474. doi: 10.1046/j.1540-8167.2002.00465.x
Keywords: mouse atrial, action potential, mathematical model, action potential, Ca2+ handling, atrial heterogeneity, signaling pathway
Citation: Zhang H, Zhang S, Wang W, Wang K and Shen W (2020) A Mathematical Model of the Mouse Atrial Myocyte With Inter-Atrial Electrophysiological Heterogeneity. Front. Physiol. 11:972. doi: 10.3389/fphys.2020.00972
Received: 16 April 2020; Accepted: 16 July 2020;
Published: 06 August 2020.
Edited by:
T. Alexander Quinn, Dalhousie University, CanadaReviewed by:
Jordi Heijman, Maastricht University, NetherlandsThomas Hund, The Ohio State University, United States
Copyright © 2020 Zhang, Zhang, Wang, Wang and Shen. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Henggui Zhang, henggui.zhang@manchester.ac.uk; Weijian Shen, weijian@shen.in
†These authors have contributed equally to this work