Neural Interactions in Developing Rhythmogenic Spinal Networks: Insights From Computational Modeling

The mechanisms involved in generation of rhythmic locomotor activity in the mammalian spinal cord remain poorly understood. These mechanisms supposedly rely on both intrinsic properties of constituting neurons and interactions between them. A subset of Shox2 neurons was suggested to contribute to generation of spinal locomotor activity, but the possible cellular basis for rhythmic bursting in these neurons remains unknown. Ha and Dougherty (2018) recently revealed the presence of bidirectional electrical coupling between Shox2 neurons in neonatal spinal cords, which can be critically involved in neuronal synchronization and generation of populational bursting. Gap junctional connections found between functionally-related Shox2 interneurons decrease with age, possibly being replaced by increasing interactions through chemical synapses. Here, we developed a computational model of a heterogeneous population of neurons sparsely connected by electrical or/and chemical synapses and investigated the dependence of frequency of populational bursting on the type and strength of neuronal interconnections. The model proposes a mechanistic explanation that can account for the emergence of a synchronized rhythmic activity in the neuronal population and provides insights into the possible role of gap junctional coupling between Shox2 neurons in the spinal mechanisms for locomotor rhythm generation.

Introduction

The mammalian spinal cord contains neuronal circuits that can generate locomotor-like oscillations in the absence of supraspinal and afferent inputs. These circuits include rhythm-generating kernels capable of producing, maintaining, and coordinating populational rhythmic activity. A rhythmogenic kernel on each side of the cord is thought to be comprised of mutually connected excitatory neurons with ipsilateral projections (Kjaerulff and Kiehn, 1996; Kiehn, 2006; Hägglund et al., 2010). Although no single genetically identified neuron type has been found to be solely responsible for rhythm generation in the spinal cord, the genetically identified Shox2 interneurons have been suggested to be involved in locomotor rhythmogenesis (Dougherty et al., 2013). Later investigations provided evidence that a subpopulation of Shox2 neurons, specifically the Shox2 non-V2a neuron type, may belong to the rhythm-generating kernel in the neonatal rodent spinal cord and be involved in generation of locomotor rhythmic activity (Dougherty et al., 2013; Kiehn, 2016).

Although the existing data on the possible cellular basis for rhythmic bursting in the Shox2 neurons remains limited, there is evidence for the expression of several types of potentially rhythmogenic currents in these neurons, including the persistent inward current (Ha et al., 2019). Based on indirect data from multiple groups (Tazerart et al., 2007, 2008; Zhong et al., 2007; Ziskind-Conhaim et al., 2008; Brocard et al., 2010, 2013; Tong and McDearmid, 2012) a previous series of computational models, primarily constructed by Rybak and colleagues (Rybak et al., 2006a,b, 2013, 2015; McCrea and Rybak, 2007; Sherwood et al., 2011; Zhong et al., 2012; Brocard et al., 2013; Shevtsova et al., 2015; Danner et al., 2016, 2017, 2019; Shevtsova and Rybak, 2016; and Ausborn et al., 2019) suggested that the generation of locomotor rhythmic activity in the spinal cord relies on the persistent sodium current (I_NaP).

Regardless of the intrinsic bursting properties of single neurons involved in rhythmogenesis, mutual interactions between the neurons are critical to allow synchronization of neuronal activity leading to the generation of populational rhythmic activity. These excitatory interactions in the spinal cord are mediated by both chemical and electrical synapses (Tresch and Kiehn, 2000), the importance and role of which appear to be age-dependent (Walton and Navarrete, 1991; Chang et al., 1999). A recent study has shown that functionally related Shox2 interneurons in neonatal mice are locally connected with each other bidirectionally via electrical synapses or gap junctions (Ha and Dougherty, 2018). This bidirectional electrical coupling may represent a mechanism for neuronal synchronization and for populational oscillations in the neonatal spinal cord. Bidirectional gap junctional coupling within each population of genetically identified neurons, i.e., Shox2 non-V2a neurons, is found with high incidence in the neonatal spinal cord. With age, this electrical coupling decreases in both incidence and strength, which may occur in parallel (and possibly in connection) with an increase of the unilateral connections via excitatory chemical synapses (Ha and Dougherty, 2018). This change in the type, probability and strength of neuronal interactions may lead to specific changes in some characteristics of locomotor rhythm and pattern and the range of generated frequencies.

To theoretically investigate the potential roles of different (electrical and excitatory chemical) connections between Shox2 interneurons and their interaction with I_NaP in generation of rhythmic populational activity, we developed a model of a population of neurons incorporating I_NaP in a subset of cells that are sparsely connected by electrical and/or excitatory chemical synapses and investigated the dependence of frequency of populational bursting on the connection type and strength, the maximal conductance and distribution of I_NaP, and other neuronal characteristics. We demonstrated that the frequency of populational bursting activity depends on the relative expression of particular types of connections between neurons. The model proposes a mechanistic explanation for emergence of a synchronized rhythmic activity in a population of Shox2 neurons and provides insights into the mechanisms of the locomotor rhythm generation.

Materials and Methods

Experimental Methods

All experimental procedures followed NIH guidelines and were reviewed and approved by the Institutional Animal Care and Use Committee at Drexel University.

We performed further analysis of ventral root recordings during fictive locomotion (evoked with NMDA and 5-HT) prior to and after application of 100 μM carbenoxolone, collected for Figure 6 in Ha and Dougherty (2018) with detailed methods contained therein. Briefly, spinal cords were isolated from neonatal mice (P0-5) in a cold dissecting solution containing in mM: 111 NaCl, 3 KCl, 11 glucose, 25 NaHCO₃, 3.7 MgSO₄, 1.1 KH₂PO₄, and 0.25 CaCl₂, bubbled with 95%/5% O₂/CO₂. Cords were then transferred to a room temperature recording solution containing in mM: 111 NaCl, 3 KCl, 11 glucose, 25 NaHCO₃, 1.3 MgSO₄, 1.1 KH₂PO₄, and 2.5 CaCl₂, bubbled with 95%/5% O₂/CO₂. Glass suction electrodes were used to record activity from 2 to 3 lumbar (L) 2–5 ventral roots, band-pass filtered at 10–1,000 Hz. Locomotion was evoked by bath application of N-Methyl-D-aspartic acid (NMDA, 7 μM, Sigma) and serotonin creatinine sulfate monohydrate (5-HT, 8 μM, Sigma). The gap junction blocker carbenoxolone (Sigma) was applied at 100 μM. Here, 25 bursts immediately before carbenoxolone application and 30 min after the start of carbenoxolone application were subjected to analysis previously described (Dougherty et al., 2013). Briefly, ventral root activity was rectified and smoothed (time constant 0.2 s) in Spike2 (Cambridge Electronic Design). Recordings included either 2 ventral roots or 3 ventral roots. The cycle and amplitude values were taken from a rostral lumbar root (L1, L2, or L3) ipsilateral to a recorded caudal lumbar root (L4 or L5). The frequency was calculated from the mean cycle period, each measured from burst onset to onset. Amplitude was the measured peak amplitude of each burst. The mean values were normalized to the no drug condition to control for inter-experiment variation including root size, suction electrode size and seal.

Computational Methods

Neuron Model

All neurons were simulated in the Hodgkin-Huxley style as single-compartment neuron models. The membrane potential, V, in neurons with I_NaP was described by the following differential equation:

$\begin{array}{ll}C\times \frac{dV}{dt}=-I_{\text{Na}}-I_{\text{NaP}}-I_{\text{K}}-I_{\text{L}}-I_{\text{Gap}}-I_{\text{Syn}},& (1)\end{array}$

where C is the membrane capacitance and t is time.

The membrane potential, V, in neurons without I_NaP was described as follows:

$\begin{array}{ll}C\times \frac{dV}{dt}=-I_{\text{Na}}-I_{\text{K}}-I_{\text{L}}-I_{\text{Gap}}-I_{\text{Syn}}.& (2)\end{array}$

The ionic currents in Equations (1) and (2) were described as follows:

$\begin{array}{l}{I}_{\text{Na}}={ḡ}_{\text{Na}}\times m_{\text{Na}}^3\times h_{\text{Na}}\times (V-E_{\text{Na}});\\ I_{\text{NaP}}={ḡ}_{\text{NaP}}\times m_{\text{NaP}}\times h_{\text{NaP}}\times (V-E_{\text{Na}});\\ I_{\text{K}}={ḡ}_{\text{K}}\times m_{\text{K}}^4\times (V-E_{\text{K}});\\ I_{\text{L}}=g_{\text{L}}\times (V-E_{\text{L}}),& (3)\end{array}$

where I_Na is the fast sodium current with maximal conductance ḡ_Na; I_NaP is the persistent (slowly inactivating) sodium current with maximal conductance ḡ_NaP; I_K is the delayed-rectifier potassium current with maximal conductance ḡ_K; I_L is the leakage current with constant conductance g_L. E_Na, E_K, and E_L are the reversal potentials for sodium, potassium and leakage currents, respectively; m and h with indexes indicating currents are, respectively, the activation and inactivation variables of the corresponding ionic channels.

Dynamics of activation and inactivation variables of voltage-dependent ionic channels (Na, NaP, and K) in Equation (3) was described by the following differential equations:

$\begin{array}{l}{\tau }_{mi}(V)\times \frac{d}{dt}{m}_i=m_{\infty i}(V)-m_i;\\ {\tau }_{hi}(V)\times \frac{d}{dt}{h}_i=h_{\infty i}(V)-h_i,& (4)\end{array}$

where m_∞i(V) and h_∞i(V) define the voltage-dependent steady-state activation and inactivation of the channel i, respectively, and τ_mi(V) and τ_hi(V) define the corresponding time constants. Activation of the fast and persistent sodium channels was instantaneous. The expressions for channel kinetics in Equation (4) are described as follows:

