HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) All Versions of this Article: 287/3/H1378    most recent 00185.2003v1 Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in ISI Web of Science Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via ISI Web of Science (33) Citing Articles via Google Scholar Google Scholar Articles by Bondarenko, V. E. Articles by Rasmusson, R. L. Search for Related Content PubMed PubMed Citation Articles by Bondarenko, V. E. Articles by Rasmusson, R. L.

Computer model of action potential of mouse ventricular myocytes

¹Department of Physiology and Biophysics, School of Medicine and Biomedical Sciences, University at Buffalo, State University of New York, Buffalo, New York 14214-3078; and ²Cell Biology and Molecular Medicine Cardiovascular Research Institute, University of Medicine and Dentistry of New Jersey-New Jersey Medical School, Newark, New Jersey 07103

Submitted 28 February 2003 ; accepted in final form 11 May 2004

	ABSTRACT

TOP ABSTRACT METHODS RESULTS DISCUSSION APPENDIX: MOUSE MODEL SUMMARY GRANTS REFERENCES

 
We have developed a mathematical model of the mouse ventricular myocyte action potential (AP) from voltage-clamp data of the underlying currents and Ca²⁺ transients. Wherever possible, we used Markov models to represent the molecular structure and function of ion channels. The model includes detailed intracellular Ca²⁺ dynamics, with simulations of localized events such as sarcoplasmic Ca²⁺ release into a small intracellular volume bounded by the sarcolemma and sarcoplasmic reticulum. Transporter-mediated Ca²⁺ fluxes from the bulk cytosol are closely matched to the experimentally reported values and predict stimulation rate-dependent changes in Ca²⁺ transients. Our model reproduces the properties of cardiac myocytes from two different regions of the heart: the apex and the septum. The septum has a relatively prolonged AP, which reflects a relatively small contribution from the rapid transient outward K⁺ current in the septum. The attribution of putative molecular bases for several of the component currents enables our mouse model to be used to simulate the behavior of genetically modified transgenic mice.

cardiac myocytes; computer modeling

The pioneering model of DiFrancesco and Noble (27) was the first to incorporate active transport and secondary active transport mechanisms. These were included in response to the growing experimental evidence that these electrogenic processes were important in repolarization. Development of the whole cell patch-clamp technique gave rise to a number of animal- and region-specific models, based on detailed descriptions of ionic currents from isolated single-myocyte experiments (41, 55, 62, 63, 76, 77, 91, 98). These models began to address the basis for the diverse AP behavior observed in different species and regions of the heart. Models of bullfrog atrial and cardiac pacemaker cells (76, 77) replicated the lack of a significant intracellular Ca²⁺ release contribution in these cells. Similarly, models of the guinea pig ventricular myocyte (62, 63, 98) reflected the absence of the transient outward current in these cells. Lindblad et al. (55) developed a model of the rabbit atrial cell with prominent transient outward (I_to) and rapid delayed rectifier (I_Kr) currents. Human atrial models were developed by Nygren et al. (72) and Courtemanche et al. (25), based partially on measurements from atrial appendage tissue obtained from bypass surgery. I_to had a prominent role in these AP models, reflecting the experimentally obtained data. Recently, a rat ventricular myocyte model was published by Pandit et al. (73). The rat AP has a short AP duration (APD) and relatively dominant transient outward current.

APs of cardiac cells from different animals and different heart regions have different shapes, durations, and sets of ionic currents. These variations reflect corresponding differences in the molecular basis of repolarization and hence pharmacological responses (32, 37, 56). For example, APD at 50% of repolarization (APD₅₀) in guinea pig and canine ventricular myocytes is 300 ms (41, 62, 63, 91), whereas APD₅₀ for rat ventricular myocytes from the epicardial region is only 13 ms (73). The guinea pig ventricular myocyte AP has an elevated plateau (41, 62, 63, 91), whereas rat and mouse ventricular APs are relatively short, with a "triangular" shape (11, 93). The human atrial AP has a relatively sharp, short peak followed by a slow repolarization phase, which takes up to 250 ms (33, 72). In general, APD seems to be correlated to the size of the heart or the heart rate for a particular species, with smaller, faster-beating hearts having shorter APDs than larger, slower-beating hearts. The shape and underlying currents for the various APs have considerable variability. This variability is less related to heart size than it is to the species and region of the heart. Presumably, the restitution patterns, pharmacology, and arrhythmogenic potential of each waveform reflects the different molecular bases and timing of underlying currents. Understanding the timing of these currents under different situations can be critical. For example, one of the arrhythmogenic defects in human ether-à-go-go-related gene (HERG) (which encodes the rapid delayed rectifier current I_Kr) may be spontaneously triggered arrhythmias resulting from a complex interaction of I_Kr with the L-type Ca²⁺ current and subcellular Ca²⁺-handling mechanisms during repolarization (67). Clearly, understanding abnormalities in cardiac repolarization requires an understanding of the functioning of the cell as an integrated system. The many complex potential interactions of transmembrane currents and intracellular ions (particularly Ca²⁺) leave computer modeling as the only method currently available to synthesize and understand these multiple nonlinear interactions.

In general, most models simulate sarcolemmal currents, pumps, and exchangers and combine them with intracellular ionic homeostasis to reproduce an AP. One limitation of many older models is their description of intracellular Ca²⁺ dynamics (55, 62, 63, 76, 77, 98). These older models did not take into account the importance of spatial localization in Ca²⁺ cellular dynamics, which results in generation of relatively high-amplitude Ca²⁺ sparks (see, for example, Ref. 11) and can only roughly approximate the relationship between Ca²⁺ release and Ca²⁺ channel inactivation. More recent models (12, 41, 91) have begun to address this issue. Our model builds on the strengths of these previous models and was developed to reproduce a mouse myocyte with a characteristically short AP.

