On the Organization of the Locomotor CPG: Insights From Split-Belt Locomotion and Mathematical Modeling

Rhythmic limb movements during locomotion are controlled by central pattern generator (CPG) circuits located in the spinal cord. It is considered that these circuits are composed of individual rhythm generators (RGs) for each limb interacting with each other through multiple commissural and long propriospinal circuits. The organization and operation of each RG are not fully understood, and different competing theories exist about interactions between its flexor and extensor components, as well as about left–right commissural interactions between the RGs. The central idea of circuit organization proposed in this study is that with an increase of excitatory input to each RG (or an increase in locomotor speed) the rhythmogenic mechanism of the RGs changes from “flexor-driven” rhythmicity to a “classical half-center” mechanism. We test this hypothesis using our experimental data on changes in duration of stance and swing phases in the intact and spinal cats walking on the ground or tied-belt treadmills (symmetric conditions) or split-belt treadmills with different left and right belt speeds (asymmetric conditions). We compare these experimental data with the results of mathematical modeling, in which simulated CPG circuits operate in similar symmetric and asymmetric conditions with matching or differing control drives to the left and right RGs. The obtained results support the proposed concept of state-dependent changes in RG operation and specific commissural interactions between the RGs. The performed simulations and mathematical analysis of model operation under different conditions provide new insights into CPG network organization and limb coordination during locomotion.


INTRODUCTION
It is commonly accepted that the spinal locomotor central pattern generator (CPG) circuits include separate rhythm generators (RGs) that each control a single limb and interact with each other via multiple commissural and long propriospinal pathways. These connections set up phase relationships between the RGs and thus coordinate limb movements and locomotor gait Danner et al., 2017). Each RG is thought to contain two excitatory neuron populations representing flexor and extensor half-centers connected by reciprocal inhibition, whose activity defines the flexor and extensor phases of limb movements, respectively. According to the classical half-center concept (Brown, 1914), switching between the flexor and extensor activity phases (for review see McCrea and Rybak, 2008;Stuart and Hultborn, 2008) occurs through a so-called release mechanism (Wang and Rinzel, 1992) based on an adapting (decrementing) activity of each half-center and mutual inhibition between them. This mechanism does not necessarily require the ability of each half-center to intrinsically generate rhythmic activity, and the resultant RG pattern is usually flexorextensor balanced, so that the durations of both phases are approximately equal.
The other potential mechanism is based on the intrinsic ability of one or both half-centers to generate rhythmic bursting (Wang and Rinzel, 1992;Skinner et al., 1994;Marder and Calabrese, 1996;Marder and Bucher, 2001). Optogenetic studies in the isolated spinal cord have demonstrated that rhythmic flexor and extensor activities can be evoked in certain conditions independent of each other (Hägglund et al., 2013), confirming that both flexor and extensor half-centers are conditional intrinsic oscillators, i.e., capable of endogenous generation of rhythmic bursting activity. Pearson and Duysens (1976) and Duysens (2006) have previously proposed a flexor-driven concept (so called swing generator model (for review see Duysens et al., 2013), in which only the flexor half-center is intrinsically rhythmic, hence representing a true RG, while the extensor half-center shows sustained activity if uncoupled and exhibits anti-phase oscillations due to rhythmic inhibition from the flexor half-center.
To meet both concepts, we previously suggested that both half-centers are conditional oscillators, whose ability to intrinsically generate rhythmic bursting depends on the level of excitation Danner et al., 2016Danner et al., , 2017Shevtsova and Rybak, 2016). In this case, a relatively strong excitation of the extensor half-center keeps it in the mode of sustained activity (if uncoupled), whereas a relatively weak excitatory drive to the flexor half-center allows generation of intrinsic oscillations. Therefore, the mechanism for rhythm generation in the RG may vary and, depending on external drives to its half-centers or their level of excitation, it can operate according to the classical half-center or the flexordriven scenario as was previously demonstrated and analyzed by Ausborn et al. (2018).
In the present study, we extend the RG model of Ausborn et al. (2018) by assuming that increased activation of the flexor halfcenter is accompanied by a corresponding decrease in the activity of the extensor half-center. To implement this assumption in the model, we suggested that the external excitatory drive to the flexor half-center simultaneously provides inhibition to the extensor half-center, thus reducing the level of its excitation and directing its operation toward intrinsic rhythmicity. To this end, with an increase of the drive to the RG (with the corresponding increase of oscillation frequency) the operating rhythmogenic mechanism changes from the flexor-driven rhythmicity to classical half-center oscillations with a quasi-balanced flexorextensor pattern.
To study the behaviors of the proposed RGs in the context of left-right interactions and limb coordination, we incorporated these new RG implementations in the model of left-right circuit interactions in the spinal cord previously described by Danner et al. (2019). The resultant model included two (left and right) RGs interacting via several commissural pathways presumably mediated by genetically identified V0 V , V0 D , and V3 interneurons. The main goal of this study was to investigate leftright interactions and coordination under different symmetric and asymmetric conditions, which were defined by the same or different drives to left and right RGs, respectively. We assumed these conditions to be, at first approximation, similar to overground or regular tied-belt treadmill locomotion in cats (symmetric conditions) and their stepping on split-belt treadmills with different speeds of the left and right belts (asymmetric conditions). The experimental data were collected from intact and spinal cats in previously published Kuczynski et al., 2017) and new experiments. These experimental data were compared with the results of our simulations, in which the modeled circuits operated in similar symmetric and asymmetric conditions. We used these comparisons to evaluate the plausibility of our model and, thus, to formulate important insights into the organization of spinal CPG circuits and their role in limb coordination during locomotion.

