Sodium channel subpopulations with distinct biophysical properties and subcellular localization enhance cardiac conduction

Sodium (Na⁺) ionic current carried by voltage-gated Na⁺ channels (Nav) is primarily responsible for driving the propagation of electrical activity in the heart. Na⁺ current (I_Na) drives the rapid depolarization of myocytes that trigger the upstroke of the cardiac action potential (AP). I_Na dysfunction can result in pathological conditions associated with both pathological conduction and repolarization, including Brugada syndrome, long QT syndrome, atrial fibrillation, sick sinus syndrome, and heart failure (Abriel, 2010; Marchal and Remme, 2022; Rivaud et al., 2020). As such, the proper expression and regulation of voltage-gated Na⁺ channels in the heart is essential for healthy cardiac function. The primary voltage-gated Na⁺ channel in cardiac myocytes is Nav1.5. A wide range of proteins has been identified to interact and regulate Nav1.5 (Abriel, 2010; Marchal and Remme, 2022; Rivaud et al., 2020). Different interacting proteins have been found to impact the expression, trafficking, internalization, phosphorylation, and modulation of channel biophysical properties. Additionally, Nav1.5 is known to associate with potentially four different β subunits, which can also impact channel surface expression and biophysical properties, in addition to playing a role in cellular adhesion (Angsutararux et al., 2021; Veeraraghavan et al., 2018; Watanabe et al., 2009; Zhu et al., 2017).

Over the past two decades, multiple studies have provided evidence supporting the existence of so-called “multiple pools” of Nav1.5 in cardiac myocytes, with distinct subcellular localization, specifically referring to different subpopulations of Na⁺ channels localized at the intercalated disk (ID), i.e., at the cell–cell junctions, and at the lateral membrane (Balycheva et al., 2015; Marchal et al., 2021; Marchal and Remme, 2022; Petitprez et al., 2011; Shy et al., 2013). Note that evidence has also suggested a third distinct subpopulation of Na⁺ channels in the transverse tubules; however, in this study, we primarily focus on ID- and lateral membrane-localized Na⁺ channels. In addition to distinct subcellular localization, these subpopulations associate with different interacting and scaffolding proteins, which support the potential for these subpopulations to exhibit distinct biophysical properties. Indeed, macropatch measurements from Lin and colleagues showed that the ID-localized I_Na steady-state activation (SSA) and inactivation (SSI) voltage dependence are shifted compared with lateral membrane-localized I_Na (Lin et al., 2011). Subsequent studies also demonstrated differences in lateral and ID-localized I_Na biophysical properties in ventricular myocytes from a transverse aortic constriction (TAC) model of heart failure (Rivaud et al., 2017) and a connexin43 (Cx43) mutant model (Agullo-Pascual et al., 2014). Interestingly, while these two experimental models exhibit changes in the I_Na current magnitude at the lateral vs. ID membrane, the relative shifts in SSA and SSI are quite comparable with the shifts observed in myocytes from the corresponding sham or control mice. However, while the multiple Nav1.5 pool conceptual model is well-demonstrated across different experimental measurements, including imaging studies, protein pull-down experiments, and patch clamp (Agullo-Pascual et al., 2014; Gillet et al., 2014; Lin et al., 2011; Marchal et al., 2021; Petitprez et al., 2011; Rivaud et al., 2017; Veeraraghavan et al., 2015), the functional roles of these distinct subpopulations are not fully appreciated.

Although it is difficult to quantify the relative fraction of ID vs. lateral membrane-localized Na⁺ channels or I_Na, measurements from Lin et al. (2011) and others (Agullo-Pascual et al., 2014; Marchal et al., 2021; Rivaud et al., 2017) suggest that the majority of I_Na (relative to the entire cell) is localized at the ID for normal conditions. Computational modeling studies have provided strong evidence for the functional role of these ID-localized Na⁺ channels in supporting and modulating electrical conduction (Kucera et al., 2002; Lin and Keener, 2010; Mori et al., 2008). At the ID, Nav1.5 has been shown to cluster with electrical and mechanical junctional proteins, primarily Cx43 in the ventricles and N-cadherin, respectively, which results in the close apposition of Na⁺ channel-dense ID membranes, separated by a narrow (on the order of 5–50 nm) intercellular cleft space. These conditions, specifically the clustering of interacting channels adjacent to a narrow cleft space, support ephaptic coupling, a cell–cell coupling mechanism predicted to occur in cardiac myocytes and other excitable cells (Kucera et al., 2002; Lin and Keener, 2010; Mori et al., 2008; Rhett et al., 2013), such as neurons (Anastassiou et al., 2011).

In brief, this mechanism of ephaptic coupling can be described as follows: consider electrical propagation between two coupled cells, cell 1 and cell 2, for which cell 1 is “upstream” of cell 2, i.e., propagation proceeds from cell 1 to cell 2. In cell 1, rapid inward I_Na at the ID membrane during the AP upstroke increases the cell 1 intracellular potential; concurrently, this same prejunctional I_Na results in a hyperpolarization of the intercellular cleft space between the two cells. The cleft hyperpolarization in turn increases the transmembrane potential across the postjunctional ID membrane in cell 2, which can promote activation of the cell 2 postjunctional I_Na. This cell–cell coupling occurs in parallel with direct electrical current flowing from cell 1 to cell 2 via gap junctions. The role of ephaptic coupling in cardiac electrical propagation was proposed over 50 years ago (Pertsov and Medvinskiĭ, 1976; Sperelakis, 1983; Sperelakis and Mann, 1977; Sperelakis and McConnell, 2002), and subsequent computational studies have predicted that conduction depends on key properties governing ephaptic coupling, including the width of the intercellular cleft space and the relative fractions of I_Na at the ID and the lateral membrane (Greer-Short et al., 2017; Hand and Peskin, 2010; Hichri et al., 2018; Ivanovic and Kucera, 2021; Ivanovic and Kucera, 2022; Kucera et al., 2002; Lin and Keener, 2010, 2013, 2014; Ly and Weinberg, 2022; Moise et al., 2021; Mori et al., 2008; Nowak et al., 2020; Nowak et al., 2021a; Nowak et al., 2021b; Poelzing et al., 2021; Veeraraghavan et al., 2015, 2016; Wei et al., 2016; Wei and Tolkacheva, 2022; Weinberg, 2017; Yu et al., 2022). As the clustering of Na⁺ channels at the ID occurs in conjunction with gap junctions, it is difficult to experimentally separate the relative contribution of gap junction and ephaptic coupling mechanisms, motivating computational modeling studies.

In a seminal study, Kucera and colleagues demonstrated that the impact of ephaptic coupling depended on the relative strength of gap junctional coupling (Kucera et al., 2002). Specifically, a narrow cleft space (thus promoting ephaptic coupling) enhances conduction for weaker gap junctional coupling due to the earlier activation of postjunctional ID I_Na (which was later termed “self-activation”; Hichri et al., 2018). However, for a stronger gap junction, a narrow cleft can slow conduction as the cleft hyperpolarization and postjunctional ID transmembrane potential depolarization can ultimately result in a smaller I_Na driving force, which in turn reduces peak I_Na and thus slows conduction, which was termed “self-attenuation” (Kucera et al., 2002). Subsequent studies have found that conduction due to ephaptic coupling can further depend on additional properties, such as the relative fraction of ID-localized I_Na and total I_Na conductance (across the whole cell; Nowak et al., 2021a), and structural properties, including heterogeneous structure within different regions of the ID (Hichri et al., 2018; Ivanovic and Kucera, 2021; Ivanovic and Kucera, 2022; Moise et al., 2021), the separation between Na⁺ channels and gap junctions on the ID membrane (Hichri et al., 2018; Ivanovic and Kucera, 2021), and cell size (Nowak et al., 2021b).

