Tapping Into the Human Spinal Locomotor Centres With Transspinal Stimulation

Human locomotion is controlled by spinal neuronal networks of similar properties, function, and organization to those described in animals. Transspinal stimulation affects the spinal locomotor networks and is used to improve standing and walking ability in paralyzed people. However, the function of locomotor centers during transspinal stimulation at different frequencies and intensities is not known. Here, we document the 3D joint kinematics and spatiotemporal gait characteristics during transspinal stimulation at 15, 30, and 50 Hz at sub-threshold and supra-threshold stimulation intensities. We document the temporal structure of gait patterns, dynamic stability of joint movements over stride-to-stride fluctuations, and limb coordination during walking at a self-selected speed in healthy subjects. We found that transspinal stimulation 1) affects the kinematics of the hip, knee, and ankle joints, 2) promotes a more stable coordination at the left ankle, 3) improves interlimb coordination of the thighs, 4) improves intralimb coordination between thigh and foot, 5) promotes greater dynamic stability of the hips, and lastly 6) affects the mechanical stability of the joints. These results support that transspinal stimulation is an important neuromodulatory strategy that directly affects gait symmetry and dynamic stability. The conservation of main effects at different frequencies and intensities calls for systematic investigation of stimulation protocols for clinical applications.


Introduction
Spinal neuronal networks can generate exible and adaptable rhythmic motor activity in absence of descending control and movement mediated afferent inputs [1][2][3] .These spinal neuronal networks are known as central pattern generators (CPGs) and have been extensively studied in animals under different preparations including genetic transcription factors and computational modeling [4][5][6][7] .Seminal works on the existence of mammal CPGs were the rhythmic motor discharges in a decerebrate, spinalized and deafferented cat postulated by Thomas Graham Brown (1882-1965) and his proposal of "half-centres" on each side of the spinal cord followed by Anders Lundberg (1920-2009) and Elzbieta Jankowska (1930-)   works on the physiology and function of the spinal interneuronal networks and their control on stepping [8][9][10] .Human CPGs have been a black box for decades, giving rise to great debates regarding their existence and function 11 .Direct evidence in humans is di cult, but the similarities between the long-latency exor re ex in l-3,4-dihydroxyphenylalanine treated spinalized animals and humans 12,13 , functional organization of interneurons mediating sensory feedback [14][15][16] , and spontaneous rhythmic leg movements in paralyzed persons 17,18 support for a spinal control of stepping in humans.The deep similar modulation pattern of the soleus stretch and H-re ex in humans and mesencephalic cats during walking 19,20 , suggest that quadrupedal and bipedal locomotion share common neuronal networks regardless of the evolution of gait.
In the human isolated spinal cord, step-like phasic motor activity is generated following tonic epidural stimulation at 25-60 Hz 21 .Epidural stimulation produces either locomotor-like patterns, synchronous activation of all limb muscles, or different combinations of both 22,23 .Epidural stimulation at 5 to 15 Hz and 25 to 50 Hz generates either co-contraction or alternated activity between antagonists supporting for activation of neuronal networks for standing and walking based on frequency 23 .Similarly, transcutaneous spinal cord (termed here transspinal) stimulation over the thoracolumbar enlargement can provide tonic excitatory inputs to spinal locomotor circuits.Transspinal stimulation at 30 Hz produces step-like movements based on muscle activation patterns in people with intact or partially injured spinal cord when legs were held in a gravity neutral position or during walking [24][25][26] .Transspinal stimulation with a single 1 ms pulse at supra-threshold intensities has direct actions on human spinal locomotor networks based on the profound suppression of both exor and extensor re exes and electromyographic (EMG) activity of knee and ankle muscles during walking [27][28][29] .
The function of locomotor centers during transspinal stimulation at different frequencies and intensities has never been systematically investigated before.This lack of knowledge limits our understanding on the underlying mechanisms of actions of transspinal stimulation.In this study we established the function of locomotor networks during transspinal stimulation at 15, 30, and 50 Hz at sub-threshold and supra-threshold transspinal stimulation intensities via 3D joint kinematics.Speci cally, we determined 1) the temporal structure of gait patterns variability, 2) dynamic stability of joint movements over stride-tostride uctuations, 3) forward walking mechanical stability, and 4) intralimb-interlimb coordination computed as the spatiotemporal pattern between limb segments when transspinal stimulation was delivered at 15, 30, and 50 Hz at sub-and supra-threshold stimulation intensities during walking at a selfselected speed in healthy subjects.We hypothesized that transspinal stimulation, regardless of frequency and intensity, promotes physiological step progression and joint movements.

Participants
Ten physically active subjects (5 female; mean age 26.7 ± 3.7 years; mean height 1.70 ± 0.2 m; and mean mass 72.7 ± 24.2 kg) with no history of neuromusculoskeletal disorder participated in the study.The experimental protocol was approved by the City University of New York (CUNY) Institutional Review Board Committee (IRB No. 2022-0003-CSI).All experimental procedures were conducted in compliance with the Declaration of Helsinki and the CUNY IRB-wide regulations and guidelines.All participants provided written informed consent before study enrolment and participation.

Transspinal stimulation
With subjects seated, the Thoracic 10 spinous process was identi ed via palpation of spinal processes and anatomical landmarks.A single reusable self-adhered cathode electrode (10.2 ´ 5.1 cm 2 , Uni-Patch, Massachusetts, USA) was placed along the spinal processes, and covered from Thoracic 10 to Lumbar 1-2 vertebral levels.A pair of interconnected anode electrodes (same type as the cathode) was placed on either side of the iliac crests [30][31][32] .With subjects' supine, hip/knee joints exed at 30° and legs maintained in midline via external support, single 1-ms pulses were delivered by a constant current stimulator (DS7A, Digitimer, Welwyn Garden City, UK) triggered by Spike 2 scripts (1401 plus, Cambridge Electronics Design Ltd., Cambridge, UK).The position of the cathodal stimulating electrode was based on the latency, shape, and amplitude of the right and left soleus transspinal evoked potential (TEP) at increasing intensities observed on a digital oscilloscope (Tektronix, USA).After the optimal location was identi ed, the electrode was a xed to the skin via Tegaderm transparent lm (3M Healthcare, Minnesota, USA), and maintained under pressure via a custom-made pad.For each subject while supine, the stimulation intensity that the right soleus TEP was rst noted on the oscilloscope was termed as TEP threshold 30,33 .

