The experimental protocol was approved by the local medical research ethics committee, and written informed consent was obtained from each participant. The study followed a longitudinal design (four time points: ≤7 days, 1 month, 3 months, and >6 months), and the dataset (44 patients and 10 healthy controls) was used to explore the differences in the evolution patterns of specific CSTs and M1–M1 rsFC changes in patients with subcortical stroke.

The inclusion criteria for patients with stroke were as follows: (a) first-onset acute ischemic stroke; (b) a single lesion in the basal ganglia and neighboring regions; and (c) right-handedness before stroke onset (Oldfield, 1971). The exclusion criteria were as follows: (a) recurrent stroke defined by clinical history and magnetic resonance imaging (MRI) evaluation; (b) any other brain abnormalities on MR images; (c) modified Fazekas scale for white matter hyperintensities greater than 1 (Fazekas et al., 1987); and (d) history of any other neurological and psychiatric disorders. During the collection of the dataset, 58 patients with subcortical stroke were initially recruited for this study. Fourteen patients were excluded due to loss of follow-up after inclusion (n = 9), recurrent stroke (n = 3), and other brain abnormalities (n = 2). Finally, 44 patients with data on the four time points were included in this study. All patients with subcortical ischemic stroke were recruited from the Tianjin Medical University General Hospital (n = 7) and the First Affiliated Hospital of Zhengzhou University (n = 37).

Multimodal MRI data were obtained using two 3.0-Tesla Discovery MR750 MR scanners (General Electric, Milwaukee, WI) from two hospitals. Tight but comfortable foam padding was used to minimize head movement, and earplugs were used to reduce the scanner noise. Diffusion-weighted imaging (DWI), sagittal 3D T1- and T2-weighted images, and T2 fluid-attenuated inversion recovery images were acquired to identify the stroke lesion, recurrent stroke, white matter hyperintensity, and other brain abnormalities. DTI was used to extract the FA values of affected CSTs with different origins. The rs-fMRI was used to calculate the M1–M1 rsFC.

Diffusion tensor imaging data were acquired using a spin-echo single-shot echo-planar imaging (EPI) sequence. Diffusion-sensitized gradients were applied along 64 non-collinear directions with a b-value of 1000 s/mm². In addition, three sets of b = 0 images were obtained. Using an integrated parallel acquisition technique (iPAT) with an acceleration factor of 2 allowed us to obtain images with less distortion from susceptibility artifacts. We collected 50 slices from each participant. The scan parameters were repetition time (TR)/echo time (TE) = 11000/77.6 ms; field of view (FOV) = 256 × 256 mm; matrix = 128 × 128; flip angle (FA) = 90°; and slice thickness = 3 mm without gap. Resting-state fMRI data were obtained using a gradient-echo single-shot EPI sequence with the following imaging parameters: TR/TE = 2000/30 ms; FOV = 240 × 240 mm; matrix = 64 × 64; FA = 90°; slice thickness = 3 mm; gap = 1 mm; 38 interleaved transversal slices; and 180 volumes. During the resting-state fMRI scans, all subjects were instructed to keep their eyes closed, stay as still as possible, think of nothing, and not fall asleep. Sagittal 3D T1WI was acquired by brain volume with the following imaging parameters: TR/TE = 8.2/3.2 ms; inversion time = 450 ms; FA = 11°; FOV = 256 mm × 256 mm; matrix = 256 × 256; slice thickness = 1 mm, no gap; and 188 slices. DWI parameters: TR/TE, matrix, FOV, slices, and slice thickness were 3,000 ms/61 ms, 160 × 160, 240 × 240 mm, 20, and 6 mm, respectively, with b-value = 1,000 s/mm². T2-FLAIR parameters: TR/TE, matrix, slices, and slice thickness were 8,500 ms/158 ms, 256 × 256, 20, and 5 mm, respectively. T2WI was acquired from clinically used sequences.

For longitudinal analysis of patients with stroke (four time points: ≤7 days, 1 month, 3 months, and >6 months), global functional deficits were assessed using the National Institutes of Health Stroke Scale (NIHSS), and motor outcomes were evaluated by the Fugl–Meyer Assessment of the whole extremity (WE_FM) at four time points.

The DTI data were preprocessed using the FMRIB's Diffusion Toolbox (FSL 5.0; http://www.fmrib.ox.ac.uk/fsl). All diffusion-weighted images were visually inspected by two radiologists for apparent artifacts due to subject motion and instrument malfunction. For each subject, the diffusion-weighted images were registered to the corresponding b = 0 images with an affine transformation to correct for eddy-current distortion and motion displacement. Then, skulls in the images were removed using the brain extract toolbox. The diffusion tensor was reconstructed using the linear least-square fitting algorithm, which was used for calculating the diffusion indices. Then, we coregistered the individual diffusion indices into the Montreal Neurological Institute (MNI) space using a two-step method. First, we coregistered the brain-extracted b = 0 images of each subject with his/her T1-weighted images using an affine method (12 parameters), and then, the T1-weighted images were affinely coregistered into the T1 template of the MNI space. Finally, the diffusion indices were written into the MNI space using the affine parameters generated from the above steps and were resliced to 2 × 2 × 2 mm³.