While these computational studies have provided important insights and predictions as to how the subcellular distribution of Na⁺ channels can influence cardiac conduction via ephaptic and gap junctional coupling, all of these prior studies have assumed that the biophysical properties of Na⁺ channels localized at the ID and lateral membrane are the same. In this study, we perform a series of simulations in both single cells and tissue to predict the impact of subpopulations of Na⁺ channels with distinct biophysical properties on cardiac conduction. In tissue simulations, we consider several different subcellular distributions to identify the contributions of Na⁺ channels with distinct spatial localization, Na⁺ channels with shifted biophysical properties, and the contribution of both occurring in tandem. Collectively, our results support the conclusion that ID-localized Na⁺ channels with shifted biophysical properties enhance cardiac conduction, with enhanced sensitivity to cleft width changes (and thus greater adaptability to ID structural perturbations). Further, while simulations predict that lateral membrane-localized Na⁺ channels do contribute to the total cell I_Na, ID-localized Na⁺ channels provide the greatest contribution to the total Na⁺ charge, a property that is enhanced for weaker gap junctional coupling.

To account for a subpopulation of Na⁺ channels with different biophysical properties, we define a fraction of I_Na current governed by “shifted” steady-state and time constant voltage-dependence. That is, we define steady-state relationships $m∞shift(V)=m∞base(V−Vactshift)$⁠, $h∞shift(V)=h∞base(V−Vinactshift)$⁠, and $j∞shift(V)=j∞base(V−Vinactshift);$ and time constants $τmshift(V)=τmbase(V−Vactshift)$⁠, $τhshift(V)=τhbase(V−Vinactshift),$ and $τjshift(V)=τjbase(V−Vinactshift).$ Specifically, we define the activation rates $αmshift=αm(V−Vactshift)$ and $βmshift=βm(V−Vactshift),$ and the inactivation rates $αhshift=αh(V−Vinactshift),$$βhshift=βh(V−Vinactshift)$⁠, $αjshift=αj(V−Vinactshift),$ and $βjshift=βj(V−Vinactshift).$ Based on measurements of the steady-state activation and inactivation curves from Lin et al. (2011), we define $Vactshift=−5.5mV$ and $Vinactshift=7.8mV,$ unless otherwise stated. The shifted SSA $m∞shift(V)$ and SSI $h∞shift(V)$ curves are shown for the shifted model parameters in Fig. 1 A (red). Note that to account for differences in species and other experimental conditions, we incorporated the relative shift in SSA and SSI between ID and lateral membrane I_Na, not the exact values for half-activation or -inactivation, from Lin et al. (2011).

Here, g_Na is the total sodium conductance (16 mS/cm²) and E_Na is the sodium reversal potential. Importantly, we highlight that, while these single-cell simulations account for varying proportions of the two subpopulations with baseline and shifted I_Na biophysical properties, they do not account for different spatial localization of the subpopulations.

To investigate the impact of Na⁺ channel subpopulations with distinct biophysical properties and subcellular spatial locations, we simulate a tissue model of a one-dimensional chain of 50 cells, previously used by us and others (Kucera et al., 2002; Wei et al., 2016; Greer-Short et al., 2017; Weinberg, 2017; Nowak et al., 2021a), which accounts for gap junction and ephaptic coupling (Fig. 1 B). We account for the non-uniform spatial localization of Na⁺ channels by spatially discretizing each cell into three membrane “patches,” specifically the postjunctional, lateral, and prejunctional membranes, with total membrane currents denoted as I_post, I_lat, and I_pre, respectively. Thus, each cell is represented by three intracellular nodes, with intracellular resistance modeled by two resistors, R_i = ρ_iL/(2πr²), where L is the cell length of 100 μm, r is the cell radius of 11 μm, and intracellular resistivity $ρi=150Ω⋅cm$⁠. Capacitance of each membrane is proportional to surface area, such that the lateral membrane capacitance C_lat = 2πrLc_m and pre- and postjunctional intercalated disk membrane capacitance C_d = πr²c_m, where capacitance per unit area c_m = 1 μF/cm².

The dynamics of the intercellular cleft voltages ϕ_cleft (between adjacent cells) are governed by the corresponding pre- and postjunctional membrane currents I_pre and I_post, respectively, the cleft radial resistance (Kucera et al., 2002) R_c = ρ_cl/(8πw), where w is the cleft width (varied between 5 and 60 nm; Greer-Short et al., 2017; Veeraraghavan et al., 2015; Veeraraghavan et al., 2018; Wu et al., 2021), and cleft resistivity $ρcl=150Ω⋅cm$⁠. Thus, cleft resistivity and width are inversely related. Gap junctional coupling is represented by a resistor, with conductance g_gap, directly connecting pre- and postjunctional intracellular spaces. In this study, we vary g_gap from 7.67 to 767 nS based on the range of values previously measured experimentally (Desplantez et al., 2007; Kwak and Jongsma, 1996; McCain et al., 2012; Moreno et al., 1994; Nielsen et al., 2012; Rudisuli and Weingart, 1989; Valiunas et al., 2002; Verheule et al., 1997; Weingart, 1986; White et al., 1990; Wittenberg et al., 1986).

Estimates based on macropatch measurements suggest that between 50 and 90% of Na⁺ current is localized at the ID (Lin et al., 2011); therefore, in this study, we consider different values for the fraction of current with shifted I_Na and for the fraction of current localized at the ID junctional membranes. Parameter f_Na defines the fraction of shifted I_Na and/or the fraction of junctional membrane-localized I_Na, depending on the Na⁺ channel distribution, as described below and in Table 1.

We consider seven distinct Na⁺ channel distributions, which collectively enable a thorough investigation of the role of different I_Na biophysical properties and spatial localization (Fig. 1 C), where we denote the fraction of the total Na⁺ current conductance of each I_Na subpopulation on each membrane: (1) baseline/lateral—all baseline I_Na localized on the lateral membrane; (2) shifted/lateral—all shifted I_Na localized on the lateral membrane; (3) mix/lateral—a mix of f_Na shifted I_Na and (1 − f_Na) baseline I_Na subpopulations, all localized on the lateral membrane; (4) baseline/polarized—(1 − f_Na) baseline I_Na localized on the lateral membrane and f_Na baseline I_Na localized on junctional membranes; (5) shifted/polarized—(1 − f_Na) shifted I_Na localized on the lateral membrane and f_Na shifted I_Na localized on junctional membranes; (6) mix/heterogeneous-polarized—(1 − f_Na)² baseline I_Na and f_Na(1 − f_Na) shifted I_Na localized on the lateral membrane, and f_Na(1 − f_Na) baseline I_Na and f_Na² localized on the junctional membrane; and (7) mix/homogeneous-polarized—(1 − f_Na) baseline I_Na localized on the lateral membrane and f_Na shifted I_Na localized on the junctional membrane. Note the junctional I_Na is evenly split between the pre- and postjunctional membranes. The fraction of the total Na⁺ current conductance of each I_Na subpopulation on each membrane is shown in Table 1, with the specific f_Na = 0.7 shown for clarity.