There is a great diversity of anatomic and electrophysiological features in myocytes from different species, e.g., transverse tubules, Ca²⁺-handling systems, as well as the specific properties of the major ionic currents. Wherever possible, we developed our model of the mouse ventricular myocytes with experimental data from adult mouse ventricular cardiac cells. This has been possible largely because of the increase in detailed examination of currents in various transgenic mice.

Although a majority of the most important currents and parameters are well constrained by this large body of data, there is still uncertainty in the exact kinetic parameters of some current systems. Where there was uncertainty or an absence of mouse data, we used data from less closely related species and heterologously expressed clones and adjusted them to match related mouse data.

Our mouse model is specifically designed to make use of the information from recent advances in molecular biology. The mouse is the most widely used animal in genetic research. There is now an abundance of transgenic mice with genetic defects associated with human diseases. Where possible, the component currents of the model have been assigned putative molecular bases, so the model can be used to test predictions about the consequences of transgenic experiments with relative ease. Our model mouse myocyte has a diversity of K⁺ channels, not all of which may be present in the same cell, depending on the region from which the cells were isolated. The model has seven distinct K⁺ currents: a rapid transient outward K⁺ current (I_Kto,f), a slow transient outward K⁺ current (I_Kto,s), a time-independent K⁺ current (I_K1), an ultrarapidly activating delayed rectifier K⁺ current (I_Kur), a noninactivating steady-state K⁺ current (I_Kss), a rapid delayed rectifier K⁺ current (I_Kr), and a slow delayed rectifier K⁺ current (I_Ks).

Molecular biology has provided considerable information on the structure and function of ion channels, e.g., mechanistic information concerning ion channel gating and its modification by either genetic disease or pharmacological intervention. To be able to assess the consequences of these modifications, we used Markov models for several important currents: the Na⁺ current (I_Na), the L-type Ca²⁺ current (I_CaL), and I_Kr. Many genetic manipulations and diseases alter the kinetic properties of a single ion channel and can be related to changes in specific rate constants in the Markov models for these channels. Our mouse model can therefore be used to predict the kinetic consequences of many transgenic or pharmacological modifications on channel properties.

The main goal of this paper was to develop a computer model of the mouse ventricular AP that will provide insight into the ionic basis of AP behavior. We modeled the behavior of myocytes from the apex and the septum regions of the heart (36, 37, 69, 92) to demonstrate that the different regional expression of I_Kto,f, I_Kto,s, I_Kur, and I_Kss can account for regional differences in myocyte repolarization in the mouse heart. The model contains transmembrane pumps, currents, and exchangers and a detailed intracellular Ca²⁺-handling system based on data from voltage-clamp experiments on individual currents (e.g., steady-state inactivation, recovery from inactivation). We developed Markov models for the fast Na⁺ current (I_Na), the L-type Ca²⁺ current (I_CaL) and I_Kr. These formalisms allow for future investigation of the cellular mechanisms of arrhythmia generation due to genetic mutations or state-dependent binding that alter functional electrical properties such as restitution.

	METHODS

TOP ABSTRACT METHODS RESULTS DISCUSSION APPENDIX: MOUSE MODEL SUMMARY GRANTS REFERENCES

 
All animal procedures conformed to the "Guiding Principles for Research Involving Animals and Human Beings" of the American Physiological Society.

Preparation of Left Ventricular Myocytes

Mouse cardiac myocytes were prepared as described previously (48). Briefly, the heart was rapidly excised and submerged in Ca²⁺-free Tyrode solution containing (in mM) 140 NaCl, 5.4 KCl, 1 MgCl₂, 0.33 NaH₂PO₄, 10 glucose, and 5 HEPES adjusted to pH 7.4. The aorta was cannulated with a blunt-tip needle (20 gauge) placed on a perfusion apparatus. The heart was retrogradely perfused for 3 min with Tyrode solution and then perfused for 18 min with Tyrode solution with 2% calf serum (Sigma-Aldrich) and 75 U/ml each of collagenase type I and type II (Worthington). All solutions were continuously bubbled with 95% O₂-5% CO₂ at 37°C. Isolated myocytes were stored at room temperature in low-Cl^–, high-K⁺ solution (in mM: 70 KOH, 147 L-glutamic acid, 40 KCl, 20 taurine, 20 KH₂PO₄, 10 glucose, 10 HEPES, and 0.5 EGTA, adjusted to pH 7.3 with Tris base) before experiments. Calcium-tolerant, rod-shaped ventricular myocytes with clear striations were randomly selected for electrophysiological studies.

Cellular Electrophysiological Studies

Cell-attached patch-clamp macroscopic currents were recorded with an Axopatch 1D amplifier (Axon Instruments). A DigiData 1200 (Axon Instruments) controlled by pCLAMP software (Axon Instruments) was used to generate command pulses and acquire data. Electrodes (Garner Glass) had tip resistances of 1–4 M. For Na⁺ current measurement, the following extracellular solution was used (in mM): 40 NaCl, 5 CsCl, 20 TEA-Cl, 2.5 MgCl₂, 5 CoCl₂, 5 4-aminopyridine (4-AP), 10 glucose, 80 sucrose, and 5 HEPES adjusted to pH 7.3. The Na⁺ current intracellular solution contained (in mM) 110 CsCl, 20 TEA-Cl, 100 aspartic acid, 5 EGTA, 2 MgCl₂, 5 HEPES, 2 MgATP, and 0.1 Na₃GTP adjusted to pH 7.3. For Ca²⁺ current measurement, the extracellular solution contained (in mM) 130 TEA-Cl, 2 CaCl₂, 1 MgCl₂, 10 HEPES, 5 4-AP, and 10 sucrose adjusted to pH 7.3. The intracellular solution was the same as that used for I_Na measurement. For K⁺ current measurement, the pipette was filled with a standard internal solution containing (in mM) 100 L-aspartic acid, 110 KOH, 20 KCl, 2 MgCl₂, 5 EGTA, 5 HEPES, 2 Mg₂ATP, and 0.3 Na₃GTP, with pH adjusted to 7.2 with Tris base. Mouse ventricular myocytes were superfused at room temperature (20–22°C) with a HEPES-buffered Tyrode solution containing (in mM) 135 NaCl, 5.4 KCl, 1 MgCl, 2 CaCl₂, 5 HEPES, and 10 glucose, with pH adjusted to 7.3 with NaOH. I_CaL was blocked by 1 µM nifedipine, and Mg₂ATP in the pipettes suppressed the ATP-sensitive K⁺ current.