Experimental data collection
With subjects standing upright on a treadmill, 15 retro-re ective markers of 12.5 mm diameter were placed on the anterior-superior iliac crest, sacrum, medial and lateral femoral epicondyle, medial and lateral malleolus, calcaneal tuberosity, and mid-point between the second and third metatarsals.Markers were used to de ne a mechanical model of pelvis, right and left thigh, shank, and foot rigid body segments (Figure 1A).Rigid clusters of 3 or 4 retro-re ective markers were placed on the thigh, shank, and foot segments away from the anatomical landmarks to track the movement of these segments during walking 34 .The three pelvic anatomical markers were used for tracking the pelvis segment.A static calibration trial while standing on the treadmill was collected while the anatomical markers and rigid clusters were attached.Then, we removed the anatomical retro-re ective markers that were not needed, and each subject walked on a motorized treadmill (Gait Trainer, Biodex Medical Systems, Shirley, NY, USA) at their self-selected speed (3.9 ± 0.5 km/h).For each subject, 7 walking bouts of 10-min duration each were completed.Transspinal stimulation was delivered with a DS8R stimulator (Digitimer Ltd., UK).The rst session was walking on the motorized treadmill without transspinal stimulation, followed by randomized sessions during which transspinal stimulation was delivered at 15, 30 and 50 Hz with 10 kHz carry over frequency at sub-threshold (80.2 ± 27.6 mA) and supra-threshold (169 ± 47.8 mA) stimulation intensities.Intensities and frequencies for transspinal stimulation were randomly assigned within and across subjects.Subjects were free to move their arms while walking on the motorized treadmill.

3-Dimensional (3D) Kinematics
The coordinates of the markers at 150 Hz sampling frequency were collected with an 8-camera optoelectronic system (Qualisys Miqus M1, Gothenburg, Sweden).Calibration residuals were on the order of < 1 mm for an active volume of approximately 4 × 2 × 2 m 3 .The 3D coordinates were reconstructed, interpolated, and labelled in Qualisys Track Manager (Qualisys, Gothenburg, Sweden).The static calibration and walking trials were then exported as C3D les and further processed using Visual3D (C-Motion Inc., Rockville, MD).The generalized cross-validation quintic splines were used for smoothing marker coordinates 35 .A mechanical model with 6 degrees-of-freedom at the hip, knee and ankle joints was created in Visual3D using the static calibration trial.The hip joint center was computed from the pelvic markers using regression equations 36 .The knee and ankle joint centers were estimated based on the kinematic analysis of marker clusters located on each joints' distal segment during a walking trial 37 .
The knee and ankle joint centers were located at the midpoint between the medial and lateral femoral epicondyle, and medial and lateral malleolus markers, respectively, both aligned with the functional rotation axis of the joints.Segment coordinate systems were constructed according to the International Society of Biomechanics (ISB) standards, with the exception that we used the projected anatomical markers to the functional axis to de ne the frontal plane 38 .The ankle joint was zero degrees during the static calibration trial, referenced to the shank coordinate system.The joint movements were described with Cardan's XYZ angular convention in accordance with the ISB standards 39 , which is equivalent to the Joint Coordinate System 40 .The medio-lateral (X) axis of the proximal segment speci ed the exion and extension at the joint, the oating (Y) axis speci ed the adduction and abduction and was de ned as the cross-product of the exion-extension axis with the longitudinal axis, and the longitudinal (Z) axis of the distal segment speci ed the axial rotation.Heel contact and toe-off during walking were determined based on the kinematics data 41 .

Detrended uctuation analysis
Gait spatiotemporal time series were analyzed using fractal analysis techniques to establish to what extent the temporal structure of gait patterns was affected by transspinal stimulation.Step length and step width times series (Figure 1B) during treadmill walking were calculated as the anteroposterior or mediolateral distance between two successive footprints, measured from the heel contact of the right foot to the heel contact of the left foot.Standard deviation quanti ed the amount of step width and step length variability.The temporal structure in step width and step length variability in terms of persistence in uctuations was quanti ed using the fractal scaling exponent (α-DFA; Figure 1C), which was obtained from the detrended uctuation analysis (DFA) using the average evenly spaced algorithm 42 .The minimum length of the analyzed time series (N = 884 steps) and the window size range [16, N / 9] was consistent with previous recommendations for DFA calculation 43 .
The α-DFA value quanti ed the degree of statistical persistence (0.5 < α-DFA < 1.0) or anti-persistence (0 < α-DFA < 0.5) in the variability of step length and step width intervals.In a persistent correlation, an increase or decrease of step length or step width is more likely to be followed by another increase or decrease in the same direction, which indicates a long-range correlation.In contrast, in an anti-persistent correlation, a decrease is more likely to follow an increase and vice versa.A time series with α-DFA = 0.5 indicates that the uctuations in the time series are uncorrelated, implying that the pattern of the time series cannot be predicted from its past values.Even if the time series is randomly permutated, it will still produce an α-DFA ≈ 0.5.
To con rm the temporal structure in the step length or step width time series, a surrogate data analysis 44 based on permutation test was conducted to test the null hypothesis that the step time or step length uctuations were indistinguishable from independent and identically distributed (IID) noise.Rejection of the null hypothesis would indicate that uctuations exhibited a deterministic temporal structure.The surrogate time series were generated by random permutation of the original data and the α-DFA value was calculated from each surrogate to construct the distribution of the α-DFA value under the IID H 0 .The α-DFA values from the original time series and from the surrogates were compared with a two-sided rankorder test at signi cance level of α = 0.05.To compensate for multiple comparisons the signi cance level was set at α adjusted = 0.005 (Bonferroni correction).For each original time series β = (2 / α adjusted ) -1 = 399 surrogates were generated.The p-values for the surrogate test were calculated by ranking the α-DFA values (399 surrogates + 1 original = 400 values) in numerical order and determining the rank of the DFA value of the original time series (p-value = rank / 400).

