Individual Structural Covariance Network Predicts Long-Term Motor Improvement in Parkinson Disease with Subthalamic Nucleus Deep Brain Stimulation

Study Design and Participants

The study retrospectively included patients from 2016 to 2022 at Beijing Tiantan Hospital, Capital Medical University, China. Inclusion criteria were: 1) diagnosed with idiopathic PD,²⁶ 2) underwent bilateral STN-DBS, 3) completed motor and nonmotor assessments at baseline and at least 1 year postsurgery, 4) met quality control standards on neuroimaging examination. Age-matched health controls (HCs) were individuals aged between 40–80 years without any neurologic disorders. Ultimately, a total of 138 patients with PD were retrospectively included, with 40 HCs. Baseline information is summarized in Table 1.

Clinical Examinations

Patients underwent detailed clinical assessments included Movement Disorder Society–sponsored Unified Parkinson's Disease Rating Scale (MDS-UPDRS), especially Part III of the MDS-UPDRS (MDS-UPDRS-III), Berg Balance Scale (Berg), Freezing of Gait Questionnaire (FOG-Q), Montreal Cognitive Assessment (MoCA), Hamilton Anxiety Rating Scale (HAMA) and Hamilton Depression Rating Scale (HAMD). These clinical assessments were conducted by 2 movement disorder experts within the first 2 weeks preceding DBS, during both the Med-OFF and on-medication (Med-ON) states. The Med-OFF state was defined as abstaining from antiparkinsonian drugs for at least 12 hours, while the Med-ON state was defined as having taken antiparkinsonian drugs within the previous hour. The levodopa equivalent daily dose for each patient's antiparkinsonian drugs is detailed in Table 1.²⁷ Follow-up assessments were conducted at 1–3 years postsurgery in both stimulation-on (Stim-ON)/Med-ON and Stim-ON/Med-OFF states. The specific surgical procedure and follow-up protocol can be found in the Online Supplemental Data.

The impact of DBS on clinical motor improvement was quantified as a percentage improvement on the MDS-UPDRS-III. MDS-UPDRS-III consists of 33 scores derived from 18 items, with each item scored on a scale of 0 to 4. A score of 0 indicates normalcy, while a score of 4 signifies severe impairment.²⁸ The improvement percentage of DBS was further used to categorize patients into a good improvement group (GIG) and a moderate improvement group (MIG). The formula for calculating the improvement percentage is as follows:

Previous studies have indicated that patients undergoing STN-DBS with an improvement rate exceeding 30% are considered to exhibit a beneficial response to stimulation.¹¹ Therefore, we defined the MIG as patients with an improvement rate of 30% or less. Additionally, the confidence interval for the overall maximum improvement in MDS-UPDRS-III ranges from 69.8 to 45.8.²⁹ Therefore, patients with an improvement rate of 70% or higher were categorized as the GIG. In the group analysis, only patients with DBS improvement rates below 30% and above 70% were included. This was done to enhance sensitivity in the analysis. Detailed information on patient grouping based on DBS improvement rates is shown in Fig 2 and the Online Supplemental Data.

All postoperative images of the patients were reviewed to confirm the absence of electrode displacement. The final follow-up stimulation parameters for both patient groups were comparable (Online Supplemental Data ).

Individual Differential Structural Covariance Network Measures

As the conventional group-level SCN analysis tends to lose individualized network information, we adopted a recently established method, NTP, as outlined by Liu et al,²⁰ to construct individual-specific SCNs. Details of T1 data acquisition and preprocessing can be found in the Online Supplemental Data.

Specifically, we followed the steps depicted in Fig 1. Initially, a reference SCN was constructed within the HC group (n = 40). This network was generated by calculating the partial Pearson correlation coefficient between gray matter volumes of each pair of brain regions, while considering total intracranial volume (TIV) as a covariate, and denoted this value as Pearson correlation coefficient for the nth patient (PCCn). Subsequently, every patient was introduced into the HC group, resulting in n + 1 subjects (n controls and 1 patient), and a new structural covariance network termed the perturbed network PCCn + 1 was constructed. Calculation of the difference between the perturbed network and the reference network, ΔPCCn = PCCn + 1 – PCCn, followed. Lastly, given that ΔPCCn exhibited a novel symmetrical distribution known as the “volcano distribution,”³⁰ we computed the Z score for ΔPCCn by using a Z-test as follows:

• Download figure
• Open in new tab
• Download powerpoint

FIG 1.

Schematic outline of the study. A, Gray matter volumes were computed by using the CAT12 toolbox, and gray matter volumes were extracted based on the AAL3 atlas for all HC and patients with PD. B, Individual structural covariance network computation process for patients. PCC indicates Pearson correlation coefficient.