Simulation Methods and Simulated Voltage-Clamp Protocols

The model is based on a set of 40 ordinary differential equations solved by a fourth-order Runge-Kutta method, with a time step of 0.0001 ms. A summary of the equations, model parameters, and initial conditions is given in the APPENDIX. The model is nominally adjusted for room temperature of 25°C (298 K). Steady-state initial conditions were obtained by running the model until changes in all variables did not exceed 0.01%. Ryanodine receptor modulation factor P_RyR was set to 0 at the beginning of each AP.

We used two types of simulated voltage-clamp protocols: a double-pulse steady-state inactivation protocol and a variable-interval gapped pulse protocol. Usually, the double-pulse steady-state inactivation protocol involved a pulse (P1) to a voltage between the holding potential and + 50 mV (in 10-mV intervals), followed immediately by a second pulse (P2) to a holding potential specific to the current under investigation. The durations of the two pulses were set appropriately for the gating time constants of the given current. Other protocols are described in the text and figures. For I_CaL, P1 and P2 were separated by a 2-ms return to the holding potential, to match experimental protocols. APs were simulated for at least 100 cycles until a steady state was established, except for simulation of the negative staircase of Ca²⁺ transients, which began from a quiescent steady state.

The variable-interval gapped pulse protocol was used to determine the rate of recovery from inactivation. A P1 pulse was applied to activate and inactivate the channel, followed by a second pulse after a variable time interval. The degree of recovery from inactivation was calculated by comparing the peak P2 current to the peak P1 current. The depolarization, duration, and range of interpulse intervals used were specific to the current being studied. Unstimulated ventricular myocytes are quiescent, so we used an 0.5-ms 80 pA/pF stimulus current (I_stim) applied at a frequency of 0.5–6.7 Hz to trigger simulated APs.

Data Analysis

The electrophysiological data were analyzed by using Clampfit (Axon Instruments), Excel (Microsoft), and Origin (OriginLab). Clampfit was used for the exponential fitting of inactivation kinetics. Other nonlinear curve fittings were performed in Origin. The quality of the fit was evaluated by visual inspection and comparison of ²-values.

	RESULTS

TOP ABSTRACT METHODS RESULTS DISCUSSION APPENDIX: MOUSE MODEL SUMMARY GRANTS REFERENCES

 
We modeled cellular electrical activity in the mouse by using a nondistributed equivalent electrical circuit with subcellular compartmental spaces. This assumes that there are no electrical gradients within the cell itself, i.e., the membrane potential is spatially homogeneous and all subcellular compartments are uniform or "well stirred." A schematic representation of the currents, fluxes, and physical compartments of the model is shown in Fig. 1. In generating APs, the model mouse membrane potential (V) was determined by the following differential equation:

(1)

where C_m is membrane capacitance, I_Na is the fast Na⁺ current, I_CaL is the L-type Ca²⁺ current, I_Kto,f is the rapidly recovering transient outward K⁺ current, I_Kto,s is the slowly recovering transient outward K⁺ current, I_Kr is the rapid delayed rectifier K⁺ current (mERG), I_Kur is the ultrarapidly activating delayed rectifier K⁺ current, I_Kss is the noninactivating steady-state voltage-activated K⁺ current, I_K1 is the time-independent inwardly rectifying K⁺ current, I_Ks is the slow delayed rectifier K⁺ current, I_NaCa is the Na⁺/Ca²⁺ exchange current, I_p(Ca) is the Ca²⁺ pump current, I_NaK is the Na⁺/K⁺ pump current, I_Cl,Ca is the Ca²⁺-activated Cl^– current, I_Cab and I_Nab are the background Ca²⁺ and Na⁺ currents, and I_stim is the externally applied stimulation current. Our model has differential expression of I_Kto,f, I_Kto,s, I_Kur, and I_Kss in myocytes from the apex and the septum region of the mouse heart.

View larger version (28K): [in this window] [in a new window]   Fig. 1. Schematic diagram of the mouse model ionic currents and Ca2+ fluxes. Transmembrane currents are the fast Na+ current (INa), the L-type Ca2+ current (ICaL), the sarcolemmal Ca2+ pump [Ip(Ca)], the Na+/Ca2+ exchanger (INaCa), the rapidly recovering transient outward K+ current (IKto,f), the slowly recovering transient outward K+ current (IKto,s), the rapid delayed rectifier K+ current (IKr), the ultrarapidly activating delayed rectifier K+ current (IKur), the noninactivating steady-state voltage activated K+ current (IKss), the time-independent K+ current (IK1), the slow delayed rectifier K+ current (IKs), the Na+/K+ pump (INaK), the Ca2+-activated Cl– current (ICl,Ca), and the Ca2+ and Na+ background currents (ICab and INab). Istim is the external stimulation current. The Ca2+ fluxes within the cell are uptake of Ca2+ from the cytosol to the network sarcoplasmic reticulum (SR) (Jup), Ca2+ release from the junctional SR (Jrel), Ca2+ flux from the network SR (NSR) to junctional SR (JSR) (Jtr), Ca2+ leak from the SR to the cytosol (Jleak), Ca2+ flux from the subspace volume to the bulk myoplasm (Jxfer), and Ca2+ flux to troponin (Jtrpn). The model includes Ca2+ buffering by troponin and calmodulin in the cytosol and by calsequestrin in the SR. [Ca2+]i, [Na+]i, [K+]i, intracellular Ca2+, Na+, and K+ concentrations; [Ca2+]o, [Na+]o, [K+]o, extracellular Ca2+, Na+, and K+ concentrations.