$\begin{array}{l}{m}_{\infty \text{Na}}(V)={(1+exp(-(V+42.5)/6.5))}^{-1};{\tau }_{m\text{Na}}=0;\\ h_{\infty \text{Na}}(V)={(1+exp((V+65.5)/10.2))}^{-1};\\ {\tau }_{hNa}(V)=35.2/cosh((V+65.5)/12.8);\\ m_{\infty \text{NaP}}(V)={(1+exp(-(V+52)/3.2))}^{-1};{\tau }_{m\text{NaP}}=0;\\ h_{\infty \text{NaP}}(V)={(1+exp((V+57)/5))}^{-1};\\ {\tau }_{h\text{NaP}}(V)=T_{hmax}/cosh((V+57)/8);\\ m_{\infty \text{K}}(V)=(1+exp{(-(V+34.5)/5)}^{-1};\\ {\tau }_{m\text{K}}(V)=10/cosh((V+34.5)/10);h_{\text{K}}=1.& (5)\end{array}$

Parameters for channel kinetics were derived from Brocard et al. (2013). The maximal conductances for the fast sodium, delayed-rectifier potassium, leakage currents were, respectively: ḡ_Na = 80 nS, ḡ_K = 100 nS, and g_L = 1 nS. The maximal conductance and inactivation time constant for I_NaP, ḡ_NaP and T_{h max}, and the leakage reversal potential, E_L, in different simulation were varied.

The synaptic current (I_Syn with conductance g_Syn and reversal potential E_Syn) for the j-th neuron was described as follows:

$\begin{array}{ll}{I}_{\text{Syn},j}=g_{\text{Syn}}·({\displaystyle \sum _{i=1}^N}{w}_{\text{Syn},ij}·s_i)·(V_j-E_{\text{Syn}}),j=1,...,N,& (6)\end{array}$

where g_Syn = 1 nS, E_Syn = 0, and w_Syn,ij = w_Syn if the i-th and j- neurons were synaptically coupled and equal to 0 otherwise. The value of w_Syn was varied.

The synaptic variable s was governed by first order activation scheme (Compte et al., 2003):

$\begin{array}{l}\frac{d}{dt}s={\alpha }_s·s_{\infty }·(1-s)-\frac{s}{{\tau }_s},\\ s_{\infty }=1/(1+exp(-(V-V_{s12})/k_s)),& (7)\end{array}$

where α_s = 1; τ_s = 15 ms; V_s12 = −20; and k_s = 2.

The electrical coupling current I_Gap for the j-th neuron was described as follows:

$\begin{array}{ll}{I}_{\text{Gap},j}={\displaystyle \sum _{i=1}^N}{g}_{\text{Gap},ij}·(V_i-V_j),j=1,...,N,& (8)\end{array}$

where the coupling strength g_Gap,ij = g_Gap if the i-th and j-th neurons share a gap junction and is equal to 0 otherwise. The value of g_Gap was varied.

The following general neuronal parameters were assigned: C = 40 pF; E_Na = 55 mV; E_K = – 80 mV.

Two-Cell Model

Each neuron in the two-cell model represented either a bursting neuron incorporating the persistent sodium current, I_NaP, or simple tonically spiking neuron without I_NaP. Dynamics of the membrane potential of the neurons was described by Equations (1)–(5) where ḡ_NaP = 4 nS, if not otherwise stated. These neurons were coupled either bidirectionally by electrical synapses as described by Equation (8) (g_Gap = {0.05, 0.1, 0.2}) or unidirectionally by AMPA synapse [see Equations (6) and (7), w_Syn = {0.5, 1, 2}]. The initial operational regime of each neuron was defined by the value of its leakage reversal potential, E_L, indicated in Tables 1, 2.

TABLE 1

Table 1. Effect of electrical coupling on neuron activity in two-cell model.

TABLE 2

Table 2. Effect of chemical coupling on neuron activity in two-cell model.

Neuron Population

The neuron population contained 100 neurons and included two subpopulations: neurons incorporating I_NaP (subpopulation 1) and neurons non-incorporating I_NaP (subpopulation 2). In our simulations, N_NaP, number of neurons incorporating I_NaP, varied from 10 to 40 and the observed results were qualitatively similar.

Heterogeneity of neurons within the population was provided by random distributions of the baseline values of the mean conductance of the leakage channel, g_L, and the leakage reversal potentials, E_L0i (i = 1, 2) for two subpopulations, mean maximal conductance of the persistent sodium channel, ḡ_NaP, in the subpopulation of neurons incorporating I_NaP, and initial conditions for the values of membrane potential and channel kinetics variables. The baseline values of distributed parameters and initial conditions were assigned prior to simulations from their defined average values and variances using a random number generator (normal distribution), and a settling period of 10–40 s was allowed in each simulation. In our simulations we varied E_L01 = {−73, −74, −75} for subpopulation of neurons with I_NaP; E_L02 = {−68, −70, −72} for subpopulation of neurons without I_NaP; parameter ḡ_NaP was varied in the range of [2, 5]. We also run simulation with different values of the maximal time constant of the persistent sodium current inactivation, T_{h max} = {7, 9, 10} s. In simulations for the populational model shown in this paper: N_NaP = 40; E_L01 = −74 ± 14.8 mV; E_L02 = −70 ± 14 mV; g_L = 1 ± 0.2; ḡ_NaP = 4 ± 0.8; T_{h max} = 9 s, if not stated otherwise.

Random connections between neurons in the population were assigned prior to each simulation based on probability of connections (p). A random number generator was used to define the existence of each synaptic connection. Typically, we used p_Gap = 0.3 for gap junctional connections (Ha and Dougherty, 2018) and p_Syn = 0.1 (Dougherty et al., 2013) for excitatory chemical connections. Similar probabilities of excitatory synaptic interconnections were estimated for neurons in the lumbar region of adult turtles (Radosevic et al., 2019) and in the pre-Bötzinger Complex involved in the generation of rhythmic inspiratory activity (Carroll and Ramirez, 2013). After the distribution of connections was set, the connection strengths of electrical and/or excitatory chemical synapses (g_Gap and w_Syn) were altered. To estimate dependence of the model behavior on coupling strength and type when both types of connections were present on the model, we introduced two variables G and W (G = N × p_Gap × g_GaP and W = N × p_Syn × w_Syn), characterizing the average total network electrical or chemical synaptic input to each neuron in the population. To define electrical and chemical synapses in the populational model we used standard definition of these synapses by Brain2 (Stimberg et al., 2014, 2019). The coupling strengths for electrical and chemical connection (g_Gap and w_Syn or G and W) are specified for each simulation.

Computer Simulations

The simulations for the two-cell model were run by using a custom Matlab script (The Mathworks, Inc., Matlab 2020a). The simulations for the populational model were performed using the custom script written in Python and implemented in Brian2 environment (Stimberg et al., 2014, 2019). Differential equations were solved using the second order Runge-Kutta integration method. Simulation results were saved as the ASCI file containing the values of g_Gap and w_Syn and calculated mean values and standard deviation over simulation time for the period and amplitude of populational bursting. We also saved the figures showing populational activity, raster plot of neuron spiking, and voltage traces of sample neurons over simulation time to visually estimate the results of simulations. The scripts and model configuration files to create the simulations presented in the paper are available at https://github.com/RybakLab/Neural-interactions-in-spinal-networks.

Data Analysis in Computer Simulations

The results of simulations were processed by custom Matlab scripts (The Mathworks, Inc., Matlab 2020a). For each neuron in the two-cell model, period of bursting was estimated as the difference between two consecutive burst onsets and averaged for the time course of simulation. Frequency was calculated as reciprocal to the period. To assess the populational model behavior, the averaged integrated activity of neurons in the population (average number of spikes per neuron per second) was used to calculate the oscillation period. The timing of burst onsets was determined at a threshold level equal to 30% of the average populational burst amplitude in the current simulation. For some simulation with less stable populational activity we determined the threshold level by visual inspection of the results. The bursting period was defined as the duration between two consecutive burst onsets. Duration of individual simulations depended on the values of parameters g_Gap and w_Syn, and to robustly estimate average value of the oscillation period, the first 10–20 transitional cycles were omitted to allow stabilization of model variables, and the values of the bursting period and amplitude were averaged for the next 10–20 consecutive cycles. To build 2D diagrams we estimated stability of the populational bursting activity in each simulation by calculating values of the period and amplitude of populational bursts and their standard deviations. All simulation results were divided in four groups: no populational activity, low amplitude (<10 spikes per neuron per second) or highly unstable populational activity (SD of calculated period ≥50%), stable populational bursting, and tonic populational activity. For the stable populational bursting, in each simulation we calculated amplitude of populational bursts as the mean of maximal values of the averaged integrated activity of neurons during the burst. The frequency of populational oscillations was determined as reciprocal to the mean bursting period. The calculated frequency and amplitude were color-coded to build 2D diagrams in G x W parameter space. Re-initialization of the randomized parameters within the same ranges resulted in qualitatively similar simulation results.

Results

The major focus of the present study was on the specific roles of, and possible cooperation between, neural interactions (through gap junction coupling and/or excitatory chemical synapses) and persistent (slowly inactivating) sodium current (I_NaP) in neuronal synchronization and generation of populational rhythmic activity related to the locomotor-like oscillations in the mammalian spinal cord. We started by simulating a pair of neurons operating in different regimes and connected by bidirectional gap junctions or unidirectional excitatory chemical synapses in order to investigate the effect of the type and strength of connections and the presence/absence of I_NaP conductance on the activity of both neurons. Then, we considered a heterogeneous population of 100 neurons which incorporated a subpopulation of neurons containing I_NaP and hence could intrinsically generate rhythmic bursting depending on neuronal excitability. The remaining neurons were simple spiking neurons. The neurons had parameters randomly distributed over the population and were sparsely connected by electrical and excitatory chemical synapses. The probabilities of connections between neurons were set based on the data from Shox2 neurons (Dougherty et al., 2013; Ha and Dougherty, 2018). Our modeling study focused on neuronal interactions, synchronization of neuronal activity, and control of frequency of populational bursting depending on the relative expression of different types of neuron connectivity and connection strength.

Gap Junctional Coupling in Two-Cell Model

We first simulated a pair of interneurons (In1 and In2) modeled in the Hodgkin-Huxley style and coupled bidirectionally by gap junctions with conductance g_Gap. Each neuron represented either a conditionally bursting neuron with I_NaP or a simple tonically spiking neuron without I_NaP. The equations describing neuron dynamics and neuron parameters are specified in section Materials and Methods.

The baseline operating regime of each neuron (prior to coupling) was dependent on its resting membrane potential, which in turn was defined by the value of the leakage reversal potential, E_L (see Table 1). With an increase in E_L, the isolated conditionally bursting neuron changed its operating regime from silence to bursting and then to tonic spiking. The isolated “simple” neuron without I_NaP could switch from silence to tonic spiking when E_L increased. The strength of gap junctional coupling between the neurons (g_Gap) was varied. The time course of the neuron membrane potentials for various baseline operating regimes and the strength of coupling is illustrated in Figure 1. The frequencies of bursting, as well as neuron resting potentials at various values of g_Gap are shown in Table 1.

FIGURE 1

Figure 1. Effects of electrical coupling on activity of two neurons (In1 and In2). (A1–A3) Both neurons incorporate the persistent sodium current (I_NaP). In1 is in a bursting mode and In2 is silent if uncoupled (A1), operates in a bursing regime (A2), or is tonically active (A3). (B1–B3) Only In1 incorporates I_NaP and, if uncoupled, is silent (B1) or operates in a busting regime (B2,B3). In2, if uncoupled, is either silent (B1,B2) or tonically active (B3). The strength of electrical coupling increases from left to right.

In our simulations shown in Figures 1A1–A3, both neurons (In1 and In2) incorporate I_NaP. In all three panels, In1, if uncoupled, operates in a bursting mode, whereas In2 is silent (at low level of excitation, panel A1), exhibits bursting (at intermediate level of excitation, panel A2), or is tonic spiking (at a highest level of excitation, panel A3). The effect of gap junctional coupling depends on the difference between the membrane potentials of coupled neurons and the strength of coupling, and hence can be depolarizing for one neuron and hyperpolarizing for the other one. In all three cases, gap junctional coupling leads to synchronization of neuronal bursting from the burst ratio M:1, M>1 at a lower g_Gap (0.375 < g_Gap < 0.66 nS) to the ratio of 1:1 at higher g_Gap (g_Gap ≥ 0.66 nS) values. Note that coupling may reduce bursting frequency relative to the uncoupled case (see panel A1) or increase it (panels A2 and A3, see also Table 1) depending on the relative values of the resting membrane potentials of coupled neurons. This is seen in Figure 1A1 when In1 depolarizes In2 and induces bursting in this neuron with common bursting frequency at higher g_Gap. At the same time, In2 hyperpolarizes In1 resulting in a general decrease of common bursting frequency relative to the uncoupled case with increased g_Gap (see Table 1). In Figures 1A2,A3, we see the opposite situation. In these cases, both neurons operate in bursting regimes. However, the level of excitation of In1 is lower than that of In2, and prior to coupling, In2 generates bursts with a higher frequency than In1 in Figure 1A2 and is tonically active in Figure 1A3. Electrical coupling of these neurons leads to depolarization of In1 and hyperpolarization of In2. As a result, burst frequency of In1 increases while burst frequency of In2 decreases. In Figure 1A3, at lower g_Gap In2 switches from tonic to bursting mode, when coupled, with the burst ratio 1:2. Similar to the simulation shown in Figure 1A1, at higher g_Gap the two-cell model shows a synchronized bursting activity; however, the resulting frequency is higher than the initial In1 frequency (see Table 1). Hence gap junctional coupling between neurons with I_NaP leads to synchronized rhythmic bursting with bursting frequency dependent on neuronal excitability and the strength of coupling.

In simulations shown in Figures 1B1–B3, only one of the two electrically coupled neurons (In1) has I_NaP. Prior the coupling, this neuron can be silent (panel B1) or operate in a bursting regime (panels B2 and B3), whereas the other, “simple” neuron (In2), depending on the level of its excitation, can be silent (panels B1 and B2) or exhibit tonic spiking (panel B3). In all these simulations, the In2 neuron has higher resting potential if uncoupled and depolarizes In1 if the two neurons are connected by gap junctions (see Table 1). This results in the emergence of bursting activity (as in panel B1) or in an increase of the bursting frequency with increasing g_Gap (panels B2 and B3). In all three cases, In1 hyperpolarizes In2 during interburst intervals but induces high-frequency spiking in In2 during bursts, resulting in synchronized activity of the two-cell model replicating bursting activity of In1.

The two neurons in Figure 2 are similar to those shown in Figure 1B1 and only one of the two electrically coupled neurons (In1) has I_NaP. Here, instead of varying the strength of gap junctional coupling, it was set at 0.1 nS and the effects of varying the strength of I_NaP are illustrated. At ḡ_NaP below the level shown in Figure 2, both neurons remain silent. However, with increasing ḡ_NaP, bursting emerges and increases in frequency. The results demonstrate that the presence of I_NaP in one neuron can initiate synchronized bursting in electrically coupled neurons even if both neurons when uncoupled are silent.

FIGURE 2

Figure 2. Effects of increasing I_NaP conductance (ḡ_NaP) on activity of two neurons mutually connected via gap junction. In these simulations, E_L1 = −71.4 mV, E_L2 = −58.4 mV, and the strength of the electrical connection (ḡ_NaP) is 0.1 nS.

Overall, gap junctional coupling leads to equilibration of membrane potentials between the coupled neurons, regardless of the presence/absence of I_NaP. Depending on the direction of change (depolarization or hyperpolarization), this coupling will either increase or decrease bursting frequency of coupled neurons with I_NaP and bring neurons without I_NaP from silent or tonic to bursting mode in synchrony with neurons with I_NaP. In turn, I_NaP in some neurons can amplify the effect of electrical coupling of these neurons with other neurons. Finally, both the presence of I_NaP in some neurons and electrical coupling between neurons promotes neuronal synchronization and both are critical for initiation and support of populational bursting.

Excitatory Chemical Synaptic Connections in Two-Cell Model

In the case of excitatory connections by chemical synapses, we considered only unidirectional connections because reported chemical synapses between spinal interneurons have been unidirectional (Dougherty et al., 2013; Chopek et al., 2018; Ha and Dougherty, 2018). Figure 3 and Table 2 illustrate effects of excitatory chemical synaptic connections in the two-cell model when both neurons (source neuron In1 and target neuron In2) are bursting neurons incorporating I_NaP (Figures 3A1–A3), or when I_NaP is incorporated in only In1 (Figures 3B1,B2) or only In2 neuron (Figures 3C1,C2).

FIGURE 3

Figure 3. Effects of the excitatory chemical coupling on activity of two neurons (In1 and In2). (A1–A3) Both neurons incorporate the persistent sodium current (I_NaP). In1 is in a bursting mode and In2, if uncoupled, is silent (A1), or operates in a bursing regime (A2), or is tonically active (A3). (B1,B2) Only In1 incorporates I_NaP and operates in a busting regime and In2, if uncoupled, is either silent (B1) or tonically active (B2). (C1,C2) In1 does not incorporate I_NaP and and is tonically active. In2 incorporates I_NaP and if uncoupled, is either silent (C1) or operates in a bursting regime (C2). The weight of synaptic connection increases from left to right.

In Figures 3A1–A3, the target neuron (In2) that incorporates I_NaP operates in various regimes without connection: silent (Figure 3A1), bursting (Figure 3A2), and tonic spiking (Figure 3A3). The synaptic input from In1 to In2 leads to bursting activity in In2 with a burst ratio of 2:1 at a low w_Syn (0.5) or synchronized bursting of both neurons (at higher w_Syn values in Figures 3A1–A3, see also Table 2).

In Figures 3B1,B2, the target neuron In2 does not have I_NaP and when uncoupled is either silent (Figure 3B1) or generates tonic spiking (Figure 3B2). In these cases, when In2 receives rhythmic excitatory synaptic input from the source neuron In1, it follows In1 activity and generates bursts synchronously with the source. In Figure 2B2, the bursting activity in In2 occurs on the background of baseline spiking activity.

In Figures 3C1,C2, the source neuron In1 is a “simple” neuron without I_NaP that generates a sustained spiking activity providing constant excitatory synaptic drive to the target In2 neuron with I_NaP. This drive depolarizes In2 and this depolarization increases with the strength of synaptic input (w_Syn). Therefore, depending on the baseline regime of the target In2 neuron, this input results either in the emergence of bursting activity (panel C1) or the increasing of bursting frequency (panel C2, see also Table 2) in the target In2 neuron.

The above simulations have shown that chemical connections in the two-cell model promote neuronal synchronization and have variable effects on the resultant common bursting frequency.

Populational Models

To study the synchronization and frequency control in a heterogeneous population of neurons connected by electrical and/or excitatory chemical synapses, we developed a model of a population of 100 neurons consisting of neurons with the intrinsic bursting properties (based on the incorporated I_NaP) and “simple” neurons without such properties (with base ratio of 40/60%). All neurons were described in the Hodgkin-Huxley style as in the two-cell models described above. The neurons in the population were sparsely connected by electrical and/or excitatory chemical synapses with probabilities of p_Gap = 0.3 (Ha and Dougherty, 2018) and p_Syn = 0.1 (Dougherty et al., 2013; see also Carroll and Ramirez, 2013; Radosevic et al., 2019). Following Butera et al. (1999b), we chose the leakage reversal potential, E_L, the leakage conductance, g_Leak, and the maximum conductance of the persistent sodium current, ḡ_NaP, to be randomly distributed over population to provide natural heterogeneity between single neuron states. These parameters are the ones that define the basic level of neuronal excitation and their ability to intrinsically generate bursting activity (Butera et al., 1999a,b; Del Negro et al., 2002a; Koizumi and Smith, 2008). E_L, ḡ_NaP, and g_Leak for neurons were distributed normally with a priory assigned means and standard deviations (see section Materials and Methods). The major focus of our studies was on the generation of populational activity and frequency control of populational bursting for different types of connections and depending on the relative expression of different types of neuron connectivity and connection strength.

Performance of Populational Model With Electrical Coupling

To investigate the role of gap junction in populational bursting, we simulated a population of 100 neurons with bidirectional electrical coupling with coupling probability of 0.3 (Ha and Dougherty, 2018). Figure 4 shows the model performance for different values of coupling (g_Gap) in the population. In each panel, the upper row shows the averaged populational activity, the middle row shows raster plots of spikes elicited by all neurons in the population, and the bottom row shows traces of membrane potentials of several sample neurons with different excitabilities. Neurons with I_NaP are shown by red and neurons without I_NaP are shown by black. In Figure 4A, neurons in the population are uncoupled (g_Gap = 0). Because of the random distribution of neuronal parameters, neurons operate in different basic regimes: some neurons with I_NaP are silent, some generate bursting activity with different burst frequencies, and some are tonically active. Similarly, some neurons without I_NaP are silent and some demonstrate tonic spiking activity with different frequencies. The integrated neuron activity is sustained (upper row in Figure 4A). With increased electrical coupling (g_Gap = 0.033 nS, Figure 4B), most of neurons incorporating I_NaP synchronize their bursting activity and some neurons without I_NaP become involved in populational bursting. Further increase of g_Gap (g_Gap = 0.066 nS, Figure 4C) results in synchronization of bursting activity of all neurons incorporating I_NaP and involvement of the majority of neurons without I_NaP in populational bursting, which increases frequency and amplitude of populational bursting activity (see upper row in Figure 4C). However, after all neurons are involved in populational bursting, further increase of g_GaP does not lead to increasing frequency of oscillations. Instead, burst duration starts to increase and interburst intervals to decrease and eventually all neurons switch to sustained tonic activity (not shown).

FIGURE 4

Figure 4. Model performance at different values of electrical connection strength in the population (g_Gap). g_Gap = 0 in (A), 0.033 in (B), and 0.066 in (C). The probability of connection is 0.3. The upper row in each panel shows populational activity as an averaged histogram of neuron spiking [spikes/(N × s), where N = 100 is the number of neurons in population; bin = 100 ms]. The middle rows show raster plots of spikes elicited by neurons in the population. The lower rows show traces of the membrane potential of sample neurons. Forty percent of neurons in the population incorporate I_NaP and are shown by red and neurons without I_NaP are shown by black.

Interestingly, the frequency and amplitude of populational bursting depends on expression of I_NaP in neurons of the population. Figure 5 shows the dependence of the frequency and amplitude of populational bursting on the average value of the I_NaP conductance (ḡ_NaP) for three different values of N_NaP, defining the number of neurons with I_NaP in the population. As N_NaP increases, the populational bursting starts at lower values of ḡ_NaP. The burst frequency increases with increasing values of both ḡ_NaP and N_NaP (Figure 5A), and amplitude of populational bursts is higher at greater N_NaP values but depends on ḡ_NaP only at lower ḡ_NaP values (Figure 5B).

FIGURE 5

Figure 5. Dependence of bursting frequency (A) and amplitude (B) on the average value of the I_NaP conductance (ḡ_NaP) for three values of N_NaP, the number of neurons incorporating I_NaP in the population. In these simulations, g_Gap = 0.06 nS.

Gap Junctional Coupling and Locomotor-Like Activity in the Spinal Cord

Previous experimental studies have shown that blocking gap junctions with carbenoxolone reduces the frequency of ventral root bursts of locomotor-like activity in the isolated rodent cord evoked by the mixture of NMDA and 5-HT (Ha and Dougherty, 2018). Figures 6A1,A2 show an example of ventral root recordings during drug-evoked locomotion prior to and 30 min after carbenoxolone application begins. Burst frequency and amplitude are quantified in Figures 6B1,B2 demonstrating that blocking gap junctions in these preparations reduces the frequency of locomotor oscillations and, in most cases, the amplitude of integrated locomotor activity. To simulate the effect of blocking electrical synapses, we reduced the strength of gap junctional coupling between neurons by 50% in our model. Figures 6C1,C2 demonstrate the average populational activity for a sample simulation when g_Gap was reduced by 50%. To generally estimate the effect of this reduction on frequency and amplitude of synchronized bursting, we ran several simulations for populations with different randomization of neuron parameters. As can be seen in Figures 6D1,D2, regardless of parameter distribution, our model exhibits changes in the frequency and amplitude of oscillations with a decrease of g_Gap similar to experimental results shown in Figures 6A1,A2,B1,B2.

FIGURE 6

Figure 6. Blocking gap junctions decreases locomotor frequency. (A1,A2) Extracellular recordings from ventral roots at lumbar level 2 (L2), flexor dominant root) and level 5 (L5), extensor dominant root) on the right (r) and left (l) sides of the spinal cord after application of NMDA (7 mM) and serotonin (8 mM) in control case and after addition of carbenoxolone (100 μM). (B1,B2) Effect of carbenoxolone application on frequency (B1) and normalized amplitude (B2) of locomotor activity for 6 spinal cord preparations. Each circle represents the mean value of estimated parameter for one recording. Thinner lines connecting two circles represent parameter change in one cord. Thick red lines represent the means. The black lines in (B1,B2) correspond to the experimental recordings shown in (A1,A2). (C1,C2) Effect of decreased strength of electrical coupling (g_Gap) on populational activity. In (C2), g_Gap is reduced by 50%. (D1,D2) Effect of g_Gap reduction on frequency (D1) and normalized amplitude (D2) of locomotor activity for 5 simulations with different randomizations of model parameters (simulation time 60 s). Each circle represents the average value of estimated parameter for one simulation. Light colored lines connecting two circles represent parameter change in simulation. Thick red lines represent the means. The black lines in (D1,D2) correspond to the simulations shown in (C1,C2).