Ethics Statement
All procedures were approved by the Animal Care Committee of the Université de Sherbrooke in accordance with policies and directives of the Canadian Council on Animal Care . The current dataset of treadmill locomotion was obtained from eight intact (four females and four males) and six spinal (four females and two males) adult cats weighing between 3.5 and 5.0 kg. These data were taken from previous studies and were either reanalyzed (from Frigon et al., 2015) or reused (from Frigon et al., 2017) (see Table 1). One cat was studied in both the intact and spinal states (IN-3 and SP-4 in Table 1). Only one new cat was used to get new data during overground locomotion . Before and after experiments, cats were housed and fed in a dedicated room within the animal care facility of the Faculty of Medicine and Health Sciences at the Université de Sherbrooke. As part of our effort to maximize the scientific output of each animal, 10 of 11 animals were used in other studies to answer different scientific questions Thibaudier et al., 2013Thibaudier et al., , 2017D'Angelo et al., 2014;Thibaudier and Frigon, 2014;Dambreville et al., 2015Dambreville et al., , 2016Hurteau et al., 2015Hurteau et al., , 2017Kuczynski et al., 2017;Harnie et al., 2018;Desrochers et al., 2019). The experimental studies complied with the ARRIVE guidelines (Kilkenny et al., 2010) and principles of animal research established by the Journal of Physiology (Grundy, 2015).  (14) SP-2 (female) Tied-belt (6-15); split-belt (12-15) SP-3 (female) Tied-belt (14); split-belt (9-14) SP-4 (male) Tied-belt (12-15); split-belt (6-13) SP-5 (male) ; split-belt (14) SP-6 (female) Tied-belt (12-14); split-belt (12-15)