Largest Lyapunov exponent
The time series of the hip, knee and ankle joint exion and extension movements were analyzed using nonlinear analysis techniques to establish to what extent the dynamic stability of joint movements over stride-to-stride uctuations was affected by transspinal stimulation during walking (Figure 1D).The temporal structure of the joint movement variability in terms of dynamic stability over stride-to-stride uctuations was investigated using the largest Lyapunov exponent (LyE).The largest LyE measures the average exponential rate of divergence of the attractor's phase trajectory (orbit) in theoretic-information terms (bits/sec).A greater rate of orbital divergence, re ected by higher positive values of the largest LyE, re ects more unstable behavior of a phase trajectory as it is repeated over time.This suggests a greater likelihood of transitioning to a more stable attractor state.A time series with largest LyE = 0 bits/sec indicates a perfectly stable, repeating pattern, which is unlikely in biological repeated movements.The reconstruction of the attractor from the joint movement time series in an m-dimensional space was based on Takens embedding theorem 45 .The appropriate delay time (τ) and embedding dimension (m) were determined for each time series using the average mutual information and global false nearest neighbor methods, respectively 46,47 .The Wolf algorithm was used to calculate the largest LyE of the respective time series 48 .In agreement with previous recommendations, the analyzed time series had a xed number of data points across participants (N = 90,000 data points; minimum strides = 442) 49 .For the Lyapunov analysis, the number of points for time evolution was set to 3, the maximum angular orientation error between evolved and nearest neighbor points was set to 0.3 radians, the minimum distance for the selection of new nearest neighbor point was set to 0.0001, and the maximum distance for the selection of a new nearest neighbor point was set to 0.1 times the maximum length of the attractor 48,50 .

Continuous Relative Phase
The interlimb spatiotemporal pattern between the left and right thigh, and the intralimb spatiotemporal pattern between the thigh and shank, shank and foot, and thigh and foot segments were analyzed with the continuous relative phase (CRP) to establish to what extent interlimb and intralimb coordination was affected by transspinal stimulation (Figure 1E).CRP quanti es the time shift of two segment angles S A and S B .To model interlimb and intralimb coordination, the sagittal plane segment angles of the thigh, shank, and foot body segments were calculated relative to the right horizontal axis (Figure 1A).Segment angles were time normalized to 101 points for the gait cycle, centered around zero, and Hilbert transformed to generate analytic signals H A and H B 51 .The CRP between the two segment angles was calculated as the arctangent of the product of H A with the conjugate of H B .CRP was unwrapped to resolve any discontinuities and adjusted to ±180° for intralimb and 0º-360º for interlimb coordination.CPR from 0° or 360º to ±180° indicates an in-phase towards an out-of-phase coordination.For the intralimb coordination, a negative CRP slope indicates that the proximal segment moves faster than the distal segment in phase space, while a positive CRP indicates that the distal segment moves faster.For the interlimb coordination, a negative CRP slope indicates that the right segment moves faster than the left segment in phase space, while a positive CRP indicates the opposite.The mean absolute relative phase (MARP) quanti ed the magnitude of in-phase and out-of-phase coordination.MARP was calculated by averaging the absolute CRP across gait cycles.Deviation phase (DP) quanti ed coordination variability.DP was calculated as the standard deviation of the CRP across gait cycles 52 .Custom-written MATLAB scripts (R2021a, MathWorks, Inc., Natick, USA) were used for data analyses.

Margin of stability
The extrapolated center of mass (xCOM) was used to determine the margin of stability in the anteroposterior direction during walking at heel contact and at toe off (Figures 1F and 1G) to establish to what extend mechanical stability was affected by transspinal stimulation 53 .The margin of stability was de ned as the anteroposterior distance between the xCOM and heel position of the leading leg at heel contact and at toe off.Changes in margin of stability during walking indicate that additional motor actions are required to preserve mechanical stability.

Statistical analysis
The effects of 15, 30, and 50 Hz transspinal stimulation at sub-and supra-threshold intensities on joint angles, the amount (standard deviation) and temporal structure (α-DFA) of the step length and step width variability, the dynamic stability (LyE) of joint angles, and interlimb/intralimb coordination (MARP, DP) during walking were established.Linear mixed-effects models followed by a generalized linear hypothesis test for planned contrasts were used to compare transspinal stimulation groups with the control group, as well as their interaction (sub-threshold vs. supra-threshold, 15 vs. 30 Hz, 15 vs. 50 Hz, 30 vs. 50 Hz) 54,55 .Linear mixed-effects models were tted by restricted maximum likelihood with subject as the random effect.The Westfall method was used to control the familywise type-I error rate.Cohen's d (d < 0.20 = very small, d < 0.50 = small, d < 0.80 = moderate, otherwise interpreted as a large effect) was used to report the effect size for planned contrasts.Before model tting, the data were preprocessed using Tukey's boxplot method to identify and remove statistical outliers.The assumptions of the linear mixedeffects models were evaluated using graphical procedures.Variance homoskedasticity was assessed by examining plots of standardized residuals against tted values, while residuals normality was assessed using quantile-quantile plots.The statistical analysis was conducted in R 56 .The signi cance level for all tests was set to α = 0.05.

Joint Kinematics
Effects of transspinal stimulation on ankle joint movements Average left ankle joint kinematics from all subjects during treadmill walking for dorsi exion (Figure 2A), inversion (Figure 2B), and adduction (Figure 2C) are shown for control conditions and during transspinal stimulation.The frequency and intensity that transspinal stimulation had signi cant effects on the left ankle joint movement compared to control walking are shown in Figure 2D.Ankle eversion increased at loading response, dorsi exion increased at terminal stance and pre-swing, inversion increased at pre-swing, and plantar exion decreased at initial swing (Figure 2E).The speci c gait phases during which the ankle movement was signi cantly affected are summarized in

Coordination Stability
Transspinal stimulation at all frequencies and intensities, except 50 Hz sub-threshold, signi cantly decreased left DP Shank-Foot at initial swing compared to control walking (p adj < 0.05).Furthermore, left DP Shank-Foot decreased at mid-swing with 30 Hz and 50 Hz supra-threshold (p adj < 0.05).These results suggest that transspinal stimulation promoted a more stable coordination pattern at the left ankle joint during the swing phase compared to control walking.Additionally, DP Thigh decreased at pre-swing with 15 Hz sub-threshold (t 48 = -3.18,p adj = 0.02, d = 0.46), indicating a more stable interlimb coordination during push off compared to control walking.