Fast Na⁺ current I_Na. The fast transient Na⁺ current drives rapid depolarization at the upstroke of the ventricular AP. Most older models of I_Na were based on the Hodgkin and Huxley formalism of the neuronal AP, where inactivation is a single independent and completely absorbing state. I_Na inactivation is fast and relatively complete after depolarization, so the exact gating mechanisms of this channel were not considered in most previous AP reconstructions, despite evidence that a small noninactivated component of I_Na might play a role in repolarization (20, 22–24). Na⁺ channel inactivation has at least two components: fast and slow. Na⁺ channel mutations in the mouse heart that alter the slow component of inactivation can underlie a form of long QT syndrome (a frequently fatal genetic disease in other species). This suggests an important role for I_Na in the pathophysiology of arrhythmias (71). Consequently, more advanced Markov-type models are now beginning to be used to characterize I_Na during the AP (20, 22, 23, 29).

We used a Markov model for I_Na (22) with three closed states (C_Na1, C_Na2, and C_Na3), an open state (O_Na), a fast (IF_Na) and two intermediate (I1_Na and I2_Na) inactivated states, and two closed inactivation states (IC_Na2 and IC_Na3), as shown in Fig. 2. There is a reversible fast inactivation O_Na IF_Na transition, which we have attributed to the initial binding of the inactivation domain to the core domain. The second inactivation transition, IF_Na I1_Na, would then represent stabilization of the inactivation domain-receptor complex. Only a few channels enter the I2_Na state via slow transitions (22). A direct transition pathway between IF_Na, the fast-inactivated state, and the C_Na1 closed state was introduced to match inactivation recovery data (29). Inclusion of two closed-inactivation states, IC_Na2 and IC_Na3, allows for more accurate representation of Na⁺ channel availability (22). The transition rates were based on the model of Clancy and Rudy (22), with parameters adjusted to match mouse data. The transition rates and dynamic equations for I_Na are detailed in the APPENDIX. I_Na is described by the equation:

(2)

where G_Na is the whole cell conductance (mS/µF), E_Na is the Na⁺ reversal potential, and O_Na is the probability of the Na⁺ channel being in the open state.

View larger version (15K): [in this window] [in a new window]   Fig. 2. State diagram of the Markov model for the Na+ channel. CNa1, CNa2, and CNa3 are closed states; ONa is the open state; IFNa is the fast inactivated state; I1Na and I2Na are the intermediate inactivated states; and ICNa2 and ICNa3 are the closed-inactivation states (22). and are the transition rates between the states, as given in the APPENDIX.

View larger version (24K): [in this window] [in a new window]   Fig. 3. The fast Na+ current INa. A: a family of simulated current traces in response to a step depolarization. A 500-ms depolarizing pulse to between –130 and +50 mV was applied (in 10-mV increments) from a holding potential of –140 mV. Only the first 30 ms are shown to demonstrate details of activation and inactivation. B: peak INa-voltage relationship for simulated and experimental data. Data were simulated with 1.8 mM [Ca2+]o and 10 mM (solid line I), 52.5 mM (solid line II), or 140 mM [Na+]o (solid line III). , Our experimental data on mice (n = 20). Simulations with 10 mM (dashed line I) and 52.5 mM [Na+]o (dashed line II) have activation time constants negatively shifted by 15 mV to mimic experimental conditions with 0.5 mM [Ca2+]o [ (86), (51), (9)]. C: INa steady-state inactivation relationship. The simulated (solid line) and experimental steady-state inactivation data from double-pulse protocols are plotted against voltage. Experimental data from mice: (86), (51), (9). D: time constant of inactivation (Na) vs. voltage for simulated (solid line) and experimental data [ (86)]. Dashed line shows calculations of Na with activation time constants negatively shifted by 15 mV to simulate experimental conditions with [Ca2+]o = 0.5 mM, as predicted by surface charge effects (see Ref. 38 for a review).

Recovery from inactivation in the Na⁺ channel was simulated by using a gapped-pulse protocol with steps to –20 mV with an interstimulus holding potential of either –80 mV or –100 mV. The rate of recovery from inactivation is strongly dependent on the interstimulus holding potential, so two very different interpulse intervals (1–50 ms in 1-ms steps and 10–500 ms in 10-ms steps) were used. Fitting simulated recovery curves with a single exponential function showed that the time constant of recovery increased by a factor of 8 when the holding potential was reduced from –80 mV (time constant for recovery from inactivation _rec,Na = 71.5 ms) to –100 mV (_rec,Na = 8.6 ms). This compares well with our data from mouse myocytes (_rec,Na = 66.0 ± 1.9 at –80 mV; n = 5).

One of the advantages of using a Markov model is that transition rates can be modified to simulate the effects of mutations on gating. However, our model of I_Na gating has several important limitations. The most important of these stems from the size of I_Na and the consequential difficulty in measuring this current. I_Na is large and can cause resistive drops of 10 mV or more across the microelectrode tip, resulting in loss of voltage clamp control. Most investigators, therefore, partially block the channel or otherwise reduce conductance to maintain voltage clamp. These manipulations (e.g., reducing [Na⁺]_o or applying partial doses of blocker) are known to alter gating kinetics and equilibrium. Consequently, we also validated our model parameters for I_Na by comparing the calculated and experimental rate of maximal upstroke velocity (dV/dt_max) due to I_Na AP upstroke (see Mouse Action Potential below).