We highlight that the mix/homogeneous-polarized distribution is consistent with a physiological representation of the I_Na distribution, specifically channels with shifted biophysical properties localized at the junctional membranes, while channels with baseline properties are localized on the lateral membrane. In contrast, the baseline/lateral distribution represents the typical model assumption, i.e., consistent with the standard monodomain model, in which all I_Na is represented with baseline biophysical properties and localized on the lateral membrane. The baseline/polarized distribution is consistent with prior studies of ephaptic coupling, in which a fraction of Na⁺ channels is localized at the junctional membranes but with the same I_Na biophysical properties as channels localized on the lateral membrane. We also highlight the difference between the mix/heterogeneous-polarized and mix/homogeneous-polarized distributions. The later distribution, again representing physiological conditions, and the former distribution both have the same total fraction of I_Na with shifted biophysical properties and the same total fraction of I_Na localized at the junctional membranes. However, the heterogeneous distribution has a mix of both I_Na subpopulations on both the lateral and junctional membranes, while for the homogeneous distribution, the two I_Na subpopulations are localized in spatially distinct regions. Importantly, unless otherwise stated, the total I_Na conductance (for a given cell) is the same for all seven Na⁺ channel distributions, with the differences only arising due to the I_Na subpopulations and subcellular localization.

Numerical experiments were performed by varying key model parameters in single-cell and tissue simulations. In single cells, we vary the fraction of shifted I_Na (f_Na) and quantify the I_Na peak and total Na⁺ charged carried Q_Na. In tissue, we vary the Na⁺ channel distributions, cleft width w, gap junction conductance g_gap, and the fraction of Na⁺ channels with shifted biophysical properties and localized at the junctional membranes (f_Na), as described in Table 1. For each tissue simulation, we quantify the conduction velocity (CV), and in a subset of simulations, the total charge carried by I_Na on the lateral and junctional membranes. Unless otherwise stated, CV was measured following electrical wave initiation from resting conditions.

In single-cell simulations, electrical activity initiated an applied stimulus of 1 ms duration for a cycle length of 500 ms. For single-cell simulations, activation time was calculated as the first time the intracellular voltage increases above +10 mV for single cells (to avoid the timing of the stimulus). The ordinary differential equations for gating variables and voltage are integrated using the MATLAB (The Mathworks, Inc.) ode solver ode15s. In the tissue, a propagating wave was initiated by simulating cells 1–5. For each tissue simulation, activation time was calculated as the first time the intracellular voltage increased above −60 mV for each membrane patch. CV was calculated by measuring the difference in activation time between the middle 50% of the one-dimensional tissue (cells 13–38). Note that a tissue length of 50 cells was found to be sufficiently long to avoid boundary and stimulus artifacts on CV measurements. Tissue simulations are solved using an operator splitting method, as described in detail in our recent work (Moise et al., 2021). In brief, gating variables for each membrane patch are integrated using forward Euler or Rush–Larson methods (for Hodgkin–Huxley type variables), which are subsequently used to update membrane ionic currents. The resulting system of ordinary differential equations and algebraic expressions governed by the electrical circuit (Fig. 1 B) are integrated using the backward Euler method to update intracellular and cleft voltages. For both single cells and tissue, initial conditions are set to steady-state values consistent with the resting potential of −84.5 mV.

Fig. S1 shows that negative shifts in steady-state activation and positive shifts in steady-state inactivation promote earlier upstroke time (t_up) and a larger contribution from the shifted Na⁺ channels to the total Na⁺ charge (Q_Na) in single cells. Fig. S2 shows that larger fractions of shifted I_Na and polarized I_Na result in greater dependence on Na⁺ channel distribution. Fig. S3 shows that larger fractions of shifted I_Na and polarized I_Na result in greater dependence on Na⁺ channel distribution. Fig. S4 shows that postjunctional and lateral I_Na activate earlier for narrow clefts and the mix/homogeneous-polarized distribution. Fig. S5 shows that postjunctional I_Na contribution is enhanced for moderate gap junction coupling. Fig. S6 shows the membrane-specific I_Na block in the mix/homogeneous-polarized distribution. Fig. S7 shows membrane-specific I_Na block in the mix/homogeneous-polarized distribution. Fig. S8 shows membrane-specific I_Na block in the mix/homogeneous-polarized distribution. Fig. S9 shows that the mix/homogeneous-polarized distribution exhibits faster conduction compared with the baseline/polarized distribution for all diastolic intervals (DIs). Fig. S10 shows that outward prejunctional I_Na drives the Na⁺ transfer mechanism for narrow clefts and large total I_Na conductance. Fig. S11 shows that outward prejunctional I_Na drives the Na⁺ transfer mechanism for lower gap junctional coupling and large total I_Na conductance. Fig. S12 shows that the addition of a subpopulation of I_Na with shifted biophysical properties results in earlier AP upstroke and larger I_Na current in single-cell mouse and human models.

We first investigate the impact of mixed I_Na populations in a single cell. Critically, in the single-cell model, all ion channels are represented as localized on the same membrane patch, i.e., there is no I_Na spatial distribution. We first consider the case of a mixed 50:50 I_Na distribution (f_Na = 0.5), equally split between baseline and shifted I_Na (magenta [voltage traces], red/blue/magenta [I_Na traces]), and compare with the nominal homogeneous model, with only baseline I_Na (black). At the level of the single-cell AP, we find almost no difference between the cell with the mixed I_Na subpopulations and only baseline I_Na, except for small changes in the peak voltage (Fig. 2 A). However, closer examination reveals that the presence of the shifted I_Na results in an earlier AP upstroke (magenta) compared with the only baseline I_Na model (black). We find that the early AP upstroke is due to the subpopulation of shifted I_Na activating earlier (red), which in turn activates the subpopulation of baseline I_Na earlier (blue), compared with the only baseline I_Na model (black). Additionally, in the cell with both subpopulations, the peak of the total I_Na (magenta) is greater than the peak of the baseline I_Na (black). That is, even though the total I_Na conductances of the two cells are equal, the magnitude of the total I_Na in the cell with the shifted I_Na subpopulation is larger compared with the baseline I_Na.

In a single cell with 70% shifted I_Na (f_Na = 0.7), we find that the effects described above are magnified (Fig. 2 B). The AP upstroke is earlier, compared with the 50–50% distribution, and the shifted I_Na subpopulation (red) is greatly enhanced, such that the total I_Na (magenta) in the cell with the shifted I_Na subpopulation is much larger in magnitude compared with the only baseline I_Na (black).