Interlimb Coordination
The CRP Thigh was closer to 180º, indicating maintenance of the reciprocal out-of-phase coordination of the limbs during walking (Figure 8A).At the beginning of the terminal stance that corresponds to the start of the push off phase, the left trailing thigh started to move faster than the right thigh termed as reversal.
A reversal was also present at the beginning of the terminal swing which corresponds to 90% of the gait cycle.At this reversal point, the right (trailing) thigh moved faster than the left thigh.The CRP Thigh curve was symmetric with respect to the double support interval, whereas legs reach an out-of-phase (180º) coordination.The relative motion of the legs produced reversals during either left or right leg push off, which corresponds to the loading response and preswing phases, respectively.The timing of the reversals in the coordination dynamics of CRP Thigh (timing of maxima and minima) was similar between transspinal stimulation and control walking.Transspinal stimulation decreased CRP Thigh at the rst reversal (Figure 8A).
Statistically signi cant MARP Thigh differences between transspinal stimulation and control walking were revealed at all frequencies and intensities except with 15 Hz sub-threshold (Figure 8B).MARP Thigh decreased from initial contact to terminal swing and increased at pre-swing, indicating a less out-ofphase coordination compared to control walking -i.e., further away from 180º (Figure 8C).

Intralimb Coordination
Effects of transspinal stimulation on thigh-foot coordination The CRP Thigh-Foot was closer to 0º at early stance indicating that thigh and foot exhibited a relationship closer to being in-phase (Figures 9A-9C).During initial contact the foot was moving faster than the thigh in phase space (positive slope).A reversal occurred during the loading response at about 8% of the gait cycle, where the foot achieved its maximum lead over the thigh.At this reversal point, the foot started to lose its lead, with the thigh moving faster than the foot in phase space up to pre-swing.Another reversal occurred at pre-swing at about 55% of the gait cycle, and the foot began to regain its lead over the thigh.
The thigh-foot coordination became closer to being in-phase again in the late swing phase.The timing of the reversals in the coordination dynamics of CRP Thigh-Foot (timing of maxima and minima) was similar between transspinal stimulation and control walking.Transspinal stimulation changed the magnitude of CRP Thigh-Foot when the reversals occurred (Figures 9B-9D).
Statistically signi cant differences on MARP Thigh-Foot between transspinal stimulation and control walking were noted at all frequencies and intensities (Figure 9E).MARP Thigh-Foot decreased at initial contact and loading response indicating a more in-phase coordination compared to control walking and increased from terminal stance to mid-swing indicating a more out-of-phase coordination (Figure 9F).These results suggest that transspinal stimulation promotes an out-of-phase thigh-foot coordination during push off that is consistent with the hip and ankle joint kinematics we observed here ensuring adequate step progression and clearance of the foot to enter a synchronized swing phase.Effects of transspinal stimulation on shank-foot coordination The CRP Shank-Foot was closer to 0º indicating that the shank and foot exhibited an in-phase relationship at early stance (Figures 10A and 10C).The shank was moving faster than the foot in phase space (negative slope) up to preswing, whereas a reversal occurred at about 50% of the gait cycle.At the reversal point the foot started to move faster than the shank (positive slope) until it reached a maximum at about 80% of the gait cycle.The shank-foot coordination became closer to being in-phase again at mid-swing.
Transspinal stimulation increased the magnitude of the negative slope during which the shank was moving faster compared to control walking at about 12% of the gait cycle which lead to greater CRP Shank- Foot when the reversal occurred compared to control walking (Figures 10B and 10D).The timing of the reversal in the coordination dynamics of CRP Shank-Foot (timing of minima) was similar between transspinal stimulation and control walking.
Statistically signi cant differences on MARP Shank-Foot between transspinal stimulation and control walking were found at all frequencies and intensities (Figure 10E).MARP Shank-Foot was signi cantly increased at loading response and from terminal stance to initial swing indicating a more out-of-phase coordination compared to control walking and decreased at mid-swing and terminal swing indicating a more in-phase coordination (Figure 10F).Additionally, transspinal stimulation promoted an in-phase coordination at ankle joint during mid-and terminal swing to prepare for heel contact.

Effect of transspinal stimulation on thigh-shank coordination
The CRP Thigh-Shank was closer to 0º indicating that thigh and shank exhibited a relationship closer to being in-phase at early stance (Figures 11A and 11C).A reversal occurred near the beginning of the gait cycle at about 12 % of the gait cycle, where the shank achieved its maximum lead over the thigh.At this reversal point, the shank started to lose its lead, with the thigh moving faster than the shank in phase space up to initial swing.Another reversal occurred at the initial swing, at about 69% of the gait cycle, and the shank began to regain its lead over the thigh (Figures 11B and 11D.The thigh-shank coordination became closer to being in-phase again in the late swing phase.The timing of the reversals in the coordination dynamics of CRP Thigh-Shank (timing of maxima and minima) was similar between transspinal stimulation and control walking.Transspinal stimulation decreased CRP Thigh-Shank at the rst reversal (Figures 11A and 11C).
Statistically signi cant differences on the MARP Thigh-Shank between transspinal stimulation and control walking were found at all frequencies and intensities (Figure 11E).MARP Thigh-Shank decreased from initial contact to midstance (where the rst reversal occurred), indicating a more in-phase coordination compared to control walking, and increased at mid-swing indicating a more out-of-phase coordination (Figure 11F).These results suggest that transspinal stimulation promoted an out of phase coordination at knee joint during swing phase to increase foot clearance and advancement of the swinging leg.These results are consistent with the knee and ankle joint kinematics we observed in this study, which showed an increased knee exion at mid-swing resulting from the excessive lead of shank over foot during push off and initial swing.Additionally, transspinal stimulation promoted an in-phase coordination at knee joint during initial contact and loading response, which is required to decelerate center of mass movement (power absorption) and maintain forward propulsion.15D).These results suggest that transspinal stimulation promoted the activity required to continue walking, which is consistent with the increased out-of-phase found at thigh-foot coordination during push off (Figures 9B and 9D), the increased out-of-phase at ankle joint during push off by slowing down the foot rotation (Figures 10B and 10D), and the increased dorsi exion at terminal stance and decreased plantar exion at preswing found at ankle joint kinematics results (Figures 3 and 4).
The results suggest that mechanical stability was affected most by transspinal stimulation delivered at 30 Hz.