L-type Ca²⁺ current I_CaL. Perturbations in the amplitude and gating kinetics of I_CaL have been associated with many arrhythmias and pathological states such as failing, hypertrophic, stunned, ischemic, and hibernating myocardium (47). Activation of the L-type Ca²⁺ channel is a purely voltage-dependent process, similar to Na⁺ and K⁺ channel activation. The molecular basis of L-type Ca²⁺ channel inactivation is not fully understood. Inactivation of the L-type Ca²⁺ channel is both a time-dependent (i.e., apparently voltage dependent because of its coupling to activation) and a Ca²⁺-dependent process (43, 52). Ca²⁺-induced inactivation of the Ca²⁺ channel appears to be modulated by calmodulin (74), which may bind to the Ca²⁺-binding motif (EF hand) on the carboxy tail of the main _1C-subunit (53), thus transducing calmodulin binding into channel inactivation (75).

Describing the complex inactivation behavior of the L-type Ca²⁺ channel requires assumptions about coupling and inactivation biophysics. Jafri et al. (41) used a mode-switching Markov model for Ca²⁺ inactivation of the L-type Ca²⁺ channel to link inactivation to subspace [Ca²⁺]. The channel is modeled with four independent subunits that can close the channel and whose mode of operation is determined by intracellular Ca²⁺ binding. Transitions between the two modes are controlled by a Ca²⁺-dependent factor, . The probability of entering the Ca²⁺ mode is dependent on voltage, conformation, and subspace [Ca²⁺]. This model makes little use of the emerging structure and function information from voltage-gated channels and fails to model recovery from inactivation.

We developed a new Markov model to describe L-type Ca²⁺ channel dynamics (Fig. 4), based on homology and extrapolated structure-function relationships of Shaker B voltage-gated K⁺ channels. Our model has two inactivation states, one of which is a Markov representation of "ball and chain" inactivation in which the ball binding is Ca²⁺ dependent (74, 75). The model Ca²⁺ channel has four closed states (C₁, C₂, C₃, and C₄), one open state (O), and three inactivated states (I₁, I₂, and I₃). Activation is controlled by the voltage-dependent rate constants, and . Transitions from O to I₁ correspond to a relatively fast Ca²⁺-dependent inactivation, controlled by the Ca²⁺-dependent rate constant . Transitions from O to I₂ represent relatively slow voltage-dependent inactivation and are controlled by the rate constants and K_pcf (forward voltage-insensitive rate constant). The I₃ inactivated state represents a channel with both Ca²⁺- and voltage-dependent inactivation. Direct transitions between the closed state C₄ and the inactivated states I₁, I₂, and I₃ result in a voltage-dependent recovery from inactivation. All rate constants fulfill the condition of thermodynamic reversibility (for review see Ref. 38), i.e., the product of the rate constants in any clockwise loop is equal to that of the same loop measured in a counterclockwise direction.

View larger version (17K): [in this window] [in a new window]   Fig. 4. State diagram of the model L-type Ca2+ channel. C1, C2, C3, and C4 are closed states; O is the open state; and I1, I2, and I3 are inactivated states. The rate constants and are voltage dependent, and is calcium dependent. Kpcf and Kpcb are the forward and backward voltage-insensitive rate constants, respectively.

(3)

where G_CaL is the whole cell conductance (mS/µF), E_Ca,L is the L-type Ca²⁺ channel reversal potential, and O is the channel open probability.

Sarcoplasmic reticulum. No matter how well the Ca²⁺ channel is modeled, it cannot be considered in isolation from its surroundings. In ventricular myocytes, L-type Ca²⁺ channels are found mostly in clusters in the T-tubules directly opposed to the ryanodine receptors on the sarcoplasmic reticulum (SR) (13). Ca²⁺ that enters the cell as part of I_CaL initiates Ca²⁺-induced Ca²⁺ release from the SR, which then contributes to Ca²⁺-dependent inactivation of the Ca²⁺ channel (10, 30, 31). A more physiologically realistic model of the L-type Ca²⁺ channel must therefore include an appropriate representation of the SR release channels and the local [Ca²⁺] near the L-type Ca²⁺ channel.

Ca²⁺ fluxes in and around the SR are shown in Fig. 1. Ca²⁺ entering the cell via through I_CaL results in a localized increase in [Ca²⁺] ([Ca²⁺]_ss) in a restricted subsarcolemmal space. [Ca²⁺]_ss both inactivates the L-type Ca²⁺ channel and initiates Ca²⁺-induced Ca²⁺ release, which further adds to Ca²⁺ channel inactivation. Ca²⁺ diffuses from the subsarcolemmal space to the general cytosol, where it binds to contractile proteins and initiates contraction. Ca²⁺ is removed from the cytosol by translocation across the sarcolemmal membrane by the Na⁺/Ca²⁺ exchanger and the sarcolemmal Ca²⁺ pump and taken back into the network SR by Ca²⁺-ATPase. Ca²⁺ diffuses from the network SR to the junctional SR, where it can be released into the subsarcolemmal space once more. The equations governing the time dependence of the ryanodine receptor are based on those of Keizer and Levine (45). We also modified the Jafri et al. (41) model of Ca²⁺ dynamics to match Ca²⁺ transients reported for the mouse. A Gaussian distribution function and integral of the time course of I_CaL were used to modulate Ca²⁺ release and simulate graded release of Ca²⁺ from the intracellular store (Eq. A15, APPENDIX). Incorporation of local effects of Ca²⁺ release in our mouse model allowed simulation of Ca²⁺ sparks, the elementary units of Ca²⁺ signaling, and a precise description of Ca²⁺ handling and Ca²⁺ channel inactivation that resembles that of more complex biophysically based models (12).