We next quantified how the I_Na peak and total Na⁺ charge (Q_Na) for the two I_Na subpopulations vary as the fraction of the shifted I_Na subpopulation changes in a single cell (Fig. 3). We first measure the ratio of the shifted I_Na peak to the baseline I_Na peak as a function of f_Na (Fig. 3 A). As f_Na increases, the shifted I_Na peak ratio also increases, as expected. However, we can compare this increase with expectations based solely on the change in conductance between the two subpopulations. Based on conductance changes alone, we would expect the ratio between the two peaks to scale with f_Na/(1 − f_Na) (dashed line). Critically, the actual simulated ratio (black) is always greater than this expectation, meaning that the shifted I_Na peak is always proportionally larger than the shifted I_Na conductance.

We also quantified the total charge carried by the two I_Na subpopulations and measured the fraction of the total Na⁺ charge (Q_Na) carried by the shifted I_Na subpopulation (Fig. 3 B). Based on conductance alone, we would expect this fraction to scale with f_Na (dashed line). However, we find that the shifted Q_Na (black) is always greater than this expectation as well. That is, the shifted I_Na subpopulation carries a proportion of the total Na⁺ charge greater than its respective conductance. Plots of the I_Na and Q_Na ratios against f_Na/(1 − f_Na) similarly show that the shifted I_Na subpopulation proportionally contributes more to the total Na⁺ current (Fig. 3, C and D), compared with the baseline I_Na, than the expectation based on conductance (dashed line). In Fig. S1, we measure the action potential upstroke time and the fraction of shifted Q_Na for different values of $Vactshift$ and $Vinactshift$⁠. For all values of f_Na, a negative (left) shift in SSA and positive (right) shift in SSI both promote an earlier upstroke time and larger contribution to Q_Na from the shifted I_Na subpopulation, with greater sensitivity to $Vactshift$⁠.

The above results demonstrate that in a single cell, channels with shifted I_Na biophysical properties contribute to an earlier AP upstroke. We next systematically investigate to what extent conduction in a one-dimensional tissue depends on the presence of Na⁺ channels with shifted I_Na and their localization at the ID. Further, we investigate how these properties also depend on ID structure, specifically the intercellular cleft width and the strength of gap junctional coupling.

As described in the Materials and methods, we consider seven distinct Na⁺ channel distributions, chosen to elucidate the role of both the presence of channels with shifted I_Na biophysical properties and their localization at the ID. In Fig. 4, we plot CV as a function of the cleft width for different values of the gap junction conductance g_gap. As expected, for all Na⁺ channel distributions, CV increases with increasing g_gap. For the three distributions with only lateral Na⁺ channels (Fig. 4, A–C), CV does not depend on cleft width, also as expected, since these distributions do not engage ephaptic coupling. Of these three lateral distributions, conduction is the fastest in the tissue with only shifted I_Na (Fig. 4 B), slightly faster than in tissue with a 30:70 baseline:shifted mixed I_Na distribution (Fig. 4 C), with the only baseline I_Na the slowest (Fig. 4 A).

The presence of the shifted I_Na subpopulation in the polarized distributions similarly results in faster conduction. In contrast with the lateral distributions, the four polarized distributions (Fig. 4, D–G) exhibit a clear dependence on the cleft width. Consistent with prior studies (Kucera et al., 2002; Mori et al., 2008; Lin and Keener, 2013; Weinberg, 2017; Nowak et al., 2021b), for all polarized distributions, narrowing the cleft tends to slow conduction for stronger gap junctional coupling, while narrowing the cleft tends to enhance conduction for weaker gap junctional coupling. Further, we find that the tissue with shifted/polarized and mixed/polarized I_Na distributions exhibit greater sensitivity to the cleft width compared with the baseline/polarized distribution.

In Fig. 5, we plot the same CV measurements as a function of g_gap for different values of cleft width w. The lateral I_Na distributions similarly exhibit no dependence on cleft width (Fig. 5, A–C). From this presentation, it is also more clear to what extent CV varies for the four polarized distributions, as a function of cleft width and g_gap (Fig. 5, D–G). As noted above, cleft narrowing results in conduction slowing for high g_gap and conduction enhancing for low g_gap. Importantly, we find that the mixed/homogeneous-polarized distribution depends on cleft width for the widest range of g_gap values and demonstrates the greatest sensitivity to cleft width (for a given g_gap).

In Fig. 6, we plot CV as a function of g_gap and directly compare a subset of the Na⁺ channel distributions. For all cleft widths and most g_gap values (except the extreme high and low g_gap), conduction is fastest for the mixed/homogeneous-polarized (blue), followed by the mixed/heterogeneous-polarized (magenta), mixed/lateral (black), baseline/polarized (red), and baseline/lateral (green). We also show ratios of the CV values for several combinations of Na⁺ channel distributions to identify to what extent Na⁺ channel composition vs. localization regulates conduction. Comparing the polarized vs lateral distribution with baseline I_Na, the polarized distribution has faster CV, i.e., ratio >1, for a range of g_gap values between about 50 and 500 nS, but slower CV for larger or smaller g_gap outside this range (Fig. 6 B, solid black). Comparing the polarized vs. lateral distribution for the mixed I_Na distributions, the polarized distributions is similarly faster for a g_gap range and slower for outside the range (Fig. 6 B, dashed black). However, the g_gap range for enhanced CV is widened to about 20–500 nS. That is, polarizing the mixed I_Na subpopulation enhances conduction for a wider range of g_gap.

We can also directly compare the mixed/homogeneous-polarized vs. mixed/heterogeneous-polarized distributions, which have the same fraction of baseline vs. shifted I_Na current and the same fraction of I_Na polarized at the ID, with the only difference being whether or not the I_Na subpopulations are spatially distinct. We find that for nearly all g_gap values, this CV ratio >1, i.e., polarization specifically of the shifted I_Na subpopulation enhances conduction (Fig. 6 B, dotted line).

Comparing the mixed/lateral vs. baseline/lateral distribution, we find that the CV ratio >1 for all g_gap values and the ratio increases as g_gap is decreased (Fig. 6 C, solid line). That is, with only lateral Na⁺ channels, the mixed I_Na subpopulations always enhance conduction, relative to baseline I_Na, and, to a greater extent, for weak gap junctional coupling. Interestingly, when we make the comparable comparison in polarized tissue (mixed/homogeneous-polarized vs. baseline/polarized), we observe the same general trend with the addition of an additive bump between g_gap of 30 and 200 nS that is greater for narrower cleft width (Fig. 6 C, dashed line). Across all distributions and g_gap, the magnitude of these effects generally decreases as cleft width increases.

Collectively, these results support the following general conclusions: (1) the inclusion of a subpopulation with shifted I_Na biophysical properties speeds up conduction, with CV enhancement increasing as gap junctional coupling weakens; (2) polarization of any Na⁺ channels enhances conduction over a wide range of gap junctional coupling levels; (3) polarization of the shifted I_Na subpopulation results in an additional conduction enhancement (over the g_gap range of 30–200 nS); and (4) the influence of Na⁺ channel polarization, for all polarized distributions, on conduction depends on the cleft width, with both conduction enhancement and slowing effects diminished with wider cleft widths.