Surgical Procedures
Surgical procedures were described in detail in Frigon et al. (2015Frigon et al. ( , 2017 and also apply to the new cat used here. Briefly, we performed all surgical procedures in an operating room with sterilized equipment. Before surgery, the cat was sedated with an intramuscular (i.m.) injection of butorphanol (0.4 mg/kg), acepromazine (0.1 mg/kg), and glycopyrrolate (0.01 mg/kg). Induction was done with Ketamine/Diazepam (0.11 ml/kg in a 1:1 ratio, i.m.). The fur overlying the back, stomach, and hindlimbs was shaved. The cat was then anesthetized with isoflurane (1.5-3%) using a mask for a minimum of 5 min and then intubated with a flexible endotracheal tube. We confirmed isoflurane concentration during surgery by monitoring cardiac and respiratory rates, by applying pressure to the paw to detect limb withdrawal, and by assessing muscle tone. A rectal thermometer was used to monitor body temperature and keep it between 35 • and 37 • C using a water-filled heating pad placed under the animal and an infrared lamp positioned ∼50 cm above the cat. During each surgery, we injected an antibiotic (Convenia, 0.1 ml/kg) subcutaneously and a transdermal fentanyl patch (25 mcg/hr) was taped to the back of the animal 2-3 cm rostral to the base of the tail. During surgery and approximately 7 h later, another analgesic (Buprenorphine 0.01 mg/kg) was administered subcutaneously. After surgery, cats were placed in an incubator and closely monitored until they regained consciousness. At the conclusion of the experiments, cats received a lethal dose of pentobarbital through the left or right cephalic vein.

Spinal transection
The spinal cord was completely transected at low thoracic levels in six cats (four females and two males; see Frigon et al., 2017). A small laminectomy was performed between the junction of the 12th and 13th vertebrae. After exposing the spinal cord, lidocaine (Xylocaine, 2%) was applied topically and injected within the spinal cord. The spinal cord was then transected with surgical scissors. Hemostatic material (Spongostan) was then inserted within the gap and muscles and skin were sewn back to close the opening in anatomic layers. Following spinalization and for the remainder of the study, the bladder was manually expressed 1-2 times each day. The hindlimbs were frequently cleaned by placing the lower half of the body in a warm soapy bath. For training the recovery of hindlimb locomotion, see Frigon et al. (2017). Cats studied in the intact (IN-1 to  and spinal (SP1 to SP-6) states. Tied-belt locomotion was studied from 0.4 m/s to 1.0 m/s in the intact state and from 0.1 m/s to 1.0 m/s in the spinal state. Split-belt locomotion was studied with the slow hindlimb stepping at 0.4 m/s and the fast hindlimb from 0.5 m/s to 1.0 m/s in both states.

Implantation
All 13 cats were implanted with electrodes to chronically record muscle activity (EMG, electromyography), although EMG recordings in the present studies were only used for demonstration that flexor and extensor burst durations change in parallel with phase durations (i.e., swing and stance). Pairs of Teflon insulated multistrain fine wires (AS633; Cooner wire, Chatsworth, CA, United States) were directed subcutaneously from 1-2 head-mounted 34-pin connectors (Omnetics Connector Corporation, Minneapolis, MN, United States) and sewn into the belly of selected hindlimb muscles for bipolar recordings. We verified electrode placement by electrically stimulating each muscle through the appropriate head connector channel.

Experimental paradigms
Experiments in the 12 cats from previous studies  were performed on an animal treadmill with two independently controlled running surfaces 120 cm long and 30 cm wide (Bertec Corporation, Columbus, OH, United States). Cats performed three locomotor paradigms: (1) Tied-belt locomotion from 0.1 m/s (spinal cats) or 0.4 m/s (intact cats) up to 1.0 m/s in 0.1 m/s increments; (2) split-belt locomotion with one side (slow side) stepping at 0.4 m/s and the other side (fast side) stepping from 0.5 m/s to 1.0 m/s in 0.1 m/s increments; (3) split-belt locomotion with the slow side stepping at 0.1 m/s and the fast side stepping from 0.2 m/s to 1.0 m/s in 0.1 m/s increments (spinal cats only). In spinal cats, the forelimbs remained on a stationary platform with a Plexiglas separator placed between hindlimbs. In the cat that contributed new data, we trained the animal to step along an oval-shaped walkway at self-selected speeds. The walkway has 2.07 m straight paths (0.32 m wide) on each side and we only analyzed data during straight path stepping.

Data acquisition and analysis
Videos of the left and right sides during overground and treadmill locomotion were captured with two cameras (Basler AcA640-100 gm) at 60 frames per second with a spatial resolution of 640 by 480 pixels. A custom-made Labview program acquired images and synchronized the cameras with the EMG. Videos were analyzed off-line at 60 frames per second using custom-made software. Contact of the paw and its most caudal displacement were determined for both hindlimbs by visual inspection. We defined paw contact as the first frame where the paw made visible contact with the treadmill surface while the most caudal displacement of the limb was the frame with the most caudal displacement of the toe. We measured cycle duration from successive contacts of the same hindpaw while stance duration corresponded to the interval of time from paw contact to the most caudal displacement of the limb. Swing duration was measured as cycle duration minus stance duration. Durations from 6 to 15 cycles for each limb were averaged for an episode during treadmill locomotion. In one cat, we obtained and analyzed 44 cycles from 10 runs of overground locomotion.
The EMG was pre-amplified (×10, custom-made system), bandpass filtered (30-1000 Hz) and amplified (×100-5000) using a 16-channel amplifier (AM Systems Model 3500, Sequim, WA, United States). EMG data were digitized (2000 Hz) with a National Instruments card (NI 6032E) and acquired with custom-made acquisition software and stored on computer. The EMG data set shown came from recordings in the anterior sartorius (Srt, hip flexor/knee extensor), the vastus lateralis (VL, knee extensor) and the lateral gastrocnemius (LG, ankle plantarflexor/knee flexor).

Mathematical Modeling
We implemented a reduced mathematical model based on the work of Danner et al. (2017). Simulating flexor and extensor half-centers using activity-based neuron models describing neuron populations (Ermentrout, 1994) significantly simplifies mathematical analysis. The voltage variable of each flexor and extensor units represents the average voltage of the population of flexor and extensor neurons. Such a reduction provides an accurate description of the network dynamics in the CPG circuits controlling mammalian locomotion (Molkov et al., 2015;Ausborn et al., 2018). The CPG network controlling rhythmic locomotion is known to include both excitatory and inhibitory connections between flexor half-centers Molkov et al., 2015;Shevtsova et al., 2015;Danner et al., 2016Danner et al., , 2017Danner et al., , 2019Shevtsova and Rybak, 2016;Ausborn et al., 2019). We only included reciprocal inhibition between flexors in the model assuming a net inhibitory interaction. Flexor and extensor halfcenters comprising left and right RGs also inhibit each other. Additionally, the model included inhibition from extensors to contralateral flexors. This connection was first introduced by Danner et al. (2017) who found that inhibition of flexor halfcenters by contralateral extensor stabilize anti-phase left-right alternations in corresponding gaits. In this study, we show that this interaction is essential for symmetric left-right alternations and explain the mechanism.
All neurons were modeled using the formalism described in Rubin et al. (2009) and then used in a number of previous publications (Rubin et al., 2011;Molkov et al., 2014Molkov et al., , 2015Molkov et al., , 2016Danner et al., 2016Danner et al., , 2017Danner et al., , 2019Ausborn et al., 2018Ausborn et al., , 2019. Intrinsic bursting properties resulted from slowly inactivating sodium current dynamics. The membrane potential (V) of flexors and extensors was governed by the following equation: Here, C is the capacitance, t is time, I L is the leak current, I NaP is the slowly inactivating (persistent) sodium current, and I syn is the synaptic current that is the sum of input currents from other neurons and the excitatory drive current. The leak current and the persistent sodium current were defined in the same manner in flexors and in extensors. (2) In the expression for the leak current (2), g L is the conductance of the leak current and E L is the leak reversal potential. In the expression for the persistent sodium current (3), g NaP is the persistent sodium maximal conductance and E Na is the sodium reversal potential. m NaP∞ (V) is the voltage-dependent steady-state activation function of the persistent sodium current. Persistent sodium current activation is considered to be instantaneous. h NaP is the persistent sodium inactivation gating variable. The steady state activation functions for persistent sodium activation and inactivation are given by the following expressions: and the dynamics of the persistent sodium inactivation variable were governed by the following differential equation: Here, τ NaP (V) is the voltage-dependent time constant for the inactivation of the persistent sodium current. In the gating variable expressions, V xNaP is the half-(in)activation voltage and k xNaP is the (in)activation slope, where x ∈ {m, h, τ } . In the differential equation for the membrane potential the third current is the synaptic current I syn and is defined by the synaptic input from neurons in the network as well as external drives. For flexors, this included input from the contralateral flexor, the ipsilateral extensor, and the contralateral extensor. For extensors, the synaptic current included input from the ipsilateral flexor. In flexors and extensors, drive was implemented as the conductance of an excitatory input. The general expression for the synaptic current in neuron i is as follows: Here, d i is the excitatory drive to neuron i and V i is the voltage of neuron i. E ex is the reversal potential for the excitatory synaptic currents. I syni includes the sum over all synaptic inputs from j = 1 : 4 [Left Flexor (1), Right Flexor (2), Left Extensor (3), Right Extensor (4), see Figure 2]. d 1 and d 2 are drives to flexors. d 3 and d 4 are net drives to extensors representing the difference between constant drive to an extensor and drive to the ipsilateral flexor (Drive to E -Drive to F in Figure 1B) which implements the inhibitory effect of flexor drives on extensors. Drive values are varied as explained in the corresponding subsections of "Results." E inh is the reversal potential for the inhibitory synaptic currents. b ji is the weight of the synaptic connection from neuron j to neuron i, which represents the maximal conductance of the corresponding synaptic channel. f (V) is the activity (normalized firing rate) as a function of voltage and is defined by the following piecewise linear function.
The activity function f (V) varies from 0 to 1. Here, V min and V max define the voltages at which threshold and saturation are reached, respectively. The values of all parameters are provided in Table 2. In our simulations, the synaptic weights of commissural connections b 12 , b 21 , b 41 and b 32 were varied, while synaptic weights within each RG b 31 , b 42 , b 13 and b 24 were fixed.

Model of the Rhythm Generator (RG) Controlling a Single Limb
In the present study, we accepted the model of Ausborn et al. (2018) and their suggestion that rhythmic activity in the RG may be based on flexor-driven or classical half-center mechanisms, depending on the level of excitation of flexor and extensor half-centers, both considered conditional bursters. They independently varied flexor and extensor drives and identified parameter areas in which the above mechanisms operate. Here, we extended the model of Ausborn et al. by using the assumption that an increase in activation of the flexor half-center is accompanied by a decrease in the activity of the extensor halfcenter. Specifically, we assumed that the excitatory drive to the flexor half-center provides inhibition to the extensor half-center (through inhibitory interneurons), reducing the initial level of its excitation (Figures 1A,B). In this case, at relatively low drives to the flexor half-center, the frequency of RG oscillations (defined by flexor activity) is low, and the locomotor pattern is not balanced, i.e., has a short flexor and long extensor bursts. An increase in the drive to the flexor half-center increases the RG frequency, making the pattern more flexor-extensor balanced while concurrently reducing the level of excitation of the extensor half-center, shifting the extensor half-center's operation toward an intrinsically rhythmic state. Figure 1C shows a two-parameter frequency dependence on flexor and extensor drives similar to shown in Ausborn et al. (2018) that was calculated for a set of parameters used in the present study. According to our suggestion, with the changes in the drive to flexor half-center (Drive to F) and the net drive to extensor half-center (Drive to E minus Drive to F), the parameter point representing a state of RG operation moves along the yellow line intersecting both areas for flexor-driven and classical half-center oscillations ( Figure 1C). Specifically, with an increase of drive to flexor center, the RG Maximal conductance ( nS ) g L = 2.8, g NaP = 5 Reversal potentials (mV) Threshold and saturation voltage (mV) operation regimes shifts from flexor-driven intrinsic oscillations (with short flexor bursts and long extensor bursts, Figure 1D1) toward the classical half-center mechanism of rhythmicity with a quasi-balanced flexor-extensor pattern ( Figure 1D2).

Commissural Interactions Between RGs Controlling Left and Right Limbs
The main goal of this study was to investigate left-right coordination of limb movements under different symmetric and asymmetric conditions. Left-right limb coordination relies on neural interactions between the two RGs controlling the left and right limbs. The connectome of these interactions was drawn from the model of Danner et al. (2019). In that model, the left and right RGs interacted via three commissural pathways (Figure 2A). Two of them, mediated by genetically identified inhibitory V0 D and excitatory V0v (V2a-V0v paths, acting via the inhibitory Ini populations) populations of commissural interneurons (CINs), promoted left-right alternation (Talpalar et al., 2013) through mutual inhibition between the left and right flexor half-centers (see also Shevtsova et al., 2015). The third pathway, mediated by genetically identified V3 CINs, promoted left-right synchronization via mutual excitation between the left and right extensor half-centers and diagonal inhibition of the contralateral flexor half-centers (Danner et al., 2016; see Figure 2A. In the present study, to simplify the model and make it more mathematically tractable, all commissural interactions were replaced by functionally equivalent direct connections, as shown in Figure 2B.

Speed-Dependent Changes in Phase Durations During Left-Right Symmetric and Asymmetric Locomotion
Our objective was to evaluate the RG circuit organization proposed above by considering their operation in two cases: a symmetric case, when left and right drives vary but remain equal, and an asymmetric case, when one of two drives changes while the other maintains a constant value. We assumed that these two regimes are functionally comparable to regular overground or tied-belt treadmill locomotion (symmetric case) and split-belt treadmill locomotion with different speeds for the left and right belts (asymmetric case). We focused on the analysis of speeddependent changes in the durations of the main locomotor phases (swing and stance) using data from previous experiments during tied-belt and split-belt treadmill locomotion in intact and spinal cats  and new experiments performed during overground locomotion in an intact cat.

Speed-Dependent Changes in Phase Durations During Left-Right Symmetric Locomotion
Left-Right Symmetric Locomotion in Cats reduction of stance phase duration with small or absent changes in swing phase duration, consistent with previous studies in cats (Halbertsma, 1983;Frigon and Gossard, 2009;Frigon et al., 2013Frigon et al., , 2014Frigon et al., , 2017). An interesting difference between the three cases shown in Figure 3 is that during overground locomotion in intact cats, at a speed of ∼1.1 m/s, the swing and stance phase durations become equal and then at higher speeds, stance becomes shorter than swing (Figure 3A). Despite a similar tendency, stance did not become shorter than swing during tied-belt treadmill locomotion in intact ( Figure 3C) or spinal (Figure 3D) cats. The treadmill locomotion is not usually performed at speeds greater than 1.0 m/s in intact cats, because of safety concerns, as well as in spinal cats, in which the pattern starts to break down. Nevertheless, spinal cats reached swing-stance equality at about 1.0 m/s ( Figure 3D).

Simulation of Left-Right Symmetric Regime With the Model
The schematic of our simplified model is shown in Figure 2B. In this model there are mutual inhibitory interactions between the flexor half-centers, which combine and simplify two inhibitory pathways mediated by V0 D and V0 V CINs in Figure 2A. This inhibition is referred to as "flexor-flexor" (or F-F) inhibition. In addition, there are also inhibitory pathways from each extensor half-center to the contralateral flexor half-center (Figure 2A), which are presumably mediated by V3 CINs through inhibitory populations, such as V1 . The strength of this connection in the present model is referred to as "extensorflexor" (or E-F) inhibition. We therefore have four control parameters in the model: the drives to both flexor half-centers (which also define the inhibitory inputs to the extensor halfcenters; these drives are equal in the symmetrical case) and F-F and E-F inhibitions. First, we simulated the changes in locomotor phase durations in response to increasing drive to a single RG (Figure 4). The external drive to the RGs was increased from 0.2 to 0.8 producing progressively shorter extension at relatively constant flexion duration. With an increase of external drive, the frequency of oscillations increased from about 0.4 to about 1.4 Hz. The increase in frequency (decrease in the period of oscillations) occurred mainly by shortening the extensor phase with minor changes in the duration of the flexor phase. The predominant decrease in extensor phase qualitatively corresponds to the change in the duration of stance and swing phases observed with increasing locomotor speed in experimental studies ( Figure 3A). Note that in our simulations, the flexor and extensor phases become equal at drive values of about 0.7, after which extension becomes shorter than flexion (similar to that in Figure 3A). This reversal in flexor-extensor durations occurs in our model because flexor and extensor half-centers receive the same external excitation at a drive value of approximately 0.7 (see Figure 1C).
To explore the system's behavior in terms of leftright coordination, we simulated the model and identified synchronization patterns while varying inhibition strengths at different drive values. Figures 5A-D shows the parameter plane partitions for four representative drive values corresponding to low and high frequencies. Qualitatively, the F-F inhibition promotes alternating (anti-phase) flexor activity while the E-F inhibition contributes to synchronizing (in-phase) the flexor half-centers due to a phasic reduction in inhibition of flexors during contralateral flexion. Therefore, it is reasonable to expect that at high F-F inhibition and low E-F inhibition (an upper-left corner on Figure 5 diagrams), the left and right RGs exhibit alternating activity, and at low F-F inhibition and high E-F inhibition their activity synchronizes at all frequencies.
These regimes of exact anti-phase and in-phase oscillations are observed in the white and black parameter regions, respectively. There is an overlap between the two regions (shown in gray), which corresponds to bistability in the system, where both regimes can operate depending on the initial conditions chosen. A transition from in-phase to anti-phase oscillations occurs at the boundary between the gray and white regions, which is invariant to the drive (Figures 5A-D). An opposite transition occurs at the gray-black boundary, which moves up in terms of F-F inhibition with the drive, thus reducing the bistability area.
There are also regimes of asymmetric alternations at relatively low ( Figure 5A, orange region) and high (Figures 5C,D yellow region) drive values corresponding to low or high locomotor frequencies. At low frequencies (i.e., low drive values), this regime is observed at low values of E-F inhibition; it results from postinhibitory rebound activation of the flexor oscillator after the contralateral flexor deactivates. Slightly higher E-F inhibition strength prevents this post-inhibitory rebound by suppressing the contralateral flexor half-centers for the duration of strong extensor activity in the beginning of the extensor burst. At high locomotor frequencies, the asymmetric alternation regime is practically indistinguishable from pure anti-phase oscillations because the duty cycle is very close to 1/2.
Based on the analysis above, we found that the considered circuit produces robust anti-phase alternations of flexor activity in a certain parameter region for all locomotor frequencies. We chose the exemplary point (0.2, 0.4) that belongs to this region for subsequent simulations. However, this particular choice did not make a qualitative difference in the system's behavior as long as the parameter point chosen belonged to the region of monostable anti-phase oscillations.

Left-Right Asymmetric Locomotion in Cats Stepping on Split-Belt Treadmills
The split-belt treadmill locomotion experiments, in which animals step on belts with different speeds for the left and right sides, is a common way to study limb coordination during locomotion in cats and humans. Many previous studies in cats demonstrated that both intact and spinal animals adapt to such stepping conditions and demonstrate stable locomotion (Kulagin and Shik, 1970;Forssberg et al., 1980;Frigon et al., 2013Frigon et al., , 2015Frigon et al., , 2017D'Angelo et al., 2014;Kuczynski et al., 2017). In these studies, we can separate cat locomotion on the split-belt treadmill in two qualitatively different types of conditions: simple and extreme Kuczynski et al., 2017). In the simple condition, characterized by a relatively small speed difference between moving belts, animals maintain a 1:1 ratio between the number of steps made by left and right limbs. In extreme conditions, the animal starts taking more steps on the fast side compared to the slow side resulting in step ratios of 1:2, 1:3, 1:4, etc. (Forssberg et al., 1980;Frigon et al., 2015Frigon et al., , 2017. The changes in locomotor phase durations during split-belt locomotion of intact and spinal cats in simple conditions are shown in Figures 6A-C, respectively (see also Frigon et al., 2015Frigon et al., , 2017. In both cases, the slow hindlimb (SHL) stepped at a constant speed of 0.4 m/s, whereas the speed of the fast hindlimb (FHL) belt increased from 0.5 to 1.0 m/s. In these conditions, the important characteristics of locomotion observed are (see also Frigon et al., 2015Frigon et al., , 2017: (1) The step cycle period remains equal in both hindlimbs (FHL and SHL). (2) In the SHL, the durations of swing and stance phases do not change much. (3) In the FHL, the duration of stance decreases, whereas the duration of swing increases, allowing step cycle duration to remain relatively unchanged despite an increase in speed of the FHL. At FHL speed of 0.9 m/s in intact cats and 0.8 m/s in spinal cats, the durations of swing and stance phases become approximately equal and then the flexion duration and swing phase in spinal cats becomes longer than the stance phase at faster FHL speeds (Figures 6B,C, right).
The locomotor characteristics of both intact and spinal cats differ in extreme conditions, when the speed ratio between the slow and fast belts are set to 1:3 and more, up to 1:10 Kuczynski et al., 2017). In this case, the locomotor pattern changes in such a way that cats take more steps on the fast side than on the slow side. Specifically, at 1:3 and 1:4 speed ratios, the limbs on the fast side perform 2-3 steps for every step of the limb on the slow side (1:2 and 1:3 coordination pattern), whereas at ratios of 1:5 or higher, 1:4 and 1:5 coordination pattern were observed Kuczynski et al., 2017). Despite inter-animal variability, both intact (Kuczynski et al., 2017) and spinal  cats exhibit 1:2 + coordination patterns.
It is also important to note that the hindlimbs of cats do not show adaptation to prolonged split-belt locomotion, such as the return of symmetry in some interlimb parameters (e.g., step length) (Kuczynski et al., 2017), in contrast to humans (Reisman et al., 2005). In other words, the adjustments in cycle and phase durations observed during split-belt locomotion in intact and spinal cats remain unchanged over time.

Modeling Asymmetric CPG Operation
To simulate asymmetric conditions corresponding to different speeds of the treadmill belts, we varied drives to the left and right RGs in our model independently (Figure 2B), so that if disconnected they would produce unsynchronized flexor/extensor alternations with different frequencies. Due to commissural interactions, the model generated different synchronization patterns depending on parameters. We assumed that the left RG receives a smaller drive. This corresponds to a triangular region above the bisector in the bifurcation diagram shown in Figure 7A. The bisector of the bifurcation diagram corresponds to equal drives, where exact anti-phase left-right alternations of flexor activity are produced at the commissural connection weights chosen.
As we start changing the drives to the fast RG, both RGs remain synchronized (1:1 region in Figure 7A), however, left and right oscillations become asymmetric. Flexor bursts in the fast  Table 1, one cat was studied in both states) and averaged cycle and phase durations for each cat. Each data point is the mean ± standard deviation for the group of intact and spinal cats. RG occur at progressively shorter intervals after flexor bursts. When the drive to the fast RG becomes significantly larger than the drive to the slow RG, the flexor bursts of the fast RG start occurring immediately when flexor bursts of the slow RG end (Figure 7C). In addition, the duration of the flexor bursts of the fast RG becomes progressively longer (see below in relation to Figures 8A,B). These behaviors correspond to the simple asymmetric conditions, described above, where a 1:1 coordination pattern is maintained.
When the frequency of the slow RG is relatively low because of a low drive to the slow RG (left part of the bifurcation diagram in Figure 7A), a transition to extreme conditions (1:2 + coordination patterns) occurs as we increase the drive to the fast RG further (see above). In the 1:2 regime, one flexor burst of the slow RG corresponds to two flexor bursts of the fast RG (1:2 area in Figure 7A). In this regime, the first flexor burst of the fast RG starts immediately after the flexor burst of the slow RG ends (Figure 7D). Further increases in the drive to the left (fast) RG leads to the emergence of 1:3 + patterns (Figure 7E), similar to that observed in extreme conditions in intact and spinal cats (see above). Between 1:1 and 1:2 regions, there is an area of intermittent regimes where either one or two flexor bursts can be produced by the fast RG during the extension phase of the slow RG, which is commonly observed experimentally .

Changes in Locomotor Phase Duration in a Simple Asymmetric Regime (1:1)
Modeling and analysis of locomotor characteristic changes in the simple condition is more functionally relevant than the extreme cases because it occurs frequently during everyday locomotion, such as stepping along a circular path or when turning. Also, these changes provide an indirect test for the CPG network organization predicted by the model. Figures 8A,B show our simulation of such a simple asymmetric case, when the drive to the slow RG was kept constant at 0.5, while the drive to the fast RG increased from 0.5 to 0.8 (see the corresponding arrow in Figure 7A). Similar to the experimental studies during split-belt locomotion in a simple asymmetric case shown in Figure 6, despite the left-right asymmetry, the oscillation period remained almost constant and was largely defined by the slow side. Similarly, the durations of flexor and extensor phases were relatively constant on the slow side but changed dramatically on the fast side with increased drive (Figures 8A,B). The most important feature of the simulated behavior (which corresponded to experimental data in Figure 6) was the increased duration of flexion in the fast RG occurring with increased drive to that RG. We can qualitatively explain this phenomenon in the model as follows. On the slow side, the flexor half-center of the slow RG operates in a rhythmic mode, while its extensor half-center operates in a regime of tonic activity (if disconnected) as it receives higher excitatory drive. Therefore, the generation of flexor bursts in the slow RG occurs endogenously after a well-defined recovery period, which is almost unaffected by the synaptic inputs it receives from the other side (fast RG). On the fast side, however, once the net drive to the extensor half center is low enough (recall that based on our assumption an increase in drive to the flexor half-center is accompanied by a decrease in drive to the extensor half-center; see above), the extensor half-center goes into an intrinsically rhythmic mode, meaning that the duration of extension and its inter-burst intervals start to depend on intrinsic burst recovery mechanisms. At the same time, the flexor halfcenter of the fast RG receives increasingly more excitation, so flexor burst termination becomes more dependent on the onset of extensor half-center inhibition rather than on the flexor's endogenous deactivation. With a progressive reduction of net drive to the extensor half-center, the recovery period for extensor activity gets longer, which extends flexion duration. Therefore, the phenomenon of increasing duration of flexion in the fast RG results from changing the rhythmogenesis mechanism in the fast RG from an intrinsic generation of flexor oscillations to the classical half-center mechanism that was implemented in our RG model.
To illustrate this further, we removed inhibitory external inputs to both (left and right) extensor half-centers (that provided the above transition in the rhythmogenic properties of the extensor half-centers) and replaced them with a constant excitatory drive of 0.7 (see Figure 1C). In this case, rhythmogenesis was always based on intrinsic bursting of flexor half-centers without switching to the classical half-center mechanism. The results of these simulations are shown in Figures 8C,D. Note that (a) the duration of the flexor phase on the fast side never increases, and (b) the step-cycle duration on both sides clearly decreases with increasing drive to the fast RG, both contradicting to experimental observations (see Figure 6).

Organization and Operation of Spinal Rhythm Generators (RGs) Controlling Limb Movements During Locomotion
There are currently two major competing concepts concerning the organization and operation of spinal neuronal RGs. In the classical half-center concept (Brown, 1914), flexor and extensor half-centers do not require intrinsic rhythmic properties (for review see McCrea and Rybak, 2008;Stuart and Hultborn, 2008). Both half-centers operate in qualitatively similar conditions with phase switching defined by a release mechanism (Wang and Rinzel, 1992) that is based on adapting (decrementing) activity of each half-center and mutual inhibition between them. In the classical half-center, the durations of flexor and extensor phases are balanced (or equal). These durations and the corresponding duty cycles can be easily changed by the level of half-center activation or by external drive. At the same time, the control of RG oscillation frequency in this case is problematic as the oscillation period is not very sensitive to the external drive in half-center oscillators (Daun et al., 2009).
In contrast, with the flexor-driven concept (Pearson and Duysens, 1976;Duysens, 2006), the RG rhythm and pattern is defined by the intrinsically rhythmic flexor half-center, while the extensor half-center has sustained activity if uncoupled and only exhibits rhythmic bursting through rhythmic inhibition from the flexor half-center (for review see Duysens et al., 2013). Thus, the frequency of intrinsically generated flexor bursting explicitly depends on flexor half-center excitation. The distinctive feature of this regime is that the flexor bust duration does not change much and most previously suggested intrinsic oscillatory mechanisms, such as those based on intracellular dynamics of ionic concentrations or slow inactivation of ionic channels (Jasinski et al., 2013;Molkov et al., 2015), produce duty cycles of bursting usually less than 0.5 and are likely to operate at low frequencies with short flexor phases and long extensor bursts.
Both concepts have support in certain conditions. Ausborn et al. (2018) demonstrated that both mechanisms can operate depending on the state of half-centers defined by their level of excitation. Here, we used and refined this idea, by suggesting that (a) at low frequencies the extensor half-center is highly excited and operates in a regime of tonic activity, and (b) an increase in excitation of the flexor half-center, which initially operates in the intrinsic bursting regime, is accompanied by a decrease of excitation of the extensor half-center. Mechanistically, such a decrease of the extensor half-center activation may result from a reduction of excitatory afferent inputs to the extensor halfcenter when unloading the limb at the stance-to-swing transition (Pearson, 1995;Dietz and Duysens, 2000). With concurrent increases in flexor and extensor drives, the RG transitions from a flexor-driven mechanism (when the frequency changes mostly with extension duration while flexion duration remains relatively unchanged) to the classical half-center mechanism (when stepping is controlled by changes in the duty cycle at a relatively constant frequency).
The proposed idea combines both the above concepts on the operations of locomotor CPG circuits. In contrast to the previous CPG models based exclusively on the fundamental classical halfcenter concept (Brown, 1914;Wang and Rinzel, 1992;Rybak et al., 2006a,b;McCrea and Rybak, 2008;Daun et al., 2009) and the previous models based on the flexor-driven concept (Molkov et al., 2015;Shevtsova et al., 2015;Danner et al., 2016Danner et al., , 2017Danner et al., , 2019, we suggest that the operation of locomotor CPG circuits is state-dependent and is particularly dependent on the locomotor speed. To test this idea, we incorporated the above RGs in a model of spinal CPG circuits with reciprocal commissural interactions and used this bilateral RG model to simulate speed-dependent changes in the locomotor pattern of intact and spinal cats in symmetrical (during overground and tiedbelt locomotion) and asymmetrical (during split-belt treadmill locomotion) conditions. The experimental data from previously published Kuczynski et al., 2017) and new experiments were analyzed. The model reproduced and explained a series of experimental findings, including (a) the reversal in flexor and extensor phase durations with an increase of locomotor speed during left-right symmetric locomotion, and (b) the maintenance of step cycle period during split-belt locomotion due to adjustment of the flexor duty cycle. The results of these simulations provide strong support for the proposed organization and operation of spinal locomotor circuits.

Organization of Left-Right Commissural Interactions in the Spinal Cord: The Role of V3-Mediated Commissural Pathways
In the present model, the interactions between left and right RGs were based on the model by Danner et al. (2019). Importantly, that model was derived from experiments on symmetric (bilateral) and asymmetric (unilateral) optogenetic stimulations of commissural V3 neurons involved in leftright coordination performed in the same study. Interestingly, unilateral stimulation produced effects that were qualitatively similar to some features of split-belt locomotion. They provided strong evidence that spinal V3 CINs are involved in leftright limb coordination via two pathways: through mutual excitation between the left and right extensor half centers of the RGs and, importantly, via crossed inhibition from extensor half-centers to contralateral flexor half centers through an additional inhibitory interneuron population (presumably V1) (see Figure 2A). In the present study, we show that the commissural inhibition of flexor half-centers by the contralateral extensor half-centers (see Figure 5 and related texts) is critically important for the stability of anti-phase flexor oscillations at low frequencies in symmetric conditions, which corresponds to a normal locomotor pattern. Therefore, our study provides additional support for the important role of V3 CINs and the existence of inhibitory commissural pathways from extensor half-centers to contralateral flexor half-centers, mediated by V3 and (presumably) V1 interneurons . Although this prediction still awaits experimental testing, crossed inhibition to flexors (by afferent stimulation) has been observed in anesthetized preparations (Jankowska et al., 2005;Jankowska and Edgley, 2010) and during locomotion in intact cats  as well as in mouse (Laflamme and Akay, 2018) and human (Mrachacz-Kersting et al., 2017) studies.
In summary, our analysis of the model allowed us to evaluate the specific roles of the two types of inhibitory commissural interactions (called here flexor-flexor and extensorflexor inhibition) in left-right coordination. The flexor-flexor inhibition, presumably mediated by V0 CINs (Talpalar et al., 2013;Shevtsova et al., 2015), supports left-right alternation and its weakening may stabilize left-right in-phase synchronization. The extensor-flexor inhibition, presumably mediated by V3 CINs and V1 interneurons , ensures that left and right activities alternate in a strict out-of-phase manner in symmetric conditions.

Insights From Symmetric Locomotion
It is well known that during normal locomotion in cats and humans, an increase of speed is accompanied by a significant reduction of stance phase duration with or without a minor reduction of swing phase duration (Grillner et al., 1981;Halbertsma, 1983;Frigon and Gossard, 2009;Frigon et al., 2013Frigon et al., , 2014Frigon et al., , 2017) (see also Figure 3). This observation seems to support the flexor-driven concept of locomotor rhythm generation. However, in intact and spinal cats, increasing locomotor speed produces a more balanced pattern, with stance duration approaching and even becoming shorter than swing duration. This is clearly observed during overground locomotion in intact cats ( Figure 3A). We suggest that when approaching the point of equality between phases, rhythmogenesis shifts toward the classical half-center mechanism. The observation indirectly supporting this view is that after the point of equality, the oscillation period (and hence the frequency) saturates and does not change much, which is a typical feature of classical half-center dynamics (Daun et al., 2009;Ausborn et al., 2018).

Insights From Asymmetric Split-Belt Treadmill Locomotion
Previous experimental studies in cats using split-belt treadmill locomotion demonstrated that the mammalian spinal cord has a remarkable adaptive capacity for left-right coordination, from simple to extreme conditions (Forssberg et al., 1980;Halbertsma, 1983;Frigon et al., 2013Frigon et al., , 2015Frigon et al., , 2017. In simple conditions, with slow/fast speed ratios of up to 1:2.5 (0.4:1.0 m/s), animals maintain the period of oscillations (and frequency) almost unchanged and compensate for the reduction of stance phase duration on the fast belt by a corresponding increase of the duration of the swing phase ; see Figure 6. Our model was able to reproduce this feature specifically due to the implementation of our suggestion, that increased activation of the flexor half-center in each RG is accompanied by a reduction in the activity of the corresponding extensor half-center. This implementation leads to a switch in the rhythmogenic mechanism of the fast RG from flexor-driven oscillations to the classical half-center mechanism (Figures 8A,B). Removing this feature from the model leads to constant swing duration accompanied by a noticeable increase of oscillation frequency in both limbs (RGs) with increasing drive to the flexor half-centers (Figures 8C,D), contradicting the experimental results, shown in Figure 6.
Experimental studies of cat locomotion on split-belt treadmills in extreme conditions, with slow/fast speed ratios of 1:3 and more Kuczynski et al., 2017) showed that cats use a specific strategy to stabilize locomotion by taking multiple steps on the fast side per step on the slow side. Moreover, although there was some variability between animals, both intact (Kuczynski et al., 2017) and spinal  cats exhibit 1:2, 1:3 or 1:4 coordination patterns corresponding to 2, 3, or 4 steps on the fast side per step on the slow side, respectively. To simulate these behaviors, we applied different drives to the left and right RGs in the model, assuming that these conditions are qualitatively similar to the extreme case of split-belt locomotion. Under these conditions, the model predicts that the number of different coordination patterns depends on the value of the drive to the slow RG (Figure 7). For relatively high drives to the slow RG (>0.45), only a 1:1 coordination pattern is possible, which corresponds to simple conditions in split-belt locomotion (see above). However, if the drive to the slow RG is smaller, 1:2 + coordination patterns become possible. For example, for a slow RG drive value of 0.4, as the drive to the fast RG increases, there is a transition from 1:1 to 1:2 coordination pattern, but no 1:3 regime exists, while for a slow RG drive value of 0.25, as the fast RG drive progressively increases, the system undergoes 1:1, 1:2, 1:3, and 1:4 regimes. Qualitatively similar behavior is observed in extreme split-belt locomotion where in order to achieve higher order coordination patterns, one has to set lower speeds of the slow belt.

Limitations, Functional Considerations, and Future Directions
In this study, we show that a relatively simple functional connectome between populations of interneurons providing output to flexor and extensor motoneurons that control a pair of limbs can explain a variety of coordination patterns emerging in split-belt experiments. The mathematical model we developed allowed us to formulate a novel hypothesis about general mechanisms of locomotor phase duration control suggesting that variation of the excitatory drive to the flexor half-centers is accompanied by an opposite change in the drive to the extensor half-centers. However, our model does not provide any specifics on neuronal pathways mediating these interactions.
What would be the benefit of switching from a flexor-driven RG operation to a classic half-center mode with increasing speed? Although we can only speculate, the goal of the spinal locomotor network might be to optimize efficiency or balance (avoid falling). At slow to moderate speeds, the stance duration is long and inputs from group I/II extensor muscle afferents and paw pad cutaneous afferents have a relatively long time to regulate stance duration and adjust/correct for destabilizing perturbations. Thus, at slow to moderate speeds, a flexor-driven RG mode is less costly and more efficient. However, as speed increases, stance duration decreases and afferent inputs do not have as much time to adjust or correct for postural perturbations. As such, at high speeds, a classic half-center mode, whereby both stance and swing phase durations are balanced, becomes more efficient to avoid falling, as each phase can be more flexibly controlled.
Locomotion in mammals results from a complex interplay between spinal CPGs, descending commands from the brain and sensory feedback from the limbs and trunk (Rossignol et al., 2006). The observation that our experimental results were similar in intact and spinal cats during tied-belt and split-belt locomotion indicates that the control of cycle and phase durations is mainly mediated by spinal CPGs interacting with sensory feedback from the limbs. As the biomechanics of the limbs change throughout the step cycle (e.g., muscle stretch and contractions, contacts and liftoffs), phasic sensory feedback also changes and different inputs affect the step cycle and its structure at different time points (Rossignol et al., 2006;Frigon et al., 2017). These phasic sensory inputs strongly affect the operation of spinal locomotor CPGs.
Considering that similar coordination patterns are observed in split-belt experiments in both intact and spinal cats , it is reasonable to assume that drives controlling left and right RGs depend on sensory feedback rather than on supraspinal inputs. One obvious source of sensory feedback is muscle afferent inputs that are known to affect the dynamics of the spinal locomotor CPG circuits (see Markin et al., 2010 for review). Our model does not explicitly account for this type of feedback. Therefore, the functional interactions and intrinsic flexor and extensor half-centers' oscillatory properties can be defined in part by inputs from somatosensory afferents. Another type of sensory feedback known to influence locomotion is from the skin . Cutaneous feedback modulation by paw anesthesia alters margins of stability during split-belt cat locomotion (Park et al., 2019). It was recently suggested that this alteration occurs due to misrepresentation of the center of mass in the cat's balance control system after disrupting cutaneous feedback from the paws (Latash et al., 2020). Altogether, the balance control system (or some of its elements) and locomotor pattern generation may interact at the spinal level, which opens new ways to mathematically model these interactions and thus generate new hypotheses about neuronal pathways mapping somatosensory afferents to the spinal locomotor circuits. Decomposing the functional interactions between left and right RGs into components mediated by local commissural interneurons and spinal reflexes can be a major future research direction where mathematical modeling proves instrumental.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.