This study aimed to characterize how functional connectivity patterns differ between PwMS and healthy controls in relation to spatiotemporal gait and turning metrics using an MVPA-based approach. We hypothesized that linear and non-linear gait metrics would be associated with rs-fMRI connectivity patterns in frontal, parietal, occipital, and temporal regions, particularly within the DMN, frontoparietal network (FPN), somatomotor network (SN), and visual network (VIS) as defined by the Yeo et al. (2011) atlas17,20,21,22. Additionally, we expected that PwMS and healthy controls would exhibit distinct associations between functional connectivity and mobility performance.

Twenty-nine adults with relapsing-remitting multiple sclerosis (RRMS) and 28 healthy controls participated in the cross-sectional study, each completing two laboratory visits. Healthy controls were recruited to match the MS group in age and sex at the group level, minimizing potential confounding effects of these demographic variables on between-group comparisons. Participants were excluded if they: could not walk or stand for 10 min without an assistive device, had MRI contraindications, or had musculoskeletal or vestibular conditions. Healthy controls were also required to be free from clinically diagnosed neurological or mobility-affecting conditions. Further participant details are provided in Table 1. The study was approved by the Colorado State University Institutional Review Board (IRB#: 18-7738 H), all participants provided written informed consent, and all experiments were performed in accordance with relevant guidelines and regulations.

Walking performance was assessed using a self-paced two-minute walk test (2MWT) along a 27-meter hallway. Turning performance was evaluated through three separate 360° turn trials at a self-selected fast pace, one minute of continuous alternating 360° turns at a natural pace, and 180° turns during the 2MWT. Before each trial, foot placement was standardized using a template, and participants were instructed to stand still with a forward gaze. During the 2MWT and 180° turns, participants were asked to turn as if retrieving a forgotten item. For 360° turns, participants were encouraged to turn naturally, avoiding spinning or military-style turns. All assessments were conducted barefoot in a well-lit open space, with research staff spotting for safety.

Walking and turning metrics were collected using Opal wireless inertial sensors and quantified using validated software (Mobility Lab, V2) (Opals by Clario - APDM Wearable Technologies, Portland, OR, USA). The primary walking metrics included cadence (steps/min), double support time as a percent of total gait cycle time (%GCT), gait speed (m/s), and stride length (m). These metrics were selected based on their sensitivity to gait impairments in PwMS and their frequent use in both clinical and research settings to capture key aspects of spatiotemporal gait performance. Gait speed and stride length have shown strong discriminative power in prior MS studies, while cadence and double support time provide additional insight into rhythm and stability during walking. The primary turning measures for the 360° and 180° turns included turn angle (degrees), turn duration (s), and peak turn velocity (degrees/s), which are widely used to quantify turning dynamics and have been shown to differentiate between impaired and unimpaired mobility performance. Turn angle reflects how closely an individual matches the prescribed turn amount (e.g., a 360° turn), with greater deviations, particularly those exceeding the target, are often associated with greater mobility function, and angles closer to the target indicating a more cautious movement strategy. Turn duration captures the total time required to complete a turn, where longer durations indicate reduced turning performance. Peak turn velocity represents the highest angular velocity achieved during the turn, with higher values generally reflecting more optimal turning performance.

Resting-state image acquisition was completed with a 3T Siemens MAGNETOM Prisma (Siemens Medical Solutions, USA, Inc., Malvern, PA) MRI scanner using a 32-channel head coil and parallel imaging. The anatomical scan consisted of a high-resolution T1-weighted MPRAGE image (GRAPPA acceleration factor 2, voxel-size = 0.8 × 0.8 × 0.8mm, TR = 2400ms, TI = 1000ms, TE = 2.07ms, flip angle = 8°, FoV = 256 mm). The resting-state images were acquired with a fast echo-planar imaging sequence with BOLD contrast (TR = 460ms, TE = 27.20ms, flip angle = 44°, 56 slices, slice thickness = 3.00 mm, acceleration factor 8, number of measurements = 1,044), totaling 8 min of scan time. While in the scanner and prior to initiating the resting-state sequence, participants were instructed to remain awake, focus their gaze on a projected fixation cross, and keep a clear and relaxed mind.

Preprocessing was conducted using CONN Toolbox (version 22.a) and MATLAB (R2023a, The MathWorks Inc, Natick, MA, USA). Functional and anatomical data underwent a standard preprocessing pipeline, including realignment, slice timing correction, outlier detection, segmentation, MNI-space normalization, smoothing, and denoising. Functional data were realigned using the SPM realign and unwarp procedure, where all scans were co-registered to a reference image using a least squares approach and a 6-parameter (rigid body) transformation and resampled using b-spline interpolation to correct for motion and magnetic susceptibility interactions. Outlier scans were identified using ART if framewise displacement exceeded 0.5 mm or global BOLD signal changes exceeded three standard deviations. A total of 38 scans were flagged as outliers out of 59,470 total volumes across all participants (< 0.064%), indicating high overall data quality. A reference BOLD image was then generated by averaging non-outlier scans. Functional and anatomical data were normalized to MNI space, segmented into grey matter, white matter, and CSF, and resampled to 2 mm isotropic voxels using SPM unified segmentation and normalization algorithm with the IXI-549 tissue probability map template. Functional data were smoothed using an 8 mm full-width half maximum Gaussian kernel.

To evaluate potential group differences in motion, we compared mean head motion between groups and found no significant differences (healthy controls: 0.031 ± 0.008 mm; PwMS: 0.034 ± 0.008 mm; p = 0.26), suggesting comparable head motion and scan quality between groups.

Denoising followed a standard pipeline, regressing out potential confounds: white matter (5 CompCor components), CSF (5 CompCor components), motion parameters and derivatives (12 factors), outlier scans (< 3 factors), session effects and derivatives (2 factors), and linear trends (2 factors). A bandpass filter (0.008-0.1 Hz) was applied to the BOLD timeseries. CompCor components were estimated by averaging BOLD signals within eroded white matter and CSF masks, extracting principal components orthogonal to motion parameters and outlier scans.

MVPA was used to identify rs-fMRI patterns associated with mobility measures and group differences. This data-driven approach applies singular value decomposition to seed-based correlations, producing a low-dimensional representation of voxel-to-voxel connectivity. To reduce dimensionality, a two-step principal component analysis (PCA) was performed: (1) 64 subject-specific components were retained, and (2) five eigenpatterns were selected, maintaining a 10:1 subject-to-component ratio. Independent F-tests were conducted on the five eigenpatterns to identify significant connectivity variations between groups for the four mobility conditions (normal and fast 360° turns, normal 180° turns, and straight-ahead walking). Separate models were developed for each condition, controlling for age. Clusters meeting k ≥ 50 voxels and FDR-corrected p < 0.05 were retained as seeds for post-hoc seed-to-voxel analyses. For the seed-to-voxel analyses, each mobility condition was analyzed separately, with associated variables included as predictors (gait: cadence, double support time, gait speed, stride length; turning: turn angle, turn duration, peak velocity). Age was retained as a covariate, with voxel-level p < 0.001 and FDR-corrected cluster-level p < 0.05. The associated network membership for significant clusters was based on the 7-network Yeo et al. (2011) atlas.

Statistical analyses were conducted using JMP Pro 17 and the CONN Toolbox, with an alpha level of 0.05. Group differences in categorical demographic variables were evaluated using chi-square tests, and independent t-tests were used for continuous variables.

For all mobility outcomes, separate mixed effects models using REML estimation were conducted for each dependent variable across the three conditions: straight walking, 180° turning, and 360° in-place turning. For all straight walking and 360° in-place turning metrics, models included age and sex as fixed-effect covariates and participant ID as a random effect. For 180° turn metrics, the total number of turns completed during the 2MWT was included as a covariate due to significant group differences in turn count (Wilcoxon Two-Sample Test, p = 0.027; mean ± SD: healthy controls = 4.18 ± 0.67, PwMS = 3.65 ± 0.95). Turn metrics were calculated as the average of all completed turns per participant. To evaluate model assumptions, linearity, homoscedasticity, and normality of residuals were evaluated using residual plots and Q-Q plots.

To assess associations between resting-state functional connectivity (rs-FC) and mobility performance, z-scored connectivity values from significant MVPA-derived seed clusters were correlated with mobility metrics from post-hoc seed-to-voxel results. Correlations were conducted separately for PwMS and healthy controls. Pearson correlations were used for normally distributed data, and Spearman rank correlations were applied for non-normally distributed variables (assessed via Shapiro-Wilk tests and Q-Q plots). Correlation strength was classified for both types as follows: very strong (± 0.90-1.00), strong (± 0.70-0.89), moderate (± 0.50-0.69), weak (± 0.30-0.49), and negligible ( ± < 0.30). When appropriate, FDR correction was performed using the Benjamini-Hochberg procedure, with an adjusted p-value threshold of p < 0.05 for mixed-effect model results and seed-based associations.