Discussion
This is the rst systematic study that provides evidence on how transspinal stimulation at different frequencies and intensities affects the human locomotor rhythm.Transspinal stimulation over the thoracolumbar region at sub-threshold and supra-threshold paresthesia intensities with 15, 30, and 50 Hz frequencies 1) affected the kinematics of the hip, knee, and ankle joints, 2) promoted a more stable coordination at the left ankle, 3) improved interlimb coordination of the thighs, 4) improved intralimb coordination between thigh and foot, 5) promoted greater dynamic stability of the hips, and lastly 6) affected the mechanical stability of walking.Our novel ndings support that transspinal stimulation is an important neuromodulatory strategy that directly affects gait symmetry and dynamic stability.
Transspinal stimulation at 30 Hz suprathreshold decreased the margin of stability during walking (Figures 15A-C).The margin of stability is a measure that re ects the relationship between the center of mass position and velocity relative to foot placement and is directly related to mechanical stability of gait 53,57 .The margin of stability predicts the point on the ground that the foot will be placed before the body comes to a standstill.While walking at steady speed, a forward-stable gait without excessive spatiotemporal uctuations is achieved by placing the foot behind that predicted point ensuring negative margin of stability, for example, through constant offset control 53 .Less negative margin of stability during walking with transspinal stimulation indicates reduction in mechanical stability (less resistance to center of mass redirection) and facilitation of neural control because less volitional mechanical impulse would be required to redirect the center of mass velocity compared to control walking 58 .Redirection of center of mass velocity during steady speed walking requires negative work at the ankle and hip to slow down the center of mass at terminal stance followed by positive work at the ankle at pre-swing 59,60 .
Ankle musculature generates minimal power at terminal stance, whereas negative mechanical work is primarily performed, and elastic strain energy is stored in the tendinous tissues 61,62 .A possible way to amplify ankle power at push off is through direct effects on the functional operation of the CPG by increasing triceps surae muscle activity through locomotor Ib facilitatory pathways exerted from medial gastrocnemius onto soleus motoneurons [63][64][65] .During walking, whereas a stretch-recoil cycle is performed by the tendinous tissues in ankle extensors 66 , prolonged dorsi exion and thus tibial progression before push off can further increase the elastic recoil and thus amplify power output 67 .
Indeed, increased out-of-phase thigh-foot coordination during push off (Figures 9B, 9D) coincided concomitantly with increased out-of-phase at ankle joint (Figures 10B, 10D), slowing down the foot rotation, and increased dorsi exion at terminal stance (Figures 3, 4).Speci cally, in the left ankle transspinal stimulation mostly with 30 Hz sub-threshold and supra-threshold facilitated walking by increasing the ankle eversion at loading response, dorsi exion at terminal stance and pre-swing, and decreasing plantar exion at initial swing (Figure 2E).Accordingly, it seems that transspinal stimulation promotes ankle control strategies linked to the extensor CPG.
Push off may also be facilitated by changes in the hip joint.The biarticular nature of muscles spanning the hip and knee joints, and the monosynaptic heteronymous connections between thigh exor afferents and ankle extensor muscles 68 , support the physiological mechanism underlying the kinematics changes we observed here.The right hip exion increased at loading response and mid-stance while hip external rotation increased at initial and mid-swing with 30 Hz regardless of intensity.Moreover, the increased knee exion at early and mid-swing phases with 30 Hz suprathreshold supports further the progression of the step cycle, and that transspinal stimulation may have access to both the exor and extensor CPGs.
Interlimb (right-left thighs) and intralimb (shank-foot, thigh-foot) coordination provide signi cant information on the function of the locomotor rhythm.Coordination between joints of the same leg, and between joints of both legs is controlled by commissural interneurons with their axons projecting to the contralateral side of the spinal cord 69 , and propriospinal neurons that modulate sensory information and carry temporal information for gait control [70][71][72] .These spinal neuronal circuits organize the limb movements and thereby transform rotations in joint into the desired walking patterns.Transspinal stimulation at 30 Hz supra-threshold promoted an out-of-phase thigh-foot coordination from terminal stance to mid-swing including push-off, and out-of-phase shank-foot coordination at loading response and from terminal stance to initial swing compared to control walking (Figure 10F).Transspinal stimulation decreased the mean absolute relative phase of the thigh from initial contact to terminal swing and increased at pre-swing, indicating a less out-of-phase coordination compared to control walking.These effects were apparent with 50 Hz sub-threshold and 15 Hz supra-threshold, supporting for variable effects on intralimb coordination based on the intensity and frequency of transspinal stimulation.
Transspinal stimulation increased out-of-phase coordination at the ankle joint during push off by slowing down the foot rotation, which is consistent with an increased dorsi exion at terminal stance and decreased plantar exion at pre-swing found at ankle joint kinematics results (Figures 3 and 4).These results further support that transspinal stimulation affects locomotor gait patterns by affecting interlimb and intralimb coordination and thus may be useful in neurological disorders that limb coordination is impaired.
Important spatiotemporal gait parameters are that of step length and step width variability, which re ect the overall physiological function of the nervous system and become pathological in neurological disorders.30 Hz sub-threshold and 30/50 Hz supra-threshold transspinal stimulation decreased step length variability and increased that of step width (Figure 13).These changes suggest that transspinal stimulation contributes to greater passive stability in the forward direction, accompanied by an increased reliance on supraspinal control for lateral stability.In silico simulations, physical biped-legged models, and experiments involving human walking, revealed that foot placement in forward direction is dependent mostly on the biomechanical properties of the body and their mechanical interaction with the environment 57,[73][74][75][76][77] .On the other hand, lateral foot placement relies on an active supraspinal/spinal sensory feedback control 57,[73][74][75][76][77] .Greater step width variability indicates increased effort of active control required for lateral stabilization 73,74,78,79 .The decrease in step length variability indicates greater contribution of passive stabilization of forward dynamics involving spinal control and simple re ex actions.The decrease on Lyapunov exponents on right hip joint corroborates the greater forward stability with 50 Hz supra-threshold stimulation.The increase in step width variability with transspinal stimulation could indicate a reduction of the amount of sensory information available, probably a noisy center of mass state information [80][81][82] , which is pertinent to both push off and foot placement feedback control.
Based on this theoretical framework, a possible push off modulation did not decrease the need for active control to stabilize lateral motion through foot placement.
The push off mechanism works synergistically with foot placement to stabilize walking in both forward and lateral directions 59,78,79,83 .Lateral stability is partly modulated by changes in the ankle movement of the trailing leg observed during push off, for example with an increase in dorsi exion at mid-stance observed with transspinal stimulation.However, push off modulation of lateral stability is coupled to forward stability and to step length uctuations 83 .While walking without transspinal stimulation, the fractal scaling exponent of step length variability decreased (Figure 13), re ecting nearly uncorrelated uctuations 84,85 .By contrast, transspinal stimulation induced long-range correlations, suggesting that subjects adopted a control mechanism aimed at regulating step length similarly to what was observed for step width, which involves the persistence of an increase or decrease in step uctuations across multiple consecutive steps before reversing.The observed temporal correlation suggests that step length and step width uctuations were not controlled based only on step-by-step corrections, and probably re ects an increase on the modulation coupling between step width and step length with transspinal stimulation.The surrogation analysis revealed signi cant differences in step width and step length variability between the original data and their surrogate counterparts.This suggests deterministic properties in the uctuations of both step width and step length with transspinal stimulation.Transspinal stimulation increased the deterministic behavior for step length uctuations.The deterministic properties of the gait variability with transspinal stimulation are important because it promotes gait stability by allowing a successful adaptation of foot placement to the small inherent perturbations during walking.
Transspinal stimulation has been used at different frequencies and intensities to improve standing and walking abilities in neurological disorders, but systematic investigations to de ne the optimal frequency and intensity per targeted outcome are lacking.For example, transspinal stimulation at frequencies varying from 0.2 to 30 Hz promoted self-assisted standing and upright trunk posture with minimal or no external assistance 86,87 .Transspinal stimulation at 30 Hz changed the left hip-knee coordination in 2 out of 3 subjects with spinal cord injury, increased the hip range of movement by 8º in another subject 25 , and spatiotemporal gait characteristics in a motor incomplete person 88 .Frequencies from 5 to 15 Hz have been proposed to be optimal for recovery of standing because they produce bilateral leg extension, while frequencies from 25 to 50 Hz have been proposed to be optimal for recovery of walking because they evoke rhythmic alternated exion-extension movement of the legs following epidural stimulation (reviewed in 89 ).However, these conclusions are based on data from a few subjects and a comprehensive investigation is thus warranted.Correspondingly, in this study we found signi cant effects with all frequencies and intensities.Joint movements and joint coordination changed in all gait phases with 30 Hz transspinal stimulation (supplementary S1).This nding may be related to increases of the mean ring rate of muscle spindle afferents with increments in frequency, consistent to simulation models of epidural stimulation 90 .A systematic approach involving animal and human research as well as simulation models is needed to design stimulation protocols that can bene t gait asymmetry and impaired dynamic balance in upper motoneuron lesions.