The integrity of each CST fiber was assessed by the FA and the percentage of impairment. Based on the fine map of the CST (Liu et al., 2020), the FA values of each CST fiber were extracted from each patient at four time points. Moreover, we only extracted FA values from the cerebral peduncle, where the CST fibers have a relatively coherent arrangement. FA values can accurately estimate the white matter integrity only in the white matter regions with coherent fiber arrangement. We used the mean FA value of all voxels within each CST to represent the FA of the CST. Acute stroke lesions and the fine map (Liu et al., 2020) of the CST fibers were used to calculate the percentage of the impairment of each CST fiber for each patient in the acute stage. First, acute stroke lesions were acquired based on the DWI data. The individual's DWI data were spatially normalized to the EPI template in the MNI space and resampled into a 1 mm³ voxel. Stroke lesions were independently outlined on the normalized DWI using the MRIcron tool (https://www.nitrc.org/projects/mricron) by three radiologists with more than 9 years of experience. The intra-class correlation coefficient for the lesion volume was 0.98, and the result from the most senior radiologist was selected as the final lesion contour. Second, for each axial slice with an overlap between the stroke lesion and a given CST fiber, the impairment percentage of the CST fiber at this slice was defined as the ratio of the area of the overlap region to the area of the CST fiber. Finally, the largest percentage in these slices was defined as the impairment percentage of the CST fiber in this patient.

The left and right M1s were separately extracted from the Human Brainnetome Atlas (Fan et al., 2016). Subsequently, we separately extracted the overlapping regions of the left and right M1 with the 50% probability fine map (Liu et al., 2020) of the CST to define the left and right seed masks for the rsFC analysis (left: MNI coordinates: −12, −21, 66; cluster size: 61 voxels; right: MNI coordinates: 12, −21, 66; cluster size: 62 voxels).