Performance of Populational Model With Chemical Coupling

In the next stage of this study, we considered the same population of 100 neurons, but replaced the bidirectional gap junction connections within the population with unilateral excitatory chemical connections with probability of p_Syn = 0.1 (Carroll and Ramirez, 2013; Dougherty et al., 2013; Radosevic et al., 2019). Figure 7 shows the model performance for different weights of excitatory synaptic connections (w_Syn) in the population. As can be seen in these simulations, enhancing w_Syn increases the number of tonically active neurons both among neurons with and without I_NaP (Figures 7A,B). Simultaneously, the bursting neurons start to consolidate in clusters including synchronously active neurons with I_NaP and partially involving neurons without I_NaP in bursting activity. Further increase of w_Syn increases frequency of bursting in neurons with I_NaP and synchronizes all neurons, resulting in populational synchronous bursting of high frequency and amplitude. In contrast to the population of electrically coupled neurons, in the population with chemical coupling, an increase in w_Syn is accompanied by a simultaneous increase of the background tonic activity.

FIGURE 7

Figure 7. Model performance at different values of excitatory chemical connection strength in the population w_Syn = 0 in (A), 1 in (B), and 1.6 in (C). The probability of connection is 0.1. Arrangement and notation are the same as in Figure 4.

Contribution of Neuronal Coupling by Electrical and Chemical Synapses to Populational Activity