Simulated I_CaL voltage-clamp studies are shown in Fig. 5. A conventional double-pulse protocol was used with a 250-ms P1 pulse stepping from –80 to between –70 and +40 mV followed by a 2-ms return to –80 mV and then a 250-ms P2 pulse to +10 mV. The biexponential current decay is typical for I_CaL from mouse myocytes (see, for example, Ref. 88). The half-time of simulated I_CaL decay at +10 mV was 4.2 ms, which is smaller than the corresponding experimental value of 8 ms obtained by Wang et al. (90). However, the calculated fast inactivation time constant _1,CaL = 4.45 ms compares well to the experimental value of 4.7 ± 0.2 ms obtained by Kirchhefer et al. (49), but is somewhat smaller than 12.28 ± 1.29 ms (82) and 7.92 ± 0.74 ms (81) obtained by Santana et al. The simulated slow inactivation time constant _2,CaL = 16.1 ms is close to the value of 22.47 ± 1.89 ms obtained by Santana et al. (82); however, it is somewhat smaller than the measured values of 56.1 ± 2.6 ms observed by Kirchhefer et al. (49) and 48 ± 1.35 ms recorded by Santana et al. (81). The variability in time constants for the components of inactivation between laboratories probably reflect different experimental conditions (e.g., amount of EGTA, effectiveness of perfusion through electrode tip, isolation procedures, ionic conditions) that can effect both Ca²⁺- and voltage-dependent inactivation. The simulated traces in Fig. 5A mimic control conditions without EGTA or other nonintrinsic buffers and so are closest to the fastest time constants, which we assumed had the least perturbation from altered Ca²⁺ handling. When Ca²⁺-release flux during simulation is decreased by 20-fold to mimic the effect of strong buffering by EGTA, the rate of decay of I_CaL is reduced, with fast and slow time constants equal to 8.8 and 57.3 ms, respectively (V = +10 mV).

View larger version (20K): [in this window] [in a new window]   Fig. 5. The L-type Ca2+ current ICaL. A: a family of simulated current traces. A 250-ms depolarizing first pulse to between –70 and +40 mV (in 10-mV increments) was applied from a holding potential of –80 mV. This was followed by a brief 2-ms return to –80 mV before a second 250-ms pulse to +10 mV. B: peak ICaL voltage relationships for simulated and experimental data. The solid line shows data from simulations under normal physiological extracellular ion concentrations (see APPENDIX), and the dashed line shows data simulating the experimental conditions of Yatani et al. (96). , Our mouse data (n = 5). Other experimental results: (48), (96), (88), (15), (50). C: steady-state inactivation relationships for 500-ms P1 pulse simulations (dashed line), and 5-s P1 pulse simulations with buffered Ca2+ (solid line). , Our data from 5-s P1 pulses and [Ca2+]i buffered with 5 mM EGTA (n = 5); , experimental measurements of Yatani et al. (96), which were buffered with 10 mM BAPTA. D: voltage dependence of simulated fast and slow inactivation time constants 1,CaL and 2,CaL.

Figure 6A shows the bulk intracellular [Ca²⁺] ([Ca²⁺]_i) and subspace volume [Ca²⁺]_ss transients in response to a step depolarization to 0 mV and corresponding experimental data (42). The simulated [Ca²⁺]_i time to peak of 14.8 ms and half-time of decay of 33.6 ms compare well with the experimental values of 17 ms and 34 ms, respectively (90).

View larger version (18K): [in this window] [in a new window]   Fig. 6. Intracellular Ca2+ transients. A: simulated [Ca2+]i transients and Ca2+ subspace volume concentration ([Ca2+]ss) in response to a step depolarization to 0 mV. , Experimental data for [Ca2+]i transients from Ref. 42, normalized to peak simulated concentration. B: voltage dependence of simulated [Ca2+]i () and peak ICaL (). These data show graded release of Ca2+ from the SR.

Figure 7 shows the results of a simulated variable-gap double-pulse protocol. Two 250-ms pulses to 0 mV were applied with an interstimulus interval of between 2 and 500 ms for three different interpulse holding potentials. Recovery from inactivation was approximated by a single exponential function with time constants _rec,Ca equal to 19, 39, and 77 ms for three holding potentials, –90, –80, and –70 mV, respectively. These simulations suggest a moderate voltage dependence of recovery. However, no data exist on the voltage-dependent recovery of I_CaL in mouse myocytes. The L-type Ca²⁺ channel in bullfrog sinus venosus and atrium (17) and in rat and rabbit ventricles (97) shows voltage dependence in this range, with recovery becoming faster with hyperpolarization. Voltage-dependent recovery is also consistent with the general properties of recovery from other voltage-gated channels (see Ref. 38 for a review).

View larger version (21K): [in this window] [in a new window]   Fig. 7. ICaL recovery from inactivation. A: simulated currents from a gapped-pulse protocol used to measure recovery from inactivation. A 250-ms pulse (P1) to 0 mV from the holding potential of –80 mV was followed by a 100-ms pulse (P2) to 0 mV after an interstimulus interval of between 2 and 502 ms, in 25-ms increments. B: peak P2 current plotted against time was fitted with a single exponential function (solid line). , , and , Interpulse holding potentials of –90, –80, and –70 mV, respectively.

(4)

where G_Kto,f is the maximum whole cell conductance (mS/µF), E_K is the K⁺ reversal potential, and a_to,f and i_to,f are the activation and inactivation gating variables. Equations for a_to,f and i_to,f are given in the APPENDIX. This current is equivalent to I_to,f in mouse ventricular myocytes (92).

Mouse ventricular myocytes show substantial regional heterogeneity in repolarization characteristics and can be divided into several types depending on both the anatomic location and the number and type of K⁺ currents present (37, 92). We modeled ventricular myocytes from two regions: the apex and the septum. The fundamental differences between these cells are the APD and density of K⁺ currents present (37, 92).

Figure 8 shows the combined depolarization-activated simulated K⁺ currents (I_K,sum) from the apex and septum when depolarized from a holding potential of –80 mV to between –70 and +50 mV in 10-mV steps. Both simulations are in qualitative and quantitative agreement with the experimental results (92). In our model, there are seven types of K⁺ currents: the rapid transient outward K⁺ current (I_Kto,f), the slow transient outward K⁺ current (I_Kto,s), the rapid delayed rectifier K⁺ current (I_Kr), the ultrarapidly activating delayed rectifier current (I_Kur), the noninactivating steady-state K⁺ current (I_Kss), a very small slow delayed rectifier (I_Ks), and the time-independent K⁺ current (I_K1).