In Figs. S2 and S3, we plot these same CV ratios for different combinations of Na⁺ channel distributions for different values of f_Na. For nearly all cases, the magnitude of the relationships described above is larger as f_Na is increased. That is, the influence of the shifted I_Na subpopulation and the polarization of specific Na⁺ channel subpopulations is magnified as the subpopulation fraction is increased.

In Fig. 7, we show contour maps of CV for the baseline/polarized, mix/homogeneous-polarized, and shifted/polarized distributions as functions of cleft width and f_Na for different values of gap junction conductance g_gap. Importantly, these maps illustrate that the range of CV changes due solely to tissue structural changes (i.e., cleft width variation) and Na⁺ channel subpopulation redistribution (i.e., varying f_Na) for a fixed degree of gap junctional coupling. Several key points are apparent from these plots: For all conditions, CV is larger for the mix/homogeneous-polarized distribution, relative to baseline/polarized, but slower than the shifted/polarized distribution. Further, CV depends on both cleft width and f_Na; however, this relationship depends on gap junctional coupling. For low g_gap, CV predominantly depends on cleft width, while for high g_gap, CV predominantly depends on f_Na, with moderate g_gap exhibiting comparable dependence on both. Thus, for weak gap junctional coupling, CV is most effectively modulated by changes in cleft width, while for strong gap junctional coupling, CV is most effectively modulated by Na⁺ channel redistribution. Additionally, the number of contour regions (with 1 cm/s spacing) illustrates that CV for the mix/homogenous-polarized distribution exhibits the greatest variability compared with either the baseline/polarized or shifted/polarized distributions. That is, the physiological distribution, via cleft width changes and Na⁺ channel redistribution, exhibits the broadest range of CV values for fixed gap junctional coupling. Further, this variability in CV is wider for moderate g_gap, i.e., CV can be modulated to the greatest degree for moderate levels of gap junctional coupling.

We next sought to understand the mechanisms underlying the changes in conduction due to Na⁺ channel biophysical properties and distribution. In Fig. 8 A (top), we plot the activation time for the baseline/polarized (black) and mix/homogeneous-polarized (magenta) distributions as a function of the cell number, for which the time of postjunctional, lateral, and prejunctional membrane activation is denoted for each cell. Both simulations have moderate gap junctional coupling (g_gap = 100 nS) and I_Na polarization (f_Na = 0.7). Consistent with the faster CV measurements above, the activation time is earlier for the mix/homogeneous-polarized distribution. For direct comparison, we also shifted the curves such that the two align at the activation time of either the prejunctional membrane of cell 24 (Fig. 7 A, middle) or the postjunctional membrane of cell 25 (bottom). From these two plots, we can identify that the primary difference in activation time is due to the delay in time between the cell 24 prejunctional membrane and the cell 25 postjunctional membrane, i.e., the time for the electrical wave to propagate across the cell–cell junction. In contrast, the delay between the postjunctional and prejunctional membrane of cell 25 is near identical for the two distributions, i.e., there is no difference in the time for the electrical wave to propagate across the intracellular space.

We next investigate the voltage and I_Na changes associated with these two cases (Fig. 8, B and C). We observe a slower AP upstroke in the baseline/polarized distribution compared with the mix/homogeneous-polarized distribution (Fig. 8, B and C, row 1), consistent with the single cell simulations in Fig. 2. In both distributions, the cleft space between the cells is hyperpolarized, consistent with prior studies of ephaptic coupling mechanisms (Kucera et al., 2002; Hichri et al., 2018; Ivanovic and Kucera, 2021; Moise et al., 2021; Nowak et al., 2021b; Fig. 8, B and C, row 2). However, we find a larger magnitude hyperpolarization in the mix/homogeneous-polarized distributions (magenta). Further, while there are small differences in time course and magnitude, the gap junctional current is comparable in both distributions, highlighting that the differences in CV are directly attributable to the differences in Na⁺ channel distribution.

However, examination of the I_Na reveals critical differences between the two distributions, specifically considering the cell 24 lateral and prejunctional I_Na; cell 25 postjunctional, lateral, and prejunctional I_Na; and cell 26 postjunctional and lateral I_Na. For the baseline/polarized distribution, for a given cell (here, cell 25), we find that the postjunctional I_Na (dashed green) activates earliest, as expected. However, the lateral (solid black) and prejunctional (solid green) I_Na activate at almost the same time, with the prejunctional I_Na peak slightly earlier. Additionally, the postjunctional peak I_Na is the largest in magnitude, while the lateral and prejunctional peak I_Na are almost identical in magnitude. For the mix/homogeneous-polarized distribution, again for a given cell (considering cell 25), the postjunctional I_Na (dashed red) similarly activates earliest and has the largest peak. However, the prejunctional I_Na (solid red) activates next and has the second largest peak, while the lateral I_Na (solid blue) activates last and has the smallest peak. That is, due to the shifted I_Na biophysical properties of the junctional I_Na current, the prejunctional I_Na activates before the lateral membrane, despite being further “downstream” in the cell. In Fig. S4, we quantify the delay between the timing of the I_Na peak for the postjunctional, lateral, and prejunctional I_Na for cell 25 and the AP activation of the cell 24 postjunctional membrane, for different values of f_Na and cleft width. We find that the postjunctional I_Na peak is consistently earliest and increasingly so for larger f_Na and smaller cleft width. For nearly all cases and both distributions, the prejunctional I_Na peak is earlier than the lateral membrane. However, the difference in peak timing is greatest for the mix/homogeneous-polarized distribution, compared with the baseline/polarized distribution, and for larger f_Na and smaller cleft width.

For the value of f_Na = 0.7 (Fig. 8), the total I_Na conductance distribution is thus 35, 30, and 35% for the postjunctional, lateral, and prejunctional membranes, respectively. However, for both Na⁺ channel distributions considered in Fig. 8, the postjunctional I_Na is disproportionally larger compared to this expectation based on conductance. Further, for the mix/homogeneous-polarized distribution, the lateral I_Na is disproportionally smaller. That is, the postjunctional I_Na contributes more to the total I_Na, while the lateral I_Na contributes less, than simply the expected proportion based on the I_Na conductance distribution.

We quantify these trends by plotting the fraction of total Na⁺ charge carried by I_Na on each membrane (Q_Na) as a function of f_Na and for different cleft widths (Fig. 9). Based on conductance distribution, the post- and prejunctional membranes would be expected to scale with f_Na/2, while the lateral membrane would be expected to scale with 1 − f_Na, shown as dashed black lines. As expected, both post- and prejunctional Q_Na increase with f_Na, while lateral Q_Na decreases with f_Na. We find that for both distributions, postjunctional Q_Na is larger than expected based on the conductance distribution (Fig. 9 A), with Q_Na for the mix/homogeneous-polarized distribution consistently larger than the baseline/polarized distribution. Both of these trends are greatest for narrow cleft width.

Interestingly, the lateral Q_Na for the baseline/polarized distribution is nearly identical to the conductance expectation (dashed line), i.e., the lateral I_Na charge is proportional to its conductance (Fig. 8 B). In contrast, lateral Q_Na for the mix/homogeneous-polarized distribution is less than the conductance distribution. Due the opposing trends for the postjunctional and lateral Q_Na, the prejunctional Q_Na fraction trends are inconsistent. For smaller f_Na, the mix/homogeneous-polarized distribution prejunctional Q_Na fraction is larger than the baseline/polarized distribution, and to a greater extent for wider cleft width.