Experimental studies have shown that electrical coupling between Shox2 neurons is age-dependent. It is prevalent between functionally related Shox2 neurons in neonatal spinal cords and decreases in incidence and strength with age (Ha and Dougherty, 2018). To investigate how the activity of the population changes when both electrical and excitatory chemical connections are varied, we simulated a population in which both electrical and excitatory chemical synapses are present and sparsely distributed among neurons with probabilities of p_Gap = 0.3 and p_Syn = 0.1, and introduced two variables G and W, characterizing in average the total network electrical or chemical synaptic input to each neuron in the population, respectively (see Materials and Methods).

For each pair (G,W), we estimated stability of the populational bursting activity by calculating values of the period and amplitude of populational bursts and their standard deviations. All simulation results were divided in four groups: no populational activity, low amplitude or highly unstable populational activity, stable populational bursting, and tonic populational activity (see Materials and Methods). Figures 8A,B show regimes of populational activity and 2D diagrams of color-coded frequency and amplitude of stable populational bursting in [G, W] plane. As seen in Figure 8A, if W is fixed and G is increasing, for most values of W, frequency of populational bursting does not increase smoothly but there are distinct areas in which the frequency does not change much and transition points when frequency increases abruptly. For example, at lower values of W (W < 7) there are one or two distinct transitions from the dark blue to light blue to cyan or the dark blue to cyan color in the diagram. This frequency transition is even more sharp at intermediate values of W (i.e., W = ~10): dark blue abruptly changes to yellow. Our simulations have shown that in the areas in which frequency does not change much, the increasing G leads predominantly to synchronization of bursting neurons with I_NaP and involving in the populational bursting some neurons without I_NaP. However, at some higher values of G, the remaining tonically active neurons depolarize neurons with I_NaP that have not been initially involved in bursting, leading to abrupt transitions to a higher frequency. Then again increasing G results predominantly in synchronization of bursting neurons, this time at a higher frequency of populational bursting. Depending on particular distribution of neuronal parameters these transitions might happen up to three times when W is fixed and G is increasing. After all neurons are involved in populational bursting and synchronized, further increase of G does not produce increases in frequency of populational bursting. On the contrary, at higher G, the length of the populational bursts is increasing and the interburst interval is decreasing, leading to smooth transition to the slower populational bursting (shown by smooth transition from the cyan/yellow to blue in Figure 8A) and eventually the population switches to sustained tonic activity (shown for higher W values by black in Figures 8A,B).

FIGURE 8

Figure 8. Role of electrical vs. chemical coupling in generation of synchronized populational activity. (A,B) 2D diagrams showing the regime of populational activity and color-coded frequency (A) and amplitude (B) of stable populational bursting when the averaged total electrical (G) or chemical (W) synaptic input to each neuron in the population was varied. Simulation results are divided in four groups: no populational activity (white areas), low amplitude or highly unstable populational activity (gray), stable populational bursting (colors), and tonic populational activity (black). For the stable populational bursting, for each (G, W), the amplitude of populational bursting was determined as the mean of maximal values of the averaged integrated activity of neurons (average number of spikes per neuron per second) during the bursts. The frequency of stable populational oscillations was determined as reciprocal to the mean bursting period. (C1–C3) Examples of simulations for three pairs of G and W values indicated in (A,B) by small red circles. Arrangement and notation in (C1–C3) are the same as in Figure 4.

With the increase of W, the unilateral excitatory chemical connections allow for a higher bursting frequency. If W is high enough, (if W > 12 in our simulations) the model is capable of generating synchronous populational bursting at low values of G or even if electrical coupling in the population is absent. With a fixed G, the frequency of populational bursting is increasing with an increase of the total strength of the excitatory chemical connections (W). However, the range of G in which the model generates populational bursting activity decreases with increasing W due to the earlier transition to sustained tonic activity of the majority neurons (black area in Figures 8A,B). Figure 8B shows that the average changes in burst amplitude depending on the values G and W are quite smooth because the amplitude depends in a greater extent on how well-neuron activity in the population is synchronized and neuronal synchronization directly depends on increasing total strengths of the electrical and/or excitatory chemical connections (see also Figures 4, 7).