The resting-state fMRI data were preprocessed using the Statistical Parametric Mapping software (SPM12, http://www.fil.ion.ucl.ac.uk/spm). The first 10 volumes from each subject were discarded to allow the signal to reach equilibrium and to let the participants adapt to the scanning noise. The remaining 170 volumes were corrected for the acquisition time delay between slices. None of the 44 subjects had a maximum displacement of >2 mm or a maximum rotation of >2.0°. The fMRI dataset was spatially normalized to the MNI EPI template and resampled into 3 × 3 × 3 mm³ voxels. Thereafter, several nuisance variables were regressed out from the fMRI data, including the averaged signals of the ventricles, white matter, and the whole brain, and the Friston 24 regressors (including six head motion parameters, six head motion parameters one time point before, and the 12 corresponding squared items) (Friston et al., 1996). Then, a band-pass frequency filter (0.01–0.08 Hz) was applied to reduce the low-frequency drift and high-frequency noise (Greicius et al., 2003). Finally, the filtered BOLD images were spatially smoothed using an isotropic Gaussian kernel of 8 mm full width at half maximum (FWHM).

Using the defined seed masks as the regions of interest (ROI), ROI-based rsFC analysis was performed. For each individual dataset, the Pearson correlation coefficient between the mean time series of the left and right ROIs was computed and converted into the z value using Fisher's r-to-z transformation to improve the normality.

We used a linear mixed-effects model to investigate the evolution patterns of specific CST and M1–M1 rsFC changes in patients with stroke and healthy controls in the longitudinal dataset (44 patients and 10 healthy controls). The random intercept term accounts for the correlation due to repeated measurements within a single patient (Gibbons et al., 1988). All patients were assumed to have a common slope (fixed effect) where only the intercepts were allowed to vary (random effect). The model parameters were estimated by the restricted maximum likelihood method and considered significant if the P values were less than 0.05. In healthy controls, we characterized the trajectories of these specific CSTs and M1–M1 rsFC changes to establish references to identify the stroke-induced changes. In both the patient and control groups, we identified significant longitudinal specific CSTs (four fibers from M1, PMC, S1, and SMA) and M1–M1 rsFC changes by assessing the significance of the slopes (five comparisons for the patient group and five comparisons for the control group). For each specific CST and M1–M1 rsFC change, we investigated the differences in the evolution pattern by comparing the slopes between stroke and control groups (five comparisons between the two groups). In the above-mentioned analyses, we performed the comparisons 15 times. To reduce the possible false-positive findings, we used the Bonferroni method (P < 0.05/15 = 3.33 × 10⁻³) to correct for multiple comparisons.

We assessed the relationships between the FA values of the affected specific CSTs and M1–M1 rsFC at four time points. Then, we examined the correlations between the FA values and impairment percentage of each CST fiber in the acute stage and the WE_FM scores in the chronic stage. For all correlation analyses, we used partial correlations to factor out the age, sex, and FA values of the CST fibers with other origins. In the above-mentioned analyses, we performed the correlations 24 times (16 correlations between the FA of the affected CST and M1–M1 rsFC at four time points, and eight correlations between the integrity of the affected CST and long-term motor outcomes). To reduce the possible false-positive findings, we used the Bonferroni method (P < 0.05/24 = 2.08 × 10⁻³) to correct for multiple comparisons.

The clinical and demographic data of the stroke patients and controls are listed in Table 1. The longitudinal dataset (44 patients with subcortical stroke and 10 healthy controls) included rs-fMRI, DTI, and WE_FM data from the acute to chronic stages. The stroke lesions involved the internal capsule and surrounding structures, including the internal capsule, thalamus, basal ganglia, and corona radiata (Figure 1). A total of 24 out of 44 patients had infarct lesions in the right hemisphere and 20 in the left hemisphere. The motor function of the patients was partially or completely recovered with an FMA > 52/100 for the whole extremities.

In 44 subcortical stroke patients with longitudinal DTI data, we observed longitudinal changes in the CSTs from different origins. The trajectories of the specific CST changes are shown in Figure 2. The statistical significance of the longitudinal specific CST changes in the stroke and control groups and the slope differences in the specific CST changes between the two groups are provided in Table 2. For the M1 fiber (Figure 2A), the stroke group showed a significant longitudinal change (decline over time; P = 2.83 × 10⁻⁷, Bonferroni's correction) and had a steeper slope (P = 4.00 × 10⁻¹⁷, Bonferroni's correction) than the control group. For the PMC fiber (Figure 2B), the stroke group demonstrated a steeper slope (P = 7.08 × 10⁻⁸, Bonferroni's correction), but the longitudinal change (decline over time; P = 5.76 × 10⁻³ > 0.05/15 = 3.33 × 10⁻³, Bonferroni's correction) did not differ between the groups. For the S1 fiber (Figure 2C), the stroke group showed a significant longitudinal change (decline over time; P = 3.32 × 10⁻⁴, Bonferroni's correction) and had a steeper slope (P = 2.43 × 10⁻¹¹, Bonferroni's correction) than the control group. For the SMA fiber (Figure 2D), the stroke group showed a longitudinal change (decline over time; P = 3.27 × 10⁻³, Bonferroni's correction) and had a steeper slope (P = 8.33 × 10⁻⁹, Bonferroni's correction) than the control group.

In 44 subcortical stroke patients with longitudinal resting-state fMRI data, we observed longitudinal M1–M1 rsFC changes. The trajectories of the M1–M1 rsFC changes are shown in Figure 2. The statistical significance of the longitudinal M1–M1 rsFC changes in the stroke and control groups and the slope differences in the M1–M1 rsFC changes between the two groups are provided in Table 2. With respect to M1–M1 rsFC (Figure 2E), the stroke group exhibited a longitudinal change (an increase over time; P = 1.55 × 10⁻⁶, Bonferroni's correction) and had a steeper slope (P = 1.11 × 10⁻¹⁶, Bonferroni's correction) than the control group.

The FA values of the M1 fiber were significantly negatively correlated with the M1–M1 rsFC in the subacute stage (3 months post-stroke) (partial correlation coefficient [pr] = −0.548, P = 3.06 × 10⁻⁴) (Figure 3A) and chronic stage (> 6 months post-stroke) (pr = −0.473, P = 2.07 × 10⁻³) (Figure 3B) after subcortical stroke. Similarly, a negative correlation was observed between the FA values of the SMA fiber and M1–M1 rsFC in the subacute stage (3 months post-stroke) (pr = −0.508, P = 9.63 × 10⁻⁴) (Figure 3C) after subcortical stroke. However, none of the FA values of the affected PMC and S1 fibers were correlated with the M1–M1 rsFC at the four time points after subcortical stroke (P > 0.05/24 = 2.08 × 10⁻³, Bonferroni's correction).

In the 44 patients with subcortical stroke, the FA values of the affected M1 fiber in the acute stage were positively correlated with the WE_FM scores in the chronic stage (pr = 0.495, P = 8.55 × 10⁻⁴, Bonferroni's correction; Figure 4A). However, none of the FA values of the affected PMC (pr = 0.362, P = 1.84 × 10⁻², Bonferroni's correction; Figure 4B), S1 (pr = 0.429, P = 4.64 × 10⁻³, Bonferroni's correction; Figure 4C), and SMA (pr = 0.399, P = 8.78 × 10⁻³, Bonferroni's correction; Figure 4D) fibers was correlated with the WE_FM scores in the chronic stage (P > 0.05/24 = 2.08 × 10⁻³, Bonferroni's correction). None of the early impairments of the affected M1, PMC, S1, and SMA fibers were correlated with the WE_FM scores in the chronic stage (P > 0.05/24 = 2.08 × 10⁻³, Bonferroni's correction).