Thus, in summary, for both distributions, the postjunctional I_Na is consistently contributing the largest fraction of the total Na⁺ charge, with a greater contribution for the mix/homogeneous-polarized distribution. As a consequence, the lateral membrane I_Na contributes proportionally less to the total I_Na, more so for tissues with the mixed subpopulation.

In Fig. 10, we plot the Q_Na fractions as functions of g_gap for different values of f_Na and cleft width for the mix/homogeneous-polarized distribution. For all cleft width and g_gap, post- and prejunctional Q_Na increase with increasing f_Na, while lateral Q_Na decreases, as expected. However, Q_Na fractions exhibit an interesting dependence on g_gap. For all cases, the lateral Q_Na decreases as g_gap decreases (Fig. 10 B). However, postjunctional Q_Na exhibits a biphasic dependence on g_gap, with a positive bump between around 30 and 300 nS (Fig. 10 A), while prejunctional Q_Na exhibits a similar negative bump in the same g_gap range (Fig. 10 C). Thus, collectively these results demonstrate that as gap junctional coupling weakens, lateral Q_Na decreases, such that the two junctional membranes constitute a larger Q_Na contribution. Further, the postjunctional Q_Na is always larger than the prejunctional Q_Na. However, this disparity is greatest for g_gap values around 30–300 nS, with prejunctional Q_Na disproportionally largest, while for extreme cases of either high or low gap junctional coupling, post- and prejunctional Q_Na are much similar in magnitude. The magnitude of these trends is generally smaller for wider cleft width. In Fig. S5, we observed similar trends for the baseline/polarized distribution, with the lateral Q_Na generally contributing a larger proportion, consistent with the result in Fig. 10.

To further demonstrate the functional role of the ID-localized I_Na, we measure CV for the membrane-specific block in the mix/homogeneous-polarized distribution for f_Na = 0.7 (Fig. 11). Specifically, we reduce the I_Na conductance on either the ID (red) or lateral (blue) membranes, while maintaining the same conductance on the other membrane. We consider the full range from the nominal mix/homogeneous-polarized distribution with 70% ID-localized shifted Na⁺ channels and 30% lateral membrane-localized baseline Na⁺ channels and decrease one of these subpopulations individually until only the other subpopulation is present. We note that this results in a decrease in overall cell I_Na conductance. We find that reducing the conductance of either the ID- or lateral membrane-localized I_Na slows conduction. However, across all conditions (cleft width and gap junction conductance), we find that reducing ID-localized I_Na results in greater conduction slowing, with generally a greater disparity for larger g_gap. Further, in Fig. S6, we quantify the Q_Na fraction for the lateral membrane as a function of the I_Na conductance reduction on either the ID (red) or lateral (blue) membranes. Importantly, we find that as the ID membrane-localized I_Na conductance decreases (i.e., greater block), the Q_Na fraction for the lateral membrane increases, i.e., lateral membrane I_Na provides a larger contribution when ID-localized I_Na is reduced. We also find that, as would be expected, blocking the lateralmembrane-localized I_Na decreases the Q_Na fraction for the lateral membrane. These relationships are also consistent across different values of cleft width and gap junction conductance. Additionally, in Figs. S7 and S8, we find the same trends for f_Na = 0.5.

All of the previous simulations considered conduction from rest. We next consider how different Na⁺ channel distributions impact conduction for faster pacing rates, i.e., shorter basic cycle length (BCLs). For the baseline/polarized and mix/homogeneous-polarized distributions, we paced the tissues at different BCL values until steady state or a steady-state alternating pattern was reached and measured CV for different cleft widths and g_gap (Fig. 12). Consistent with the above results, across all BCLs, conduction is faster for the mix/homogeneous-polarized distribution (magenta) compared with the baseline/polarized distribution (black). Additionally, for all cleft widths and g_gap, tissues with the mix/homogeneous-polarized distribution could be paced at shorter BCLs, before loss of 1:1 capture or conduction block. Additionally, for faster pacing rates, the tissues exhibited spatially concordant alternans, in which action potential duration alternates on a beat-to-beat basis. This alternation occurs concurrently with a beat-to-beat alternation in CV, such that the CV vs. BCL bifurcates at this alternans onset cycle length. We find that alternans onset is consistently at a shorter BCL for the mix/homogeneous-polarized distribution, i.e., faster pacing is required to induce alternans. In Fig. S9, we plot CV against the preceding diastolic interval (DI; i.e., the CV restitution curve) for these conditions. We find that the CV restitution curve extends to shorter DI values for the mix/homogeneous-polarized distribution, and in most cases, the curve is flatter, i.e., CV is less sensitive to changes in DI.

We additionally investigate to what extent conduction differences due to the mixed distribution depend on the shift in either the SSA or SSI curves. Thus, we consider cases in which the shifted I_Na subpopulation has only either a shift in the SSA or SSI relationships. In Fig. 13 A, we plot CV as a function of g_gap for these two cases, with the baseline/polarized and mix/homogeneous-polarized (with both SSA and SSI shifts) distributions shown for comparison. CV increases with increasing g_gap, as in previous results. Importantly, for all g_gap and cleft widths, we find that the mix/homogeneous-polarized distribution (magenta), with both SSA and SSI shifts, exhibits the largest CV, while the baseline/polarized distribution (black) exhibits the slowest CV. The mix/homogeneous-polarized distribution, with only the SSA shift, has the next fastest CV, while the case with only the SSI shift is third.

Relative to the mix/homogeneous-polarized distribution, the CV ratio for tissue with only SSA-shifted I_Na is above 0.9, nearly independent of g_gap or cleft width (Fig. 13 B). In contrast, the CV ratio for the tissue with only SSI-shifted I_Na is consistently lower, varying between 0.6 and 0.9, generally decreasing as g_gap decreases. The CV ratio for the baseline/polarized distribution is the smallest and always <1 (note: this curve is the inverse of the dashed line in Fig. 6 C). Thus, we find that, while both SSA and SSI shifts result in enhancement in conduction (relative to the baseline/polarized case), the shift in the SSA curve is primarily responsible for the increase in conduction observed in the mix/homogeneous-polarized distribution.

As a final investigation illustrating how different biophysical properties impact conduction, we measure CV for different values of $Vactshift$ and $Vinactshift$ for the ID-localized Na⁺ channels in the mix/homogeneous-polarized distribution (Fig. 14). Across all cleft widths and gap junction conductances, a negative (left) shift in SSA and positive (right) shift in SSI both promote faster conduction, with shifts in SSA modulating CV with greater sensitivity, consistent with predictions of earlier action potential upstroke time in single cells (Fig. S1). Thus, the directional shifts in both SSA and SSI curves in the physiological Na⁺ channel distribution both contribute to faster conduction. Further, we note the large CV ranges observed due to shifts in the SSA and SSI curves of up to 10 mV.