The individual differential structural covariance network (IDSCN) for every patient was subsequently constructed, where the weight of each edge was determined by the Z-scores obtained from the Z-test. Moreover, we calculated the P value for each edge in the IDSCN for every patient based on the Z score. Finally, we identified edges in each patient's IDSCN that significantly deviated from the reference network, applying Bonferroni correction. The edges within the IDSCN convey how the inclusion of an additional patient altered the covariance of gray-matter volume pairs for specific brain regions compared with the reference group. Ultimately, for each patient, we constructed IDSCN consisting of 12,090 edges, based on 156 brain regions defined by the automated anatomical labelling atlas 3 (AAL3) (Online Supplemental Data). The edges in the IDSCN represent how the covariance between 2 brain regions in an individual patient deviates from the reference structural covariance network observed in HCs.

Network Analysis

We selected a range of sparsity thresholds (K = 0.14 ∼ 0.5, with an increment of 0.01) to binarize the IDSCN and computed global and nodal network topological properties. Global metrics (small-world attributes and global efficiency) and node-level metrics (degree centrality, local efficiency) were calculated by using the Graph Theoretical Network Analysis (GRETNA) software package.³¹ Definition of graph theory metrics can be found in the Online Supplemental Data.

For each sparsity threshold, network metrics were calculated, and quantitative analysis of graph theoretical metrics was performed by computing the area under the curve (AUC) within the entire sparsity range.

Patients Setting and Multilayer Perceptron Classification

We selected the top 5% of edges based on their Z-scores for subsequent analysis.²⁰ This approach was employed to ensure that such differences would not be confounded by an excessive number of nodes and edges in further intergroup analysis and predictive model construction. First, we assessed whether there existed group-level differences in network properties, which could serve as indicators of long-term prognosis. We divided the data set into training/validation 1) and testing 2) sets in a 7.5:2.5 ratio. Within a set, we ranked the edges based on their correlation with the DBS improvement and considered the absolute values. Additionally, covariates such as age, sex, TIV, and MoCA scores were taken into account. We selected the top 1% of the correlated edges from a set, resulting in 6 edges (n = 6). Subsequently, we conducted group-wise comparisons between these 6 edges in the overall groups.

Next, we employed a multilayer perceptron classifier (neural network) with default parameters from scikit-learn (https://scikit-learn.org/stable/) on the top 1% correlated edges in a set.³² This evaluation aimed to assess the model's ability to differentiate the responsiveness to STN-DBS based on combinations of network edges. We performed 5-fold cross-validation on a set to test the model's performance. Finally, the model was evaluated on testing set. This involved generating receiver operating characteristic (ROC) curves on testing set, plotting the true-positive rate against the false-positive rate. The performance of the tested model was quantified by calculating the AUC of the average ROC curve.

STN-DBS Outcome Prediction Using XGBoost

We utilized XGBoost (https://xgboost.readthedocs.io/en/stable/) to investigate whether preoperative network fingerprints can predict the long-term DBS efficacy for all patients. Initially, we divided the data set into training/validation and testing sets in a 75%:25% ratio, followed by feature selection within the training/validation set where we chose the top 1% of edges. Subsequently, data standardization was performed. Through grid search, we determined the optimal hyperparameter combination for the XGBoost regression model, including the number of base learners (n_estimators) and the weight shrinkage of each base learner (learning_rate). During the hyperparameter search, we experimented with different values to find the best combination. For n_estimators, we explored a range from 10 to 200, while for learning_rate, we explored a range from 0.01 to 0.1.

Following this, we trained the XGBoost regression model on the training data set by using the best hyperparameter combination and evaluated its performance through 5-fold cross-validation on the training/validation data set to assess hyperparameter performance. Subsequently, we assessed the model's performance by predicting long-term DBS improvement on the testing data set. To compute feature importance, we utilized the built-in feature_importances attribute of the XGBoost library, providing the relative importance of each feature in predicting the model. Evaluation metrics included the correlation coefficient between predicted and actual values, offering insights into the model's predictive accuracy. Additionally, we computed the mean squared error (MSE) as a quantifiable measure of overall prediction error.

Statistical Analyses

The comparison of changes in outcomes from baseline to long-term follow-up for all patients, which did not conform to a normal distribution, was assessed by using Wilcoxon signed-rank tests. Following grouping into GIG, MIG, and HC, an analysis of variance was employed to examine age and sex differences among the groups. Comparisons between the GIG and MIG groups underwent a normality test, and if the data followed a normal distribution, a 2-sample t test was applied. Non-normally distributed data were analyzed by using the Mann-Whitney U test. All results are presented as mean ± standard deviation, and statistical significance was defined as P < .05.

Intergroup differences in gray matter network connectivity were evaluated by using ANCOVA tests, with age, sex, TIV, and MoCA scores as covariates. Multiple comparisons were corrected by using the false discovery rate (FDR), and FDR-adjusted P-values (FDR p) below .05 were considered statistically significant. The correlation between gray matter volume network edges and long-term improvement in DBS was assessed by using partial correlation analysis with age, sex, TIV, and MoCA scores as covariates.