View larger version (14K): [in this window] [in a new window]   Fig. 8. Total depolarization-activated K+ currents (IK,sum) from different regions of mouse heart. Simulated currents were elicited by a 5-s depolarizing step to between –70 and +50 mV in 10-mV increments from a holding potential of –80 mV. A: model apical cell. B: model septal cell.

View larger version (21K): [in this window] [in a new window]   Fig. 9. The rapidly inactivating transient outward K+ current IKto,f. A: simulated double-pulse protocol current traces from a 500-ms pulse to between –100 and +50 mV (in 10-mV increments) from the holding potential of –80 mV followed by a 500-ms second pulse to +50 mV. B: peak IKto,f-voltage relationships for the apex and the septum. Simulated results from the apex and septum are shown as solid lines; , experimental data (apex; Ref. 92). Activation is shifted slightly negative to compensate for the effects of divalent ions used to block Ca2+ currents during the experiments (1). C: steady-state inactivation relationships. The simulated (solid line) relationships were calculated from a double-pulse protocol. Experimental results from mouse: (92) and (87). D: time constants of activation and inactivation (act and inact). Simulated results are shown by bold solid lines. Mouse experimental data for inactivation: (87), (100), and (1 data point only at +40 mV; Ref. 92). Experimental data for activation: (92) and (93). The dashed line shows the calculated activation time constant, shifted positive to fit the experimental conditions (92) with divalent ions.

A simulated gapped-pulse protocol was used to determine the rate of recovery of I_Kto,f from inactivation. Two pulses to +50 mV from the holding potential of –70 mV were applied with an interstimulus interval of 10–500 ms. The simulated recovery time constant was _rec,Kto,f = 27.4 ms, which is close to the experimental value of 27 ms recorded in the mouse (92). The time constants of recovery are identical for I_Kto,f from the apex and the septum (92).

Slow-inactivating transient outward K⁺ current I_Kto,s. The slow-inactivating K⁺ current, I_Kto,s, presumably encoded by Kv1.4, is present in cardiac myocytes from the septum region, but largely absent from other regions of the heart (92). The slowly recovering inactivating K⁺ current might be considered a compensatory current, because it is present in septal cells where there is a considerably smaller I_Kto,f and in ventricular myocytes from mice in which I_Kto,f has been transgenically deleted (7, 36, 69). The model formulation of I_Kto,s is:

(5)

where G_Kto,s is the maximum whole cell conductance (mS/µF) and a_to,s and i_to,s are the activation and inactivation gating variables.

Figure 10A shows simulated current traces in response to a series of step depolarizations from –100 to between –90 and +50 mV in 10-mV increments. The inactivation rate of I_Kto,s is voltage independent, with a time constant of 270 ms that corresponds to the experimental value of 258 ± 15 ms (92). The simulated and experimental peak current-voltage relationships are shown in Fig. 10B. Figure 10C shows the simulated and experimental steady-state activation and inactivation functions of I_Kto,s. Our model was developed to simulate experimental currents recorded in the absence of Ca²⁺ channel blockers (16). Data obtained in the presence of divalent ionic blockers of the Ca²⁺ channel, such as 2 mM CoCl₂, showed a positive shift in activation and inactivation due to the surface charge effect (99). We therefore adjusted our model gating parameters to fit data and account for the shift in activation and inactivation. Figure 10D shows the voltage dependence of time constants of activation for I_Kto,s. The simulated time constant for recovery from inactivation is 1.306 s, which compares well with the experimental value of 1.298 s (92).

View larger version (18K): [in this window] [in a new window]   Fig. 10. The slow inactivating transient outward K+ current IKto,s. A: families of simulated traces for IKto,s for the septum. Depolarizing pulses of 5 s were applied to between –90 and +50 mV in 10-mV increments from the holding potential of –100 mV. B: peak IKto,s-voltage relationships for the septum. Simulated results are shown as a solid line; , experimental data (92). C: steady-state activation and inactivation relationships for the IKto,s. Solid lines show simulated results, and circles with dashed lines show data from Ref. 99. Note that the experimental data (99) are shifted because of the effect of divalent ions (Co2+). D: time constant of activation for IKto,s. , Experimental data from Ref. 99; solid lines, simulated data.

(6)

I_Kr is described by the following equation:

(7)

where G_Kr is the whole cell conductance (mS/µF), O_K is the probability of the channel being in the open state, R is the ideal gas constant, F is the Faraday constant, and [K⁺]_o and [K⁺]_i are the K⁺ concentrations outside and inside the cell, respectively. The driving force for permeation through I_Kr has some permeability to Na⁺ to account for the experimentally observed reversal potential in ventricular myocytes (10 mV positive to E_K). Simulated I_Kr in response to a double-pulse voltage-clamp protocol is shown in Fig. 11. P1 was a 1-s pulse from the holding potential (–80 mV) to between –70 and +50 mV in 10-mV steps, and P2 was a second 1-s pulse to –40 mV. The simulated currents are comparable to experimental results from I_Kr (57). The calculated value of maximum I_Kr amplitude of 0.25 pA/pF was chosen to be close to experimentally observed value in adult mice (57). Our representation of I_Kr as a linear Markov model does not permit inactivation from closed states. Single-channel experiments suggest that these transitions can occur (46), and they have been included in other I_Kr models (21). Kiehn et al. (46) described a relatively rapid flicker state, which allows the channel to inactivate without opening when inactivation is sufficiently fast and complete at positive potentials. There is no experimental evidence that the flicker state is part of the activation process or that it has any consequences for drug binding or burst duration. For macroscopic currents, this flicker is readily addressed as a reduced conductance with a mean lifetime equal to burst duration. We have treated it as such in this model, so as not to complicate the model with a rapid time constant.