In this study, we demonstrate that an I_Na subpopulation with distinct biophysical properties, specifically negatively shifted SSI and positively shifted SSA voltage dependence, can modulate cardiac conduction. In a single cell, this I_Na subpopulation promotes an early AP upstroke, in a manner such that the shifted Na⁺ channels provide a proportionally greater contribution to the total Na⁺ charge during the upstroke. In a one-dimensional tissue model that accounts for the subcellular spatial localization of the different subpopulations, the shifted Na⁺ channels promote faster conduction. Further, the shifted Na⁺ channels result in greater CV sensitivity to changes in the intercellular cleft width, and this occurs across the entire range of physiological values for the gap junction conductance. The ID-localized shifted Na⁺ channels also provide a proportionally greater contribution to the total Na⁺ charge, in particular the postjunctional Na⁺ channels, while lateral membrane-localized Na⁺ channels contribute proportionally less, such that loss of the ID-localized Na⁺ channels results in greater conduction slowing. Additionally, we find that the ID-localized shifted Na⁺ channels enable conduction for faster pacing rates in a manner that results in a flatter CV restitution curve. Finally, we demonstrate in the computational model that the shift in the SSA curve is primarily responsible for the enhancement in conduction.

We highlight that our computational study is agnostic as to the specific “source” for the different biophysical properties between the Na⁺ channel subpopulations, which could arise due to differences in interacting and scaffolding proteins at the distinct subcellular locations (Abriel, 2010; Rivaud et al., 2017; Rivaud et al., 2020; Shy et al., 2013). We speculate that association with different β subunits, which have been shown to alter Nav1.5 biophysical properties (Angsutararux et al., 2021; Zhu et al., 2017), in different subcellular locations could also contribute to distinct Na⁺ subpopulations. Finally, there is increasing evidence for the presence of different Nav1.x isoforms in cardiac myocytes (Tarasov et al., 2023; Westenbroek et al., 2013), which could also result in subpopulations with distinct biophysical properties. Our computational study provides a framework for investigating the potential impact of different Na⁺ channel subpopulations. Here, we specifically considered two subpopulations with “baseline” and shifted biophysical properties based on measures from Lin et al. (2011); we highlight that the specific subpopulation biophysical properties were based on one specific set of conditions. However, these properties can and likely do vary for different conditions and disease settings. Collectively, our study illustrates that biophysical property differences can have important functional consequences for cardiac conduction and the contributions of the different subpopulations. Further, our approach could be naturally extended to consider multiple subpopulations with both distinct and mixed spatial localization, with biophysical properties based on interacting proteins, β subunits, and Nav isoforms. Future work will consider such complexities, in collaboration with experimental groups performing microscopy, patch clamp, and biophysical measurements.

We note that our computational study tests several predictions posed in the study by Lin et al. (2011). In their Discussion, the authors posit that Na⁺ channels in the middle of the cell “are mostly inactivated at a normal resting potential, leaving most of the burden of excitation to …[Na⁺] channels in the ID region.” Considering the mix/homogeneous-polarized physiological Na⁺ distribution, our results mostly agree with this prediction (Figs. 8 and 9). In particular, our results do predict that ID-localized Na⁺ channels provide the majority of the burden of excitation (quantified by the fraction of total Na⁺ charge), in particular the postjunctional ID Na⁺ channels. However, we predict that the lateral membrane-localized Na⁺ channels also contribute, albeit less than their respective proportion of the total Na⁺ conductance. Further, Lin and colleagues state that “[Middle Na⁺] channels have a minimum contribution to I_Na under control conditions but represent a functional reserve that can be upregulated by exogenous factors” (Lin et al., 2011). Our results also generally agree with these predictions; specifically, simulations predict that exogenous factors that reduce ID-localized I_Na result in a larger contribution from lateral membrane-localized I_Na (Fig. S6). We speculate that a concomitant upregulation of lateral membrane Na⁺ channels could further enhance conduction consistent with the concept of a functional reserve. Interestingly, the lateral membrane-localized Na⁺ channels contribute the least to excitation as gap junctional coupling is reduced, while prejunctional Na⁺ channels contribute the most for moderate values of gap junctional conductance (Fig. 10). While our study does not specifically test the hypothesis that lateral membrane Na⁺ channels are upregulated by exogenous factors under some conditions, our study does provide support for the notion of these channels serving as a functional reserve. That is, as f_Na decreases, the channels will provide a greater contribution to the burden of excitation. Additionally, while our study predicts that the loss of junctional Na⁺ channels results in greater conduction slowing, the loss of the lateral membrane-localized Na⁺ channels also slows conduction (Fig. 11 and Fig. S7), consistent with prior experimental studies (Shy et al., 2014).

An interesting prediction in this study is that ID-localized I_Na provide a proportionally larger contribution to total Na⁺ charge (Q_Na) entering the cell, relative to I_Na conductance, specifically postjunctional ID I_Na. This prediction has interesting implications for “Na⁺ transfer” at the cell–cell junction, a potential additional mechanism underlying ephaptic coupling recently demonstrated by Ivanovic and Kucera (2021). In brief, full depolarization of the prejunctional cell and hyperpolarization of the cleft result in a prejunctional ID transmembrane potential larger than the I_Na reversal potential (E_Na), which in turn results in a brief outward prejunctional I_Na. The timing of this outward prejunctional I_Na occurs during the inward postjunctional I_Na, suggesting a Na⁺ transfer from the pre- to postjunctional ID membranes across the intercellular cleft. For the parameter conditions considered in the above studies, however, we did not observe this mechanism, illustrated by the solely inward I_Na traces in Fig. 8. There are many differences between the model framework in Ivanovic and Kucera (2021) and this study that could potentially account for the absence of Na⁺ transfer in this study, including the level of structural detail represented at the ID. However, one relatively simple difference between the studies was the total I_Na conductance, which was slightly higher in Ivanovic and Kucera (23 mS/cm²). To assess if this increase in I_Na conductance also contributed to the presence of the Na⁺ transfer mechanism, we performed additional simulations in which the total I_Na conductance was increased (Figs. S10 and S11). Importantly, we find that increasing the total I_Na conductance results in the presence of prejunctional outward I_Na, consistent with the Na⁺ transfer mechanism as posited by Ivanovic and Kucera (2021). The magnitude of the outward I_Na component increases as the total I_Na increases and as the cleft width narrows (Fig. S10). Interestingly, the outward I_Na component decreases for stronger gap junctional coupling (Fig. S11). This arises due to the stronger gap junctional coupling shortening the delay between pre- and postjunctional I_Na, which results in a briefer cleft hyperpolarization, such that the prejunctional transmembrane potential is not elevated above E_Na. Importantly, the parameter conditions promoting Na⁺ transfer do not impact the general conclusions of the study, specifically that the inclusion of the shifted I_Na subpopulation at the ID still results in earlier and larger I_Na and thus faster condition compared with only baseline I_Na at the ID (Figs. S10 and S11). Future work will also consider how the Na⁺ transfer mechanism may impact the contributions of total Na⁺ charge across the different ID and lateral membranes.