Figures 8C1–C3 demonstrate the results of simulations for three pairs of G and W values shown by small red circles in Figures 8A,B. In Figure 8C1, W = 0 and G = 6 nS. In this case, G is high enough to generate synchronized populational bursting activity in the absence of excitatory chemical connections. In Figure 8C2, G = 3 nS and at that strength of electrical coupling the model cannot generate synchronized populational bursting in the absence of excitatory chemical connections. However, populational bursting can be re-established with increasing W, which is illustrated in Figure 8C2 for W = 7.5. On the other hand, if W is high enough, populational bursting activity emerges even without electrical coupling of neurons in the population as can be seen in Figure 8C3 at W = 15 and G = 0 nS. It should be noted that in the latter case, a higher level of background activity is present.

Altogether, our simulations indicate that both electrical and excitatory chemical connections between neurons in the heterogeneous population incorporating a subset of intrinsically bursting neurons support synchronization of activities of individual neurons resulting in the rhythmic populational bursting activity. However, the effects of increasing strength of electrical and excitatory chemical coupling of frequency of populational bursting are different.

Discussion

This study is based on a previous, experimental investigation by Ha and Dougherty (2018) that revealed the presence of electrical coupling between Shox2 non-V2a interneurons, which have been suggested to be critically involved in neuronal synchronization and generation of populational bursting in the neonatal rodent spinal cords. Here, we used computational modeling to investigate the potential roles of different neural interactions (gap junctional coupling and/or excitatory chemical synapses) and their interactions with the cellular I_NaP-dependent bursting properties in the neuronal synchronization and generation of locomotor-like oscillations. We simulated and analyzed the behavior of pair of neurons operating in different regimes as well as heterogeneous populations of 100 neurons in different states, sparsely connected by either gap junctions or excitatory chemical synapses, and studied the possible role of electrical and chemical synaptic connections in populational bursting and control of frequency and amplitude of populational oscillations. The different aspects of our simulation results in connection with previous experimental findings as well as limitations of current studies are discussed below.

The Role of Electrical Coupling Between Spinal Interneurons in Rhythmic Activity in the Isolated Spinal Cord of Neonatal Rodents

Electrical coupling between neurons is known to be widely present in different regions of the central nervous system of mammals, including rodent spinal cord (Tresch and Kiehn, 2000), especially in the early stages of the development (Walton and Navarrete, 1991; Chang et al., 1999). The role of gap junctions in neuronal synchronizations and generation of populational oscillations has been previously studied using computational neuron models of different complexities (Kepler et al., 1990; Sherman and Rinzel, 1992; Pfeuty et al., 2003; Kopell and Ermentrout, 2004). Gap junctional coupling has been demonstrated to interconnect spinal Shox2 non-V2a neurons (Ha and Dougherty, 2018) as well as Hb9 neurons (via other unidentified neurons, Wilson et al., 2007). It has been suggested that electrical coupling within each of these genetically identified populations significantly contributes to neuronal synchronization and populational bursting during spinal circuit development. Although the role of this coupling in the generation of locomotor-like oscillations in neonatal mouse cords has not been directly demonstrated for the specific spinal neuronal populations, pharmacological suppression of gap junctions in the cord in general resulted in a decrease of frequency and amplitude of locomotor-like oscillations recorded from ventral roots of the neonatal spinal cord during drug-evoked fictive locomotion (Ha and Dougherty, 2018). Our simulations of a heterogeneous population of neurons connected with gap junctions have demonstrated the dependence of frequency and amplitude of populational bursting on the strength of electrical coupling and the reduction of both characteristics with the reduction of this coupling, similar to the experimental results of pharmacological suppression of gap junctions. This provides additional support for the role of electrical coupling in the locomotor activity in the spinal cords of neonatal rodents.

The Potential Role of Persistent Sodium Current in Locomotor-Like Oscillations in the Cord and Its Interaction With Gap Junction Neuronal Coupling

Despite many years of intensive investigations, the role of intrinsic biophysical neuronal properties in the generation of locomotor-like oscillations in the mammalian spinal cords remains largely unknown. Previous studies of rhythmic oscillations in the medullary pre-Bötzinger complex (pre-BötC) have suggested the important role of I_NaP-dependent pacemaker properties in the generation of these oscillations (Butera et al., 1999a,b; Del Negro et al., 2001, 2002a; Smith et al., 2007, 2013; Koizumi and Smith, 2008). The critical role of I_NaP-dependent pacemaker oscillations in the pre-BötC has been challenged by some studies (Del Negro et al., 2002b; Feldman et al., 2013) and currently remains debated (Del Negro et al., 2018, but Jasinski et al., 2013; Bacak et al., 2016; Molkov et al., 2016; Ausborn et al., 2018; Phillips and Rubin, 2019). Based on the previous computational models on the generation of respiratory oscillations in the pre-BötC in vitro (Butera et al., 1999a,b; Rybak et al., 2004), a series of computational models of spinal locomotor circuits has been developed by suggesting the critical role of I_NaP in the generation of locomotor rhythmic activity in the spinal cord (Rybak et al., 2006a,b; McCrea and Rybak, 2007). This suggestion has been partly supported by several independent studies confirming the presence of I_NaP in the neonatal rodent spinal cords (Tazerart et al., 2007, 2008; Zhong et al., 2007; Ziskind-Conhaim et al., 2008; Brocard et al., 2010, 2013; Harris-Warrick, 2010; Tong and McDearmid, 2012). Consequently, our recent computational models of spinal locomotor circuits which incorporated I_NaP considered this current to be critically involved in the spinal cord rhythmogenesis (Zhong et al., 2012; Brocard et al., 2013; Rybak et al., 2013, 2015; Shevtsova et al., 2015; Danner et al., 2016, 2017, 2019; Shevtsova and Rybak, 2016; Ausborn et al., 2019). The persistent sodium channels, including their voltage-dependent activation, have been characterized in spinal Hb9 interneurons present in the spinal cord (Brocard et al., 2013) and presumably contribute to the generation and/or propagation of locomotor-like oscillations (Wilson et al., 2007; Ziskind-Conhaim et al., 2008; Brocard et al., 2013). Although I_NaP has not yet been characterized in Shox2 neurons, there is evidence for the expression of persistent inward currents in these neurons which likely includes I_NaP (Ha et al., 2019). In this study, we followed the previous computational models of the spinal locomotor circuits developed by Rybak et al. by suggesting that I_NaP plays critical role in the generation of locomotor-like oscillations. We, however, suggested that this current is present in, and randomly distributed over, only a subset of all neurons involved in generation of populational oscillations (with basic ratio of 40/60%). Our particular interest was in considering the behavior of neurons with and without I_NaP connected by gap junctions and the possible roles of this current and electrical coupling in initiation of synchronized bursting activity. Our modeling results show that electrical coupling can lead to equilibration of membrane potentials between the coupled neurons. Depending on the membrane potentials of neurons, the electrical coupling can increase and decrease bursting frequency and bring neurons without I_NaP from silent or tonic to bursting mode, hence amplifying neuronal synchronization. We have also shown that the presence of I_NaP in coupled neurons can amplify the effect of electrical coupling, leading to emergence of bursting activity in initially silent neurons (see Figure 2). This corresponds to conclusion made from other studies that I_NaP is a strong modulator of electrical synapses (Haas and Landisman, 2012). Finally, we have shown that the presence of I_NaP in some neurons in cooperation with electrical and/or chemical coupling between neurons promotes neuronal synchronization in heterogenous neuronal population and that this current can be critical for initiation and support of populational bursting. Since incorporation of I_NaP into only a subset of neurons in the modeled population was enough to generate a synchronized populational bursting when bidirectional connections at incidences seen electrophysiologically were incorporated, we would predict that to generate rhythmic locomotor-like oscillations in the neonatal spinal cord, it is sufficient for I_NaP to be present in a subset of rhythm generating neurons but not necessarily in all of them.

Gap Junctions vs. Chemical Synaptic Connections

In isolated cords of newborn mice, suppressing chemical synaptic interactions between neurons with lower extracellular Ca²⁺ concentrations or blocking action potentials using TTX did not eliminate rhythmic motor output (Tresch and Kiehn, 2000), hence confirming the important role of electrical coupling in neonatal cords. Similarly, in our model, oscillations are generated in the population without chemical synaptic interactions (at W = 0). Lowering extracellular Ca²⁺ concentrations not only suppresses chemical synaptic interactions but has also been shown to amplify I_NaP-dependent oscillations in spinal neurons (Brocard et al., 2013). Although not simulated by the present simplified model, this may serve as additional indirect evidence for the role of I_NaP in oscillations generated by networks of electrically coupled neurons in the neonatal spinal cord. This, of course, can only be the case under conditions that I_NaP in the spinal cord is TTX-resistant (see Cummins et al., 1999; Herzog et al., 2001; Dai and Jordan, 2011). These issues remain for future study.

Electrical transmission has been shown to be prevalent during early postnatal period but decreases as the animal matures (Walton and Navarrete, 1991; Chang et al., 1999). Although the decline is seen within the first postnatal week in motor neurons (Walton and Navarrete, 1991; Chang et al., 1999), gap junctional coupling within both Shox2 non-V2a neuron and Hb9 neuron populations remains relatively high throughout the first two postnatal weeks (Hinckley and Ziskind-Conhaim, 2006; Ha and Dougherty, 2018). The probability of coupling among neighboring neurons was found to be ~0.3 (Ha and Dougherty, 2018) and this is the value that was chosen for the modeling study. The probability value of 0.1 was set for chemical synapses and falls more in line with a prior Shox2 connectivity study (Dougherty et al., 2013). There were a few key differences between the reports of Shox2 connectivity. First, the experiments showing connectivity by chemical synapses at a rate of ~10% (Dougherty et al., 2013) was performed in mice where the Chx10 absence/presence was not distinguished, so it included the entire Shox2 population as a single sampling group. Connections were tested in a dorsal horn removed preparation and were typically within the field of view of a 60x objective but not necessarily neighboring neurons and there were no instances of electrical connections reported. In Ha and Dougherty (2018), pairs of neurons recorded were in close proximity (~50 μm on average) and were selected based on the visualization of processes crossing and in close apposition. Further, these neurons were distinguished based on the presence/absence of Chx10 expression. Therefore, the 33% incidence of electrical synapses is only among very close neighboring neurons of the same group. In this study, chemical connections were only found at 2%, indicating that neighboring neurons are less likely to be connected by chemical synapses. Both values of chemical connections are likely to be underestimated, particularly the one in which only near neighboring neurons were recorded. Interestingly, the discrepancies in the results of these studies reinforces the idea that electrical connectivity is mechanism for local synchronization.