View larger version (14K): [in this window] [in a new window]   Fig. 11. The rapid delayed rectifier K+ current IKr: a family of simulated IKr current traces. A 1-s depolarizing pulse from –70 to +60 mV was applied (in 10-mV increments) from a holding potential of –80 mV followed by a 1-s second pulse to –40 mV.

(8)

(9)

where G_Kur and G_Kss are the maximum whole cell conductances (mS/µF), a_ur and a_Kss are activation gates, and i_ur and i_Kss are inactivation gates, described in detail in the APPENDIX. I_Kur and I_Kss correspond to I_K,slow and I_ss, respectively (92). We chose to use the I_Kur terminology to emphasize the rapid activation, as opposed to the slow inactivation that is characteristic of these channels. Use of this terminology also emphasizes the similarity of this current to I_Kur in human atrium (for review, see Ref. 68). These currents are present in both the apex and the septum, but with different expression levels. The model parameters were chosen to fit mouse experimental data (92).

Figure 12 shows simulated currents in response to a series of step depolarizations from –100 to between –90 and +50 mV in 10-mV intervals. The inactivation rate of I_Kur is voltage independent and can be fit with a single exponential function (92). The simulated time constant of inactivation is 1,198 ms, which compares well with the experimental value of 1,180 ± 45 ms (92). I_Kss shows little inactivation in either simulation or experiment (92).

View larger version (14K): [in this window] [in a new window]   Fig. 12. Families of simulated traces for the ultrarapidly activating delayed rectifier current (IKur; A) and the noninactivating steady-state K+ current (IKss; B) for the septum. Depolarizing pulses of 5 s were applied to between –90 and +50 mV in 10-mV increments from the holding potential of –100 mV.

View larger version (20K): [in this window] [in a new window]   Fig. 13. The ultrarapidly activating delayed rectifier K+ current IKur and the noninactivating steady-state K+ current IKss. Simulated peak current-voltage relationships for IKur (circles) and IKss (triangles) are shown. Data from the apex are open symbols; data from the septum are filled symbols. A: simulated data. B: experimental data from Xu et al. (92). C: steady-state inactivation relationships for IKur. Solid line is simulated steady-state inactivation; squares with dashed line show experimental data from Ref. 16. D: IKur and IKss activation time constants. and , Experimental data for IKur and IKss, respectively (92); solid lines, simulated results.

(10)

The equations governing the activation gate, n_Ks, are shown in the APPENDIX. Simulated voltage-clamp data for I_Ks are shown in Fig. 14. I_Ks does not significantly affect AP shape under control conditions; however, it can be of potential importance in experiments where other K⁺ currents are knocked out.

View larger version (15K): [in this window] [in a new window]   Fig. 14. The slow delayed rectifier K+ current IKs: A family of simulated current traces. Depolarizing pulses were applied to between –70 and +50 mV (in 10-mV increments) from the holding potential of –80 mV.

(11)

where G_Cl,Ca is the maximum whole cell conductance (mS/µF), O_Cl,Ca is the probability of the channel being in the open state, K_m,Cl is the half-saturation constant, and E_Cl is the Cl^– reversal potential. This small-amplitude current was recently discovered by Xu et al. (94). Experiments (94) and simulations (not shown) suggest that the current does not affect the AP shape under control conditions. However, the current may be important under pathophysiological conditions that result in abnormal cellular Ca²⁺ loading.

Time-independent K⁺ current I_K1. The I_K1 equation in our model is based on that of DiFrancesco and Noble (27):

(12)

Figure 15 shows simulated and experimental (7, 50, 88) peak current-voltage relationships in response to a pulse from the holding potential of –80 mV to voltages from –150 to –40 mV. Simulated and experimental results give currents of similar magnitude.

View larger version (11K): [in this window] [in a new window]   Fig. 15. The time-independent potassium current IK1. The current-voltage relationship for IK1 from simulations is shown as a solid line with experimental data for comparison [ (88), (50), at –120 mV only (7), and (our data)].

(13)

(14)

where G_Cab is the background Ca²⁺ conductance, G_Nab is the background Na⁺ conductance, E_CaN is the Ca²⁺ reversal potential, and E_Na is the Na⁺ reversal potential. These currents normally maintain ionic homeostasis.

Na⁺/Ca²⁺ exchange current I_NaCa. The equation representing I_NaCa is based on that of Luo and Rudy (62):

(15)

where K_m,Na is the Na⁺ half-saturation constant for I_NaCa, K_m,Ca is the Ca²⁺ half-saturation constant, k_sat is the saturation factor at very negative potentials, and = 0.35 is the position of the energy barrier that controls the voltage dependence of I_NaCa. I_NaCa plays an important role in mouse myocyte Ca²⁺ dynamics (84). The absolute value of I_NaCa does not exceed 1.0 pA/pF within the range –120 to +40 mV. The scaling factor k_NaCa was chosen to ensure equilibrium intracellular ionic concentrations of Na⁺ and Ca²⁺ within experimental values and to match the experimental ratio of Ca²⁺ extruded by the SR ATPase to Ca²⁺ expelled from the cell via I_NaCa during a myocyte twitch.

Sarcolemmal Ca²⁺ pump current I_p(Ca). I_p(Ca) is taken from the Luo and Rudy (62) model:

(16)

where I_p(Ca)^max is the maximum Ca²⁺ pump current and K_m,p(Ca) is the Ca²⁺ half-saturation constant.

Na⁺/K⁺ pump current I_NaK. I_NaK was modeled by using the Luo and Rudy (62) formulation:

(17)

(18)

(19)

where I_NaK^max is the maximum pump current, K_m,Nai is the Na⁺ half-saturation constant, and K_m,Ko is the K⁺ half-saturation constant. This current maintains Na⁺ and K⁺ electrochemical gradients across the cell membrane. I_NaK^max was optimized to reproduce experimental values of quiescent [Na⁺]_i and [K⁺]_i for mouse myocytes.

Quiescent Cellular Properties