One interesting observation from the simulations across all seven considered Na⁺ channel distributions is that the physiological distribution (i.e., the mix/homogeneous-polarized), in general, did not support the fastest conduction (Fig. 4); rather, the shifted/polarized distribution generally had the largest CV across different cleft widths and gap junction conductances. That is, the fastest conduction occurred for tissues with only the shifted Na⁺ channels, which naturally suggests the question of why the physiological distribution may be beneficial, compared with the shifted/polarized distribution of only shifted Na⁺ channels. First, our data demonstrate that the mix/homogeneous-polarized distribution exhibits the greatest sensitivity in response to changes in cleft width (Fig. 7). Thus, if we analogize conduction as being regulated by a series of “dials,” then gap junctional conductance is the “strongest” dial, highlighted by the wide range of CV values over the physiological range for g_gap. While in the above results, the total Na⁺ channel conductance is fixed, prior work has also shown strong sensitivity to this value (Nowak et al., 2021b), suggesting overall Na⁺ channel conductance is an additional key dial modulating conduction. Our current study suggests that additional dials include the cleft width and Na⁺ channel biophysical properties and subcellular localization. However, it is important to note that these dials are interdependent, in that the strength of one dial may be context-specific and depend on the current state of the others. For example, prior work (Kucera et al., 2002) and the current study (Fig. 4) demonstrate the CV is more weakly dependent on g_gap for narrower clefts, while here we show that CV is more sensitive to cleft width for a Na⁺ channel distribution with the shifted Na⁺ channels at the ID.

Additionally, we speculate that restricting Na⁺ channels with shifted biophysical properties to the ID and not the lateral membrane or t-tubules could potentially play an anti-arrhythmogenic protective role via several mechanisms. In the setting of dysfunctional calcium (Ca²⁺) handling, prior work has shown that spontaneous calcium release events can promote an influx of depolarizing current via the Na⁺/Ca²⁺ exchanger (NCX), which in turn can promote pro-arrhythmic delay after depolarizations (DADs; Bers et al., 2006; Bögeholz et al., 2015; Bögeholz et al., 2016). Colocalization of shifted Na⁺ channels (specifically with negatively shifted SSA) with the Ca²⁺ handling proteins in t-tubules would be more likely to activate due to the NCX current and thus drive DADs, such that restricting these shifted Na⁺ channels to the ID would be anti-arrhythmogenic. However, we note that as a recent computational study predicted that the restricted extracellular space of the t-tubules is hyperpolarized in a manner similar to the ephaptic effects described in the intercellular cleft (Vermij et al., 2019), the interactions between t-tubule-localized Na⁺ and Ca²⁺ channels and NCX may in fact be quite complex during excitation–contraction coupling. In contrast, restricting the subpopulation of shifted Na⁺ channels at the ID will, by design, result in earlier activation (Figs. 8 and S4), which is ideal for supporting robust conduction. Additionally, at faster pacing rates (and thus short DIs), prior studies have shown that a steeper CV restitution curve promotes spatially discordant alternans, which can lead to conduction block and wave breaks to initiate spiral waves and arrhythmias (Choi et al., 2007; Qu et al., 2000; Weiss et al., 2006), suggesting that the flatter CV restitution curve of the physiological distribution is also anti-arrhythmic (Fig. 12 and Fig. S9).

Finally, we note the potential limitations of our study. First, we note that shifts in SSA and SSI were measured in murine ventricular myocytes in Lin et al. (2011), while the LR1 model represents guinea pig ventricular ionic current dynamics. To our knowledge, such measurements of lateral vs. ID I_Na have only been measured in rats or mice and not in guinea pigs, while the LR1 model was chosen for its low computational cost that enabled the broad parameter investigation in this study. To account for the species differences, the shifts in SSA and SSI were implemented as relative shifts between lateral and ID-localized I_Na, instead of using the experimental measures for the voltage of half activation or inactivation (V_1/2), as noted in Materials and methods. Further, while there are species differences in I_Na kinetics, the overall dynamics of fast activation and fast/slow inactivation are generally conserved across species, such that we expect our general findings with regards to the impact of shifts in SSA and SSI on an ID-localized I_Na subpopulation to also be conserved. In Fig. S12, we provide evidence for this species’ generalizability and show that including a subpopulation with shifted I_Na in single cells in a mouse (Bondarenko et al., 2004; Morotti et al., 2014) and human (O’Hara et al., 2011) ventricular myocyte model exhibits similar results as shown in Fig. 2; that is, the shifted I_Na activates earlier and with larger magnitude, producing a larger total I_Na and earlier action potential upstroke.

Additionally, our tissue simulations represent a one-dimensional “fiber” or chain of coupled cells that does not account for the complex and heterogeneous three-dimensional geometry of the heart. This limits our ability to make predictions on the potential impact of Na⁺ channel subpopulations on transverse conduction or inherently higher dimensional electrical dynamics, such as spiral waves and arrhythmia initiation. Additionally, while we vary the distribution of Na⁺ channels and their respective biophysical properties, we do not consider the role of other ionic currents. While most other ionic currents will have minimal impact on conduction, the relative expression of the inward rectifier I_K1 current has been shown to impact excitability via regulation of the resting membrane potential (Dhamoon and Jalife, 2005).

Further, we assume that the biophysical properties and localization of the subpopulations are static; however, we speculate that these characteristics may be dynamic. Indeed, our study demonstrates that variability in these properties provides additional mechanisms for regulating conduction. Thus, our study can be considered a “snapshot” of conduction for a given Na⁺ channel distribution. We also do not consider additional degrees of heterogeneity, both within the tissue and within the ID. While to our knowledge, there are no measures of the spatial heterogeneity in Na⁺ channel distributions, it highly plausible that such heterogeneity exists within the ventricles, and we hypothesize that conduction is robust to some degree of spatial heterogeneity in these properties. This question is a focus of future work. Additionally, there is significant spatial heterogeneity within the ID, specifically within plicate and interplicate regions, in membrane structure and intermembrane separation. We have recently developed a finite-element framework to model ID structure and integrate the resulting heterogeneity in ID and cleft properties into a tissue-scale model (Moise et al., 2021). Additionally, Ivanovic and Kucera recently showed in a high-resolution model of the ID that Na⁺ channel clustering within the perinexus, a region of additional cleft narrowing adjacent to the gap junction plaque, can further promote faster cell–cell electrical transmission (Ivanovic and Kucera, 2021). However, the broad parameter investigation performed here would have been computationally prohibitive with these more structurally detailed models. Future work will investigate the role of both ion channel organization and biophysical properties within the ID on cardiac conduction.

In conclusion, we investigate the role of Na⁺ channel subpopulations with distinct biophysical properties and spatial localization on cardiac conduction in a computational model of cardiac tissue. We find that ID-localized Na⁺ channels with negatively shifted SSA and positively shifted SSI support the majority of burden of excitation and contribute to faster and more robust conduction, specifically regarding tissue structural changes (i.e., cleft width), gap junctional coupling, and fast pacing rates. Our work supports the hypothesis that Na⁺ channel redistribution and biophysical properties regulation are critical mechanisms by which cardiac cells can respond to perturbations to support electrical conduction.

Simulation code for single cells and tissue is freely available in a Github repository (https://github.com/SHWeinberg/Weinberg-JGP-INa-Distribution).

David A. Eisner served as editor.

This work was supported by the National Institutes of Health (NIH) grants R01HL138003 and R01HL165751.