Gap junctions behave as low pass filters and efficiently convey membrane potential fluctuations. Additionally, action potentials in one Shox2 neuron reliably lead to spikelets in gap junctionally connected neurons. However, chemical synapses between neonatal Shox2 neurons had high failure rates (Ha and Dougherty, 2018). Therefore, it is possible that there is a higher incidence of chemical synaptic connections in older ages, not due to new synapses but increased efficacy at synapses established early in development. The experimental study considered single pairs using dual recordings and responses to current injection in the presynaptic neuron. It is unlikely that the postsynaptic cell receives input from a single presynaptic Shox2 neuron. Gap junctional coupling between multiple presynaptic neurons will synchronize their firing, as demonstrated in the current study, likely leading to a more reliable excitation of the postsynaptic neuron, even if excitatory chemical synapses are sparse and single synapses have high synaptic failure rates. This could provide a mechanism for increased efficacy of excitatory transmission not only within the population but in transmission to downstream targets.

Using our model, we investigated how the activity of the population might change with changing electrical and chemical synaptic connections. Our simulations demonstrated that both electrical and excitatory chemical connections between neurons in the heterogeneous population support synchronization of activities of individual neurons resulting in the rhythmic populational bursting activity. This implicitly confirms that the locomotor-like oscillations in spinal cord can be maintained over a continuum. Thus, the results suggest that gap junctional coupling is capable of strongly driving synchronization of the population. As bidirectional electrical coupling declines with age, other mechanisms can take over to maintain rhythmic bursting within the population. Our simulations demonstrate that the change to a dominance of unidirectional chemical synaptic interactions, particularly with I_NaP present in a subset of these neurons, is sufficient to fulfill this role and to expand the range of attainable frequencies. Interestingly, our simulation results have demonstrated that a strong, synchronous population bursting can be produced in a population when as low as 10% of neurons are connected by excitatory chemical synapses.

Our simulations have also demonstrated that the effects of increasing the strength of electrical and excitatory chemical coupling on the frequency of populational bursting are different. In both cases, a subset of neurons in the population generates spiking activity and provides an excitatory tonic input to other neurons which is enhanced with increasing connection strength and results in an increased frequency of populational bursting. In addition, in the case of bidirectional electrical coupling, tonically active neurons directly depolarize other neurons, elevating neuronal resting potentials and leading to an increase in frequency of populational bursting. However, in the case of bidirectional electrical coupling, bursting neurons have a hyperpolarizing effect on tonically firing neurons and eventually an increase of g_Gap turns spiking activity of majority of tonically firing neurons into bursting, as was shown in our simulations (see Figures 1A3,B3, 4) and in earlier studies (Sherman and Rinzel, 1992). This prevents further increase in bursting frequency. After all of the neurons become involved in populational bursting and their membrane potentials reach equilibrium, further increase of g_Gap does not increase the frequency of populational oscillations (see also Kepler et al., 1990). Instead, the frequency starts to decrease. In contrast, in populations with chemical coupling, elevation of w_Syn increases the activity of tonically active neurons which depolarize bursting neurons, hence leading to the increase the frequency of populational bursting.

Limitations of the Present Studies and Future Directions

In this modeling study we limited single neuron models by incorporating only minimal number of ionic channels necessary for generation of neuronal spiking activity, such as fast sodium (Na) and potassium rectifier (K), and persistent sodium (NaP) involved in bursting rhythmic activity. Therefore, these models remain largely generic. At the same time, we know that spinal neurons contain other ionic channels, including T-type calcium and h currents and different Ca²⁺ and Ca²⁺-dependent potassium currents that can also support pacemaker properties and have been implicated in locomotor-like activity (Wilson et al., 2005; Anderson et al., 2012; Kiehn, 2016; Brocard, 2019). In addition, there is indirect evidence of an important role of Na⁺/K⁺ pump in CPG operation (Kueh et al., 2016) and particularly in the neonatal spinal cords of mice (Picton et al., 2017). The Na⁺/K⁺ pump current has been previously included in the computational models of Hb9 neurons (Brocard et al., 2013) and can be included in Shox2 neuron models in the future. This current is known to interact with h-current (Kueh et al., 2016) and can specifically affect neuronal coupling via electrical and chemical synapses. This awaits additional experimental and computational investigations. A more realistic model of different types of spinal interneurons that includes different ionic channels with experimentally measured characteristics should be developed and investigated (see for e.g., Sharples et al., 2020), which represents one of the directions of our future studies.

The other limitation of the present study is the consideration of only excitatory interactions within the spinal network. The potential role of electrical coupling between inhibitory interneurons and the role of inhibitory network interactions in cooperation with gap junction (Bennett, 1997; Bennett and Zukin, 2004; Kopell and Ermentrout, 2004; Connors, 2017) will also be a focus of our future investigations.

The contribution of gap junctional coupling within specific populations to locomotor rhythmicity remains to be tested experimentally. The pharmacological studies are limiting in that the gap junction blockers are acting at all electrical synapses, which have been shown between several classes of interneurons (Hinckley and Ziskind-Conhaim, 2006; Wilson et al., 2007; Zhong et al., 2010; Chopek et al., 2018) and motor neurons (i.e., Walton and Navarrete, 1991; Rash et al., 1996; Chang et al., 1999; Tresch and Kiehn, 2000). Additionally, carbenoxolone, the most common gap junction blocker, has been shown to have additional actions which may affect rhythmicity (Elsen et al., 2008; Tovar et al., 2009; Connors, 2012). Therefore, conditional knockout studies are required to determine implications more precisely.

Finally, the mutual interactions between different classes of genetically identified excitatory and inhibitory interneurons present in the mammalian spinal cord also awaits further experimental and computational investigations.

Data Availability Statement

The modeling datasets for this study are publically available at https://github.com/RybakLab/Neural-interactions-in-spinal-networks.

Ethics Statement

This animal study was reviewed and approved by Drexel University Institutional Animal Care and Use Committee.

Author Contributions

NS, NH, IR, and KD: conceptualization, investigation, and writing (review and editing). NS, IR, and KD: writing (original draft). NH and KD: experiments and data curation. NS: modeling, formal analysis, and software. NS and KD: visualization. IR and KD: supervision, project administration, and funding acquisition. All authors contributed to the article and approved the submitted version.

Funding

This study was supported by grants from the National Institutes of Health (R01 NS090919; R01 NS095366; R01 NS100928; and R01 NS110550).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We are grateful to Simon Danner and Lihua Yao for technical assistance. This manuscript has been released as a pre-print at Shevtsova et al. (2020).

References

Anderson, T. M., Abbinanti, M. D., Peck, J. H., Gilmour, M., Brownstone, R. M., and Masino, M. A. (2012). Low-threshold calcium currents contribute to locomotor-like activity in neonatal mice. J. Neurophysiol. 107, 103–113. doi: 10.1152/jn.00583.2011

PubMed Abstract | CrossRef Full Text | Google Scholar

Ausborn, J., Koizumi, H., Barnett, W. H., John, T. T., Zhang, R., Molkov, Y. I., et al. (2018). Organization of the core respiratory network: insights from optogenetic and modeling studies. PLoS Comput. Biol. 14:e1006148. doi: 10.1371/journal.pcbi.1006148

PubMed Abstract | CrossRef Full Text | Google Scholar

Ausborn, J., Shevtsova, N. A., Caggiano, V., Danner, S. M., and Rybak, I. A. (2019). Computational modeling of brainstem circuits controlling locomotor frequency and gait. eLife 8:14. doi: 10.7554/eLife.43587.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Bacak, B. J., Kim, T., Smith, J. C., Rubin, J. E., and Rybak, I. A. (2016). Mixed-mode oscillations and population bursting in the pre-Bötzinger complex. eLife 5:e13403. doi: 10.7554/eLife.13403

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, M. V. (1997). Gap junctions as electrical synapses. J. Neurocytol. 26, 349–366. doi: 10.1023/A:1018560803261

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, M. V., and Zukin, R. S. (2004). Electrical coupling and neuronal synchronization in the mammalian brain. Neuron 41, 495–511. doi: 10.1016/S0896-6273(04)00043-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Brocard, F. (2019). New channel lineup in spinal circuits governing locomotion. Curr. Opin. Physiol. 8, 14–22. doi: 10.1016/j.cophys.2018.11.009

CrossRef Full Text | Google Scholar

Brocard, F., Shevtsova, N. A., Bouhadfane, M., Tazerart, S., Heinemann, U., Rybak, I. A., et al. (2013). Activity-dependent changes in extracellular Ca²⁺ and K⁺ reveal pacemakers in the spinal locomotor-related network. Neuron 77, 1047–1054. doi: 10.1016/j.neuron.2013.01.026

PubMed Abstract | CrossRef Full Text | Google Scholar

Brocard, F., Tazerart, S., and Vinay, L. (2010). Do pacemakers drive the central pattern generator for locomotion in mammals? Neuroscientist 16, 139–155. doi: 10.1177/1073858409346339

PubMed Abstract | CrossRef Full Text | Google Scholar

Butera, R. J., Rinzel, J. R., and Smith, J. C. (1999a). Models of respiratory rhythm generation in the pre-Bötzinger complex. I. bursting pacemaker neurons. J. Neurophysiol. 82, 382–397. doi: 10.1152/jn.1999.82.1.382

PubMed Abstract | CrossRef Full Text | Google Scholar