Future Perspectives and Limitations
The detailed kinematic analysis of leg movement during treadmill walking suggested that transspinal stimulation at frequencies ranging from 15 to 50 Hz and intensities from sub-threshold to supra-threshold affected spatiotemporal parameters of gait, leg coordination, and dynamic stability.The well-established phase-dependent amplitude modulation of locomotor EMG activity 91,92 and adaptation of lumbosacral motoneuron activation 27,28 suggests the need for concomitant analysis of locomotor EMG activity along with kinematic and kinetic analysis of standing and walking with different transspinal stimulation protocols.Our next step is to delineate the muscle synergies during walking, and the ground reaction forces exerted during standing following transspinal stimulation at different frequencies and intensities in healthy humans and people with spinal cord injury.This approach will enable us to develop stimulation protocols that are optimal for recovery of standing and walking.Incorporation of computational models along with experimental data will further assist in de ning the stimulation protocol that suits a speci c patient and thus develop a real targeted treatment.Con ict of Interest: The authors declare that the research was conducted in the absence of any commercial or nancial relationships that could be construed as a potential con ict of interest. Figures Figure 2F.The largest effect was observed with 30 Hz sub-threshold on dorsi exion/plantar exion at terminal stance (t 53 = 3.05, p adj = 0.02, d = 0.42), pre-swing (t 53 = 5.12, p adj < 0.001, d = 0.70), and initial swing (t 53 = 3.87, p adj = 0.001, d = 0.53).The largest effect on ankle eversion/inversion was observed with 30 Hz supra-threshold at preswing (t 53 = 2.71, p adj = 0.04, d = 0.37), and at loading response (t 53 = -3.19,p adj = 0.01, d = 0.44) and with 15 Hz supra-threshold at mid-swing (t 53 = -3.77,p adj < 0.001, d = 0.52).Average right ankle joint kinematics from all subjects during treadmill walking for dorsi exion (Figure3A), inversion (Figure3B), and adduction (Figure3C) are shown for control conditions and during transspinal stimulation.The frequency and intensity that transspinal stimulation had signi cant effects on right ankle joint movement compared to control walking are shown in Figure3D.Ankle abduction increased at initial contact, eversion increased at loading response, dorsi exion increased at pre-swing, and inversion decreased at mid-swing (Figure3E).The speci c gait phases during which the ankle movement was signi cantly affected are summarized in Figure3F.The largest effect was observed with 30 Hz supra-threshold on ankle abduction/adduction at initial contact (t 53 = 3.07, p adj = 0.03, d = 0.42), eversion/inversion at loading response (t 53 = -3.64,p adj = 0.005, d = 0.50), and mid-swing (t 53 = -4.36,p adj < 0.001, d = 0.60).The largest effect on ankle dorsi exion/plantar exion at pre-swing (t 53 = 3.75, p adj = 0.004, d = 0.52) was observed with 30 Hz sub-threshold.Effects of transspinal stimulation on knee joint movementsAverage right knee joint kinematics from all subjects during treadmill walking for exion (Figure4A), adduction (Figure4B), and internal rotation (Figure4C) are shown for control conditions and during transspinal stimulation.The frequency and intensity that transspinal stimulation had signi cant effects on right knee joint movement compared to control walking are shown in Figure4D.Transspinal stimulation increased external rotation, abduction, and exion (Figure4E).The speci c gait phases during which the right knee joint movement was signi cantly affected are summarized in Figure4F.The largest effect on internal/external rotation at initial contact (t 53 = -3.64,p adj = 0.006, d = 0.50), midstance (t 53 = -3.64,p adj = 0.006, d = 0.50), mid-swing (t 53 = -3.30,p adj = 0.02, d = 0.45), terminal swing (t 53 = -3.49,p adj = 0.009, d = 0.48), and on abduction/adduction at initial swing (t 53 = -4.60,p adj < 0.001, d = 0.63) and mid-swing (t 53 = -3.32,p adj = 0.01, d = 0.46) was observed with 30 Hz sub-threshold.The largest effect on exion/extension at initial swing (t 53 = 2.98, p adj = 0.03, d = 0.41) and mid-swing (t 53 = 4.23, p adj < 0.001, d = 0.58) was observed with 30 Hz supra-threshold.Average left knee joint kinematics from all subjects during treadmill walking for exion (Figure5A), adduction (Figure5B), and internal rotation (Figure5C) are shown for control conditions and during transspinal stimulation.Transspinal stimulation increased external rotation at initial contact and at initial swing (Figures5D-F).The largest effect on internal/external rotation at initial contact was observed with 30 Hz sub-threshold (t 53 = -3.01,p adj = 0.03, d = 0.41), and at initial swing with 15 Hz supra-threshold (t 53 = -3.70,p adj = 0.005, d = 0.51).Effects of transspinal stimulation on hip joint movementAverage left hip joint kinematics from all subjects during treadmill walking for exion (Figure6A), adduction (Figure6B), and internal rotation (Figure6C) are shown for control conditions and during transspinal stimulation.Left hip exion/extension and adduction/abduction were not signi cantly affected with transspinal stimulation.Statistically signi cant differences on left hip internal/external rotation between transspinal stimulation and control walking were revealed at all frequencies and intensities (Figure6D).Internal rotation decreased from initial contact to midstance, and external rotation increased from loading response to terminal swing (Figure6E, 6F).The speci c phases of the gait cycle during which transspinal stimulation had a signi cant effect on left hip joint movements are shown in Figure 6F.The largest effect was observed with 50 Hz supra-threshold at initial contact (t 53 = -2.88,p adj < 0.05, d = 0.40), loading response (t 53 = -2.89,p adj < 0.05, d = 0.40), and terminal swing (t 53 = -3.80,p adj = 0.002, d = 0.52), while the largest effect with 30 Hz sub-threshold was observed at midstance (t 53 = -3.26,p adj = 0.017, d = 0.45), terminal stance (t 53 = -3.83,p adj = 0.003, d = 0.53), initial swing (t 53 = -3.51,p adj = 0.008, d = 0.48), and midswing (t 53 = -4.18,p adj < 0.001, d = 0.58) phases of gait.Average right hip joint kinematics from all subjects during treadmill walking for exion (Figure 7A), adduction (Figure 7B), and internal rotation (Figure 7C) are shown for control conditions and during transspinal stimulation.The frequency and intensity that transspinal stimulation had signi cant effects on right hip joint movement compared to control walking are shown in Figure 7D.Hip exion increased at loading response (t 53 = 3.49 p adj = 0.009, d = 0.48) and midstance (t 53 = 2.89, p adj < 0.05, d = 0.40) with 30 Hz supra-threshold, and hip external rotation increased at initial swing (t 53 = -3.11,p adj = 0.03, d = 0.43) and mid-swing (t 53 = -3.25,p adj = 0.016, d = 0.45) phases of gait with 30 Hz sub-threshold (Figure 7E, 7F).
Figure8Esummarizes the reversal of thighs in phase space that occurred in stance and swing phases.
The speci c phases of the gait cycle that signi cant effects of transspinal stimulation on MARP Thigh- Shank were found are shown in Figure9G.For the right MARP Thigh-Foot , the largest effect was observed with 30 Hz supra-threshold at initial contact (t 45 = -3.65,p adj = 0.006, d = 0.54), while the largest effect with 30 Hz sub-threshold was observed at terminal stance (t 53 = 3.18, p adj = 0.022, d = 0.46), pre-swing (t 53 = 4.08, p adj = 0.001, d = 0.58), initial swing (t 53 = 4.47, p adj < 0.001, d = 0.58), and mid-swing (t 51 = 4.69, p adj < 0.001, d = 0.66) phases of gait.For the left MARP Thigh-Foot , the largest effect was observed with 50 Hz sub-threshold at initial contact (t 53 = -4.25,p adj < 0.001, d = 0.58) and loading response (t 53 = -3.40,p adj = 0.011, d = 0.47).The largest effect for the left MARP Thigh-Foot was observed with 30 Hz subthreshold at pre-swing (t 53 = 3.26, p adj = 0.017, d = 0.45) and initial swing (t 53 = 3.76, p adj = 0.003, d = 0.52), and at mid-swing with 15 Hz supra-threshold (t 53 = 3.61, p adj = 0.006, d = 0.50).At loading response, the differences on left MARP Thigh-Foot between stimulation intensities (sub-threshold vs. suprathreshold) were more pronounced when transspinal stimulation was delivered at 50 Hz compared to 30 Hz (t 52 = -2.83,p adj = 0.04, d = 0.40).A representative stride of a subject during walking is shown as stick diagrams in Figure9H.The red parts of the stick diagrams correspond to the phase when the thigh moves faster than the shank in phase space, and the blue parts when the shank moves faster than the thigh.
The speci c phases of the gait cycle whereas the signi cant effects of transspinal stimulation on MARP Shank-Foot were found are shown in Figure10G.The largest effects were observed with 30 Hz suprathreshold and/or sub-threshold (p adj < 0.01, 0.46 ≤ d ≤ 0.61) at loading response, terminal stance, preswing, and initial swing, and with 50 Hz sub-threshold at mid-swing (t 53 = -3.45,p adj < 0.01, d = 0.48) and terminal swing (t 52 = -3.48,p adj = 0.009, d = 0.48).This suggests that the coordination dynamics at the ankle joint during push off were affected most by transspinal stimulation delivered at 30 Hz.A representative stride of a subject during walking is shown as stick diagrams in Figure10H.The red part of the stick diagram corresponds to the phase where the thigh moves faster than the foot in phase space, and the blue part corresponds to the phase where the foot moves faster than the thigh.
The speci c phases of the gait cycle whereas the signi cant effects of transspinal stimulation on MARP Thigh-Shank were found are shown in Figure11G.The largest effects were observed with 50 Hz subthreshold on left MARP Thigh-Shank at initial contact (t 53 = -4.15,p adj = 0.001, d = 0.60), loading response (t 53 = -4.42,p adj < 0.001, d = 0.61), and midstance (t 53 = -2.99,p adj = 0.036, d = 0.41), and with 50 Hz supra-threshold at mid-swing (t 51 = 4.20, p adj = 0.001, d = 0.59).At initial contact, the differences on left MARP Thigh-Shank between stimulation intensities (sub-threshold vs. supra-threshold) were more pronounced when transspinal stimulation delivered at 50 Hz compared to 15 or 30 Hz (15×50 Hz: t 48 = -2.69,p adj = 0.043, d = 0.39; 30×50 Hz: t 48 = -3.24,p adj = 0.015, d = 0.47), loading response (30×50 Hz: t 52 = -3.14, p adj = 0.015, d = 0.44; 15×50 Hz: t 52 = -2.83,p adj = 0.03, d = 0.39).A representative stride of a subject during walking is shown as stick diagrams in Figure11H.The red parts of the stick gures correspond to the phase where the thigh moves faster than the shank in phase space, and the blue parts correspond to the phase where the shank moves faster than the thigh.Lyapunov exponentThe LyE of the right hip exion/extension (Figure12A) movement decreased with 15 Hz supra-threshold (t 54 = -2.94,p adj = 0.039, d = 0.40) and 50 Hz supra-threshold (t 54 = -2.80,p adj = 0.047, d = 0.38) transspinal stimulation indicating greater dynamic stability compared to control walking.Transspinal stimulation did not affect the dynamic stability of the other joint movements (Figures12B-F).Amount of gait variabilityThe amount of step length variability decreased with 30 Hz sub-threshold (t 48 = -3.43,p adj = 0.011, d = 0.49), 30 Hz supra-threshold (t 48 = -3.19,p adj = 0.014, d = 0.46), and 50 Hz supra-threshold (t 48 = -3.29,p adj = 0.014, d = 0.47) compared to control walking (Figure13A).The amount of step width variability increased at all frequencies and intensities (p adj < 0.04, 0.38 ≤ d ≤ 0.54), except with 15 Hz sub-threshold (Figure13B).Detrended uctuation analysisThe temporal structure of step length uctuations became more persistent with 30 Hz sub-threshold (t 51 = 4.24, p adj < 0.001, d = 0.59), 30 Hz supra-threshold (t 51 = 3.42, p adj = 0.007, d = 0.48), and 15 Hz suprathreshold (t 51 = 3.20, p adj = 0.017, d = 0.45) compared to control walking (Figure13C).The temporal structure of step width uctuations became less persistent with 15 Hz supra-threshold (t 49 = -2.98,p adj = 0.03, d = 0.43) compared to control walking (Figure13D).The difference on the temporal structure of step width variability between stimulation intensities (sub-threshold vs. supra-threshold) was signi cant (t 49 = 3.03, p adj = 0.033, d = 0.43).Surrogation analysisA representative step length time series is shown in Figure14A, while a random permutation surrogate obtained by randomly shu ed the time series is shown in Figure14B.The DFA analysis of step length time series and its surrogate are depicted in Figure14C.The DFA analysis of the original step length time series revealed persistent uctuations around its mean value, indicating the presence of long-range correlations (α-DFA = 0.87).On the other hand, when the step length time series is randomly permutated, the resulting surrogate shows α-DFA values close to 0.5 (speci cally α-DFA = 0.55).This indicates that the random permutation has disrupted the temporal long-range correlations observed in the original step length time series, resulting in uncorrelated noise.The distributions (i.e., null) of the 399 α-DFA values for the step width and step length surrogates for each transspinal stimulation group and subject are shown in Figures14D and 14E, respectively.The corresponding step length and step width without stimulation are indicated in Figure14Ffor each subject.For most subjects, the percentile of the α-DFA value of the original time series, denoted by a vertical line, was equal to or larger than 97.5 %, indicating a signi cant difference from uncorrelated white noise (p ≤ 0.025).This con rmed that most subjects exhibited deterministic characteristics of persistent uctuations in both step width and step length time series during walking with sub-threshold or supra-threshold transspinal stimulation delivered at the different frequencies.Margin of stabilityThe margin of stability decreased at the heel strike for the left (t 49 = 3.28, p adj = 0.016, d = 0.47) and right (t 50 = 3.73, p adj = 0.004, d = 0.53) legs with 30 Hz supra-threshold stimulation compared to control walking (Figures 15A, 15C), while no signi cant effects were observed during the toe off of the right and left legs (Figures 15B, AbbreviationCPGs: central pattern generators

Figure 1 Methods
Figure 1

Figure 2 Effects
Figure 2

Figure 3 Effects
Figure 3

Figure 4 Effect
Figure 4

Figure 5 Effect
Figure 5

Figure 7 Effect
Figure 7

Figure 9 Effects
Figure 9

Figure 10 Effect
Figure 10

Figure 11 Effect
Figure 11

Figure 12 Effect
Figure 12

Figure 13 Effect
Figure 13

Figure 15 Effect
Figure 15