Butera, R. J., Rinzel, J. R., and Smith, J. C. (1999b). Models of respiratory rhythm generation in the pre-Bötzinger complex: II. populations of coupled pacemaker neurons. J. Neurophysiol. 82, 398–415. doi: 10.1152/jn.1999.82.1.398

PubMed Abstract | CrossRef Full Text | Google Scholar

Carroll, M. S., and Ramirez, J. M. (2013). Cycle-by-cycle assembly of respiratory network activity is dynamic and stochastic. J. Neurophysiol. 109, 296–305. doi: 10.1152/jn.00830.2011

PubMed Abstract | CrossRef Full Text | Google Scholar

Chang, Q., Gonzalez, M., Pinter, M. J., and Balice-Gordon, R. J. (1999). Gap junctional coupling and patterns of connexin expression among neonatal rat lumbar spinal motor neurons. J. Neurosci. 19, 10813–10828. doi: 10.1523/JNEUROSCI.19-24-10813.1999

PubMed Abstract | CrossRef Full Text | Google Scholar

Chopek, J. W., Nascimento, F., Beato, M., Brownstone, R. M., and Zhang, Y. (2018). Sub-populations of spinal V3 interneurons form focal modules of layered pre-motor microcircuits. Cell Rep. 25, 146–156. doi: 10.1016/j.celrep.2018.08.095

PubMed Abstract | CrossRef Full Text | Google Scholar

Compte, A., Sanchez-Vives, M. V., McCormick, D. A., and Wang, X. J. (2003). Cellular and network mechanisms of slow oscillatory activity (<1 Hz) and wave propagations in a cortical network model. J. Neurophysiol. 89, 2707–2725. doi: 10.1152/jn.00845.2002

PubMed Abstract | CrossRef Full Text | Google Scholar

Connors, B. W. (2012). Tales of a dirty drug: Carbenoxolone, gap junctions, and seizures Epilepsy Curr. 12, 66–68. doi: 10.5698/1535-7511-12.2.66

PubMed Abstract | CrossRef Full Text | Google Scholar

Connors, B. W. (2017). Synchrony and so much more: Diverse roles for electrical synapses in neural circuits. Dev. Neurobiol. 77, 610–624. doi: 10.1002/dneu.22493

PubMed Abstract | CrossRef Full Text | Google Scholar

Cummins, T. R., Dib-Hajj, S. D., Black, J. A., Akopian, A. N., Wood, J. N., and Waxman, S. G. (1999). A novel persistent tetrodotoxin-resistant sodium current in SNS-null and wild-type small primarys neurons. J. Neurosci. 19:RC43. doi: 10.1523/JNEUROSCI.19-24-j0001.1999

PubMed Abstract | CrossRef Full Text | Google Scholar

Dai, Y., and Jordan, L. M. (2011). Tetrodotoxin-, dihydropyridine-, and riluzole resistant persistent inward current: novel sodium channels in rodent spinal neurons. J. Neurophysiol. 106, 1322–1340. doi: 10.1152/jn.00918.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

Danner, S. M., Shevtsova, N. A., Frigon, A., and Rybak, I. A. (2017). Computational modeling of spinal circuits controlling limb coordination and gaits in quadrupeds. eLife 6:e31050. doi: 10.7554/eLife.31050.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Danner, S. M., Wilshin, S. D., Shevtsova, N. A., and Rybak, I. A. (2016). Central control of interlimb coordination and speed- dependent gait expression in quadrupeds. J. Physiol. 594, 6947–6967. doi: 10.1113/JP272787

PubMed Abstract | CrossRef Full Text | Google Scholar

Danner, S. M., Zhang, H., Shevtsova, N. A., Borowska-Fielding, J., Deska-Gauthier, D., Rybak, I. A., et al. (2019). Modeling spinal V3 interneurons and left-right coordination in mammalian locomotion. Front. Cel. Neurosci. 13:516. doi: 10.3389/fncel.2019.00516

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Negro, C. A., Johnson, S. M., Butera, R. J., and Smith, J. C. (2001). Models of respiratory rhythm generation in the pre-Bötzinger complex. III. experimental tests of model predictions. J. Neurophysiol. 86, 59–74. doi: 10.1152/jn.2001.86.1.59

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Negro, C. A., Koshiya, N., Butera, R. J., and Smith, J. C. (2002a). Persistent sodium current, membrane properties and bursting behavior of pre-Botzinger complex inspiratory neurons in vitro. J. Neurophysiol. 88, 2242–2250. doi: 10.1152/jn.00081.2002

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Negro, C. A., Morgado-Valle, C., and Feldman, J. L. (2002b). Respiratory rhythm: an emergent network property? Neuron 34, 821–830. doi: 10.1016/S0896-6273(02)00712-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Del Negro, C. A., Funk, G. D., and Feldman, J. L. (2018). Breathing matters. Nat. Rev. Neurosci 19, 351–367. doi: 10.1038/s41583-018-0003-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Dougherty, K. J., Zagoraiou, L., Satoh, D., Rozani, I., Doobar, S., Arber, S., et al. (2013). Locomotor rhythm generation linked to the output of spinal shox2 excitatory interneurons. Neuron 80, 920–933. doi: 10.1016/j.neuron.2013.08.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Elsen, F. P., Shields, E. J., Roe, M. T., VanDam, R. J., Jonathan, D., and Kelty, J. D. (2008). Carbenoxolone induced depression of rhythmogenesis in the pre-Bötzinger complex. BMC Neurosci. 9:46. doi: 10.1186/1471-2202-9-46

PubMed Abstract | CrossRef Full Text | Google Scholar

Feldman, J. L., Del Negro, C. A., and Gray, P. A. (2013). Understanding the rhythm of breathing: so near, yet so far. Annu. Rev. Physiol. 75, 423–452. doi: 10.1146/annurev-physiol-040510-130049

PubMed Abstract | CrossRef Full Text | Google Scholar

Ha, N. T., and Dougherty, K. J. (2018). Spinal Shox2 interneuron interconnectivity related to function and development. eLife 7:e42519. doi: 10.7554/eLife.42519.023

PubMed Abstract | CrossRef Full Text | Google Scholar

Ha, N. T., Shevtsova, N. A., Rybak, I. A., and Dougherty, K. J. (2019). Data driven modeling of the rhythmogenic potential of spinal Shox2 neurons. Society for Neuroscience Annual Meeting Abstract No. 12112. Chicago, IL.

Haas, J. S., and Landisman, C. E. (2012). Bursts modify electrical synaptic strength. Brain Res. 1487, 140–149. doi: 10.1016/j.brainres.2012.05.061

PubMed Abstract | CrossRef Full Text | Google Scholar

Hägglund, M., Borgius, L., Dougherty, K. J., and Kiehn, O. (2010). Activation of groups of excitatory neurons in the mammalian spinal cord or hindbrain evokes locomotion. Nat. Neurosci. 13, 246–252. doi: 10.1038/nn.2482

PubMed Abstract | CrossRef Full Text | Google Scholar

Harris-Warrick, R. M. (2010). General principles of rhythmogenesis in central pattern generator networks. Prog. Brain Res. 187, 213–222. doi: 10.1016/B978-0-444-53613-6.00014-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Herzog, R. I., Cummins, T. R., and Waxman, S. G. (2001). Persistent TTX-resistant Na1 current affects resting potential and response to depolarization in simulated spinal sensory neurons. J. Neurophysiol. 86, 1351–1364. doi: 10.1152/jn.2001.86.3.1351

PubMed Abstract | CrossRef Full Text | Google Scholar

Hinckley, C. A., and Ziskind-Conhaim, L. (2006). Electrical coupling between locomotor-related excitatory interneurons in the mammalian spinal cord. J. Neurosci. 26, 8477–8483. doi: 10.1523/JNEUROSCI.0395-06.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Jasinski, P. E., Molkov, Y. I., Shevtsova, N. A., Smith, J. C., and Rybak, I. A. (2013). Sodium and calcium mechanisms of rhythmic bursting in excitatory neural networks of the pre-Bötzinger complex: a computational modelling study. Eur. J. Neurosci. 37, 212–230. doi: 10.1111/ejn.12042

PubMed Abstract | CrossRef Full Text | Google Scholar

Kepler, T. B., Marder, E., and Abbott, L. F. (1990). The effect of electrical coupling on the frequency of model neuronal oscillators. Science 248, 83–85. doi: 10.1126/science.2321028

PubMed Abstract | CrossRef Full Text | Google Scholar

Kiehn, O. (2006). Locomotor circuits in the mammalian spinal cord. Annu. Rev. Neurosci. 29, 279–306. doi: 10.1146/annurev.neuro.29.051605.112910

PubMed Abstract | CrossRef Full Text | Google Scholar

Kiehn, O. (2016). Decoding the organization of spinal circuits that control locomotion. Nat. Rev. Neurosci. 17, 224–238. doi: 10.1038/nrn.2016.9

PubMed Abstract | CrossRef Full Text | Google Scholar

Kjaerulff, O., and Kiehn, O. (1996). Distribution of networks generating and coordinating locomotor activity in the neonatal rat spinal cord in vitro: a lesion study. J. Neurosci. 16, 5777–5794. doi: 10.1523/JNEUROSCI.16-18-05777.1996

PubMed Abstract | CrossRef Full Text | Google Scholar

Koizumi, H., and Smith, J. C. (2008). Persistent Na⁺ and K⁺-dominated leak currents contribute to respiratory rhythm generation in the pre-Bötzinger complex in vitro. J. Neurosci. 28, 1773–1785. doi: 10.1523/JNEUROSCI.3916-07.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Kopell, N., and Ermentrout, B. (2004). Chemical and electrical synapses perform complementary roles in the synchronization of interneuronal networks. Proc. Natl. Acad. Sci. U.S.A. 101, 15482–15487. doi: 10.1073/pnas.0406343101

PubMed Abstract | CrossRef Full Text | Google Scholar

Kueh, D., Barnett, W. H., Cymbalyuk, G. S., and Calabrese, R. L. (2016). Na+/K+ pump interacts with the h-current to control bursting activity in central pattern generator neurons of leeches. eLife 5:e19322. doi: 10.7554/eLife.19322.016

PubMed Abstract | CrossRef Full Text | Google Scholar

McCrea, D. A., and Rybak, I. A. (2007). Modeling the mammalian locomotor CPG: insights from mistakes and perturbations. Prog. Brain Res. 165, 235–253. doi: 10.1016/S0079-6123(06)65015-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Molkov, Y. I., Rubin, J. E., Rybak, I. A., and Smith, J. C. (2016). Computational models of the neural control of breathing. Wiley Interdiscip. Rev. Syst. Biol. Med. 9. doi: 10.1002/wsbm.1371

CrossRef Full Text | Google Scholar

Pfeuty, B., Mato, Golomb, and Hansel, D. (2003). Electrical synapses and synchrony: the role of intrinsic currents. J. Neurosci. 23, 6280–6294. doi: 10.1523/JNEUROSCI.23-15-06280.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

Phillips, R. S., and Rubin, J. E. (2019). Effects of persistent sodium current blockade in respiratory circuits depend on the pharmacological mechanism of action and network dynamics. PLoS Comput. Biol. 15:e1006938. doi: 10.1371/journal.pcbi.1006938

PubMed Abstract | CrossRef Full Text | Google Scholar

Picton, L. D., Zhang, H. Y., and Sillar, K. T. (2017). Sodium pump regulation of locomotor control circuits. J. Neurophysiol. 118, 1070–1081. doi: 10.1152/jn.00066.2017

PubMed Abstract | CrossRef Full Text | Google Scholar

Radosevic, M., Willumsen, A., Petersen, P. C., Lindén, H., Vestergaard, M., and Berg, R. W. (2019). Decoupling of timescales reveals sparse convergent CPG network in the adult spinal cord. Nat. Commun. 10:2937. doi: 10.1038/s41467-019-10822-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Rash, J. E., Dillman, R. K., Bilhartz, B. L., Duffy, H. S., Whalen, L. R., and Yasumura, T. (1996). Mixed synapses discovered and mapped throughout mammalian spinal cord. Proc. Natl. Acad. Sci. U.S.A. 93, 4235–4239. doi: 10.1073/pnas.93.9.4235

PubMed Abstract | CrossRef Full Text | Google Scholar

Rybak, I. A., Dougherty, K. J., and Shevtsova, N. A. (2015). Organization of the mammalian locomotor CPG: review of computational model and circuit architectures based on genetically identified spinal interneurons. eNeuro 2:ENEURO.0069-15.2015. doi: 10.1523/ENEURO.0069-15.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

Rybak, I. A., Shevtsova, N. A., and Kiehn, O. (2013). Modelling genetic reorganization in the mouse spinal cord affecting left-right coordination during locomotion. J. Physiol. 591, 5491–5508. doi: 10.1113/jphysiol.2013.261115

PubMed Abstract | CrossRef Full Text | Google Scholar

Rybak, I. A., Shevtsova, N. A., Lafreniere-Roula, M., and McCrea, D. A. (2006a). Modelling spinal circuitry involved in locomotor pattern generation: insights from deletions during fictive locomotion. J. Physiol. 577, 617–639. doi: 10.1113/jphysiol.2006.118703

PubMed Abstract | CrossRef Full Text | Google Scholar

Rybak, I. A., Shevtsova, N. A., Ptak, K., and McCrimmon, D. R. (2004). Intrinsic bursting activity in the pre-Bötzinger complex: role of persistent sodium and potassium currents. Biol. Cyber. 90, 59–74. doi: 10.1007/s00422-003-0447-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Rybak, I. A., Stecina, K., Shevtsova, N. A., and McCrea, D. A. (2006b). Modelling spinal circuitry involved in locomotor pattern generation: insights from the effects of afferent stimulation. J. Physiol. 577, 641–658. doi: 10.1113/jphysiol.2006.118711

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharples, S. A., Burma, N. E., Borowska-Fielding, J., Kwok, C. H. T., Eaton, S. E. A., Baker, G. B., et al. (2020). A dynamic role for dopamine receptors in the control of mammalian spinal networks. Sci. Rep. 10:16429. doi: 10.1038/s41598-020-73230-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Sherman, A., and Rinzel, J. (1992). Rhythmogenic effects of weak electrotonic coupling in neuronal models. Proc. Natl. Acad. Sci. U.S.A. 89, 2471–2474. doi: 10.1073/pnas.89.6.2471

PubMed Abstract | CrossRef Full Text | Google Scholar

Sherwood, W. E., Harris-Warrick, R., and Guckenheimer, J. (2011). Synaptic patterning of left-right alternation in a computational model of the rodent hindlimb central pattern generator. J. Comp. Neurosci. 30, 323–360. doi: 10.1007/s10827-010-0259-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Shevtsova, N. A., and Rybak, I. A. (2016). Organization of flexor-extensor interactions in the mammalian spinal cord: insights from computational modelling. J. Physiol. 594, 6117–6131. doi: 10.1113/JP272437

PubMed Abstract | CrossRef Full Text | Google Scholar

Shevtsova, N. A., Talpalar, A. E., Markin, S. N., Harris-Warrick, R. M., Kiehn, O., and Rybak, I. A. (2015). Organization of left-right coordination of neuronal activity in the mammalian spinal cord: insights from computational modelling. J. Physiol. 593, 2403–2426. doi: 10.1113/JP270121

PubMed Abstract | CrossRef Full Text | Google Scholar

Shevtsova, N. A., Ha, N. T., Rybak, I. A., and Dougherty, K. J. (2020). Neural interactions in developing rhythmogenic spinal networks: Insights from computational modeling. Biorxiv [Preprint]. doi: 10.1101/2020.09.15.298281

CrossRef Full Text | Google Scholar

Smith, J. C., Abdala, A. P., Borgmann, A., Rybak, I. A., and Paton, J. F. (2013). Brainstem respiratory networks: building blocks and microcircuits. Trends Neurosci. 36, 152–162. doi: 10.1016/j.tins.2012.11.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, J. C., Abdala, A. P. L., Koizumi, H., Rybak, I. A., and Paton, J. F. R. (2007). Spatial and functional architecture of the mammalian brainstem respiratory network: a hierarchy of three oscillatory mechanisms. J. Neurophysiol. 98, 3370–3387. doi: 10.1152/jn.00985.2007

PubMed Abstract | CrossRef Full Text | Google Scholar

Stimberg, M., Brette, R., and Goodman, D. F. M. (2019). Brian 2, an intuitive and efficient neural simulator. eLife 8:e47314. doi: 10.7554/eLife.47314.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Stimberg, M., Goodman, D. F. M., Benichoux, V., and Brette, R. (2014). Equation-oriented specification of neural models for simulations. Front. Neuroinform. 8:6. doi: 10.3389/fninf.2014.00006

PubMed Abstract | CrossRef Full Text | Google Scholar

Tazerart, S., Viemari, J.-C., Darbon, P., Vinay, L., and Brocard, F. (2007). Contribution of persistent sodium current to locomotor pattern generation in neonatal rats. J. Neurophysiol. 98, 613–628. doi: 10.1152/jn.00316.2007

PubMed Abstract | CrossRef Full Text | Google Scholar

Tazerart, S., Vinay, L., and Brocard, F. (2008). The persistent sodium current generates pacemaker activities in the central pattern generator for locomotion and regulates the locomotor rhythm. J. Neurosci. 28, 8577–8589. doi: 10.1523/JNEUROSCI.1437-08.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Tong, H., and McDearmid, J. R. (2012). Pacemaker and plateau potentials shape output of a developing locomotor network. Cur. Biol. 22, 2285–2293. doi: 10.1016/j.cub.2012.10.025

PubMed Abstract | CrossRef Full Text | Google Scholar

Tovar, K. R., Maher, B. J., and Westbrook, G. L. (2009). Direct actions of carbenoxolone on synaptic transmission and neuronal membrane properties. J. Neurophysiol. 102, 974–978. doi: 10.1152/jn.00060.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Tresch, M. C., and Kiehn, O. (2000). Motor coordination without action potentials in the mammalian spinal cord. Nature Neurosci. 3, 593–599. doi: 10.1038/75768

PubMed Abstract | CrossRef Full Text | Google Scholar

Walton, K. D., and Navarrete, R. (1991). Postnatal changes in motoneurone electrotonic coupling studied in the in vitro rat lumbar spinal cord. J. Physiol. 433, 283–305. doi: 10.1113/jphysiol.1991.sp018426

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilson, J. M., Cowan, A. I., and Brownstone, R. M. (2007). Heterogeneous electrotonic coupling and synchronization of rhythmic bursting activity in mouse Hb9 interneurons. J. Neurophysiol. 98, 2370–2381. doi: 10.1152/jn.00338.2007

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilson, J. M., Hartley, R., Maxwell, D. J., Todd, A. J., Lieberam, I., Kaltschmidt, J. A., et al. (2005). Conditional rhythmicity of ventral spinal interneurons defined by expression of the Hb9 homeodomain protein. J. Neurosci, 25, 5710–5719. doi: 10.1523/JNEUROSCI.0274-05.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhong, G., Droho, S., Crone, S. A., Dietz, S., Kwan, A. C., Webb, W. W., et al. (2010). Electrophysiological characterization of V2a interneurons and their locomotor-related activity in the neonatal mouse spinal cord. J. Neurosci. 30, 170–182. doi: 10.1523/JNEUROSCI.4849-09.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhong, G., Masino, M. A., and Harris-Warrick, R. M. (2007). Persistent sodium currents participate in fictive locomotion generation in neonatal mouse spinal cord. J. Neurosci. 26, 6509–6517. doi: 10.1523/JNEUROSCI.1410-06.2006

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhong, G., Shevtsova, N. A., Rybak, I. A., and Harris-Warrick, R. M. (2012). Neuronal activity in the isolated mouse spinal cord during spontaneous deletions in fictive locomotion: insights into locomotor central pattern generator organization. J. Physiol. 590, 4735–4759. doi: 10.1113/jphysiol.2012.240895

PubMed Abstract | CrossRef Full Text | Google Scholar

Ziskind-Conhaim, L., Wu, L., and Wiesner, E. P. (2008). Persistent sodium current contributes to induced voltage oscillations in locomotor-related Hb9 interneurons in the mouse spinal cord. J. Neurophysiol. 100, 2254–2264. doi: 10.1152/jn.90437.2008

PubMed Abstract | CrossRef Full Text | Google Scholar