Classification of multiple sclerosis patients based on structural disconnection: A robust feature selection approach.

BACKGROUND AND
PURPOSE: Although structural disconnection represents the hallmark of multiple sclerosis (MS) pathophysiology, classification attempts based on structural connectivity have achieved low accuracy levels. Here, we set out to fill this gap, exploring the performance of supervised classifiers on features derived from microstructure informed tractography and selected applying a novel robust approach.
METHODS: Using microstructure informed tractography with diffusion MRI data, we created quantitative connectomes of 55 MS patients and 24 healthy controls. We then used a robust approach-based on two classical methods of feature selection- to select relevant features from three network representations (whole connectivity matrices, node strength, and local efficiency). Classification accuracy of the selected features was tested with five different classifiers, while their meaningfulness was tested via correlation with clinical scales. As a comparison, the same classifiers were run on features selected with the standard procedure in network analysis (thresholding).
RESULTS: Our procedure identified 11 features for the whole net, five for local efficiency, and seven for node strength. For all classifiers, the accuracy was in the range 64.5%-91.1%, with features extracted from the whole net reaching the maximum, and overcoming results obtained with the standard procedure in all cases. Correlations with clinical scales were identified across functional domains, from motor and cognitive abilities to fatigue and depression.
CONCLUSION: Applying a robust feature selection procedure to quantitative structural connectomes, we were able to classify MS patients with excellent accuracy, while providing information on the white matter connections and gray matter regions more affected by MS pathology.

INTRODUCTION

Multiple sclerosis (MS) is a chronic disease of the central nervous system characterized by inflammation and neurodegeneration. , , Both disease diagnosis and in vivo investigation of disease pathophysiology rely on the application of standard and advanced MRI techniques. Over the last years, following the developments witnessed in the artificial intelligence field, classification of MS patients via imaging‐based machine learning approaches has been pursued in order to aid MS diagnosis and/or offer new insights in MS pathophysiology. , , , Primary information source for individual classification in MS appears to be white matter (WM) lesions. Indeed, WM lesions represent the imaging hallmark of MS and the basis of the so‐called disconnection syndrome, caused by trans‐axonal degeneration of fibers connecting secluded gray matter (GM) regions. Being a pathophysiologically relevant and consistent MS feature, , brain structural connectivity seems an ideal candidate for classification tasks.

So far, however, the only study that has attempted a structural connectivity‐based classification of MS patients has reached a classification accuracy lower than 80%. This might be due to several reasons. Among these, the intrinsic limitation of fiber (streamlines) tracking in presence of WM lesions and microstructural damage, and the fact that the selection of relevant features to be fed to the classification model is based on empirically set thresholds. ,

Here, in order to overcome these limitations, we based our classification on structural metrics derived from quantitative connectomes, whose ability to characterize brain connectivity beyond impairment caused by focal lesions has been already exploited in MS. To further improve the classification accuracy, we applied a robust feature selection (FS) procedure before testing different machine learning classifiers. Finally, to confirm the relevance of the selected features, not only in terms of discrimination between MS patients and healthy controls but also in terms of clinical meaningfulness, we explored their relationships with patients’ clinical status.

METHODS

Study population

Fifty‐five MS patients (36 female [F], age 50.45 ± 11.28 years, disease duration 15.5 ± 11.6 years, 22 primary progressive, 20 secondary progressive, and 13 relapsing‐remitting, median Expanded Disability Status Scale [EDSS] 4 range: 0‐6.5) and 24 healthy controls (HC) (11 F, age 50.3 ± 8.5 years) were prospectively enrolled. Clinical examination, performed within 1 week from the MRI scan, included scales assessing neurological dysfunction (EDSS), motor performance (manual dexterity‐9 Hole Peg Test and ambulation‐25 Foot Walk Test), cognitive performance (attention‐Symbol Digit Modalities Test, visuospatial memory‐Brief Visuospatial Memory Test, verbal memory‐California Verbal Learning Test, verbal fluency‐Controlled Oral Word Association Test, executive functions‐card sorting test), fatigue (Modified Fatigue Impact Scale), and depression (Beck Depression Inventory). Table 1 contains a summary of clinical and demographical data. Written informed consent was obtained from all participants before the beginning of the study procedures, according to the Declaration of Helsinki. The protocol was approved by the Institutional Review Board of the Icahn School of Medicine at Mount Sinai.

TABLE 1

Demographic and clinical data

	Healthy controls	MS patients
Age (years)	50.3 ± 8.5	50.45 ± 11.28
Sex	11 female, 13 male	36 female, 19 male
Phenotype	–	22 primary progressive 20 secondary progressive 13 relapsing‐remitting
Disease duration (years)	–	15.5 ± 11.6
EDSS median (range)	–	4 (0‐6.5)
9HPT (seconds)	–	28.13 ± 7.84
25FWT (seconds)	–	8.44 ± 10.22
SDMT	–	48.6 ± 12.21
BVMT	–	17.96 ± 8.87
CVLT	–	52.65 ± 13.47
COWAT	–	38.45 ± 13.40
CARD SORTING	–	9.46 ± 4.74
MFIS	–	37.11 ± 20.68
BDI	–	9.16 ± 7.17

Note: Summary of demographic and clinical data for the subjects involved in this study. All measures are reported as mean ± standard deviation unless otherwise indicated.

Abbreviations: 25FWT, 25‐foot walk test; 9HPT, 9‐hole peg test; BDI, Beck depression inventory; BVMT, brief visuospatial memory test; COWAT, controlled oral word association test; CVLT, California verbal learning test; EDSS, Expanded Disability Status Scale; MFIS, modified fatigue impact scale; MS, multiple sclerosis; SDMT, symbol digit modalities test; std, standard deviation.

MRI acquisition

All subjects underwent MRI on a Siemens Skyra 3T scanner (Siemens, Erlangen, Germany) with a 32‐channels head coil. The MRI protocol included the following sequences: axial T2‐weighted 3‐dimensional (repetition time [TR]: 8000 ms, echo time [TE]: 95 ms, spatial resolution 0.5 × 0.5 × 3.0 mm3); sagittal T1‐weighted 3‐dimensional magnetization‐prepared rapid gradient echo (TR/TE: 3000/2.47 ms, inversion time: 1000 ms, spatial resolution 0.8 × 0.8 × 0.8 mm3; generalized autocalibrating partially parallel acquisitions with acceleration factor R = 2); twice‐refocused spin echo echo‐planar imaging sequence for diffusion MRI with b values of 1000 and 2000 seconds/mm2 and 30 directions each (repeated for both phase acquisitions), in addition to six b = 0 images (TR/TE: 4700/100 ms, flip angle 80°, spatial resolution 1.8 × 1.8 × 2 mm3).

MRI processing

All images were processed as described in Schiavi et al according to the pipeline illustrated in Figure 1. Briefly, after segmenting T2‐hyperintense and T1‐hypointense lesions, we used the relative masks to fill the T1‐weighted images with FMRIB software library (FSL) (https://fsl.fmrib.ox.ac.uk) and processed them with FreeSurfer version 6.0 (http://surfer.nmr.mgh.harvard.edu), to obtain an cortical parcellation comprising 85 regions of interest. The obtained parcellations were visually inspected by a neurologist with more than 5 years of experience in neuroimaging, but no additional manual editing was performed on them. After state‐of‐the‐art processing steps on diffusion MR images to build the whole brain probabilistic tractogram, we ran the convex optimization modeling for microstructure informed tractography (COMMIT) , to assign the intra‐axonal signal fraction to each entry of the connectivity matrix and derive a more quantitative structural connectome. As done in previous works investigating differences in structural connectivity between HC and MS patients, , , for each node of the whole brain connectomes, we extracted the node's strength (NS; which is the sum of weights of links connected to the node) and local efficiency (LE; which is the average inverse shortest path length in the neighborhood of the node, and is related to the clustering coefficient) using Brain Connectivity Toolbox. We then ran the FS approaches described below on the entire set of connections (Whole Net [WN] (85 × 85 – 85)/2 features), and on the set of NS and LE (85 features each) separately.

FIGURE 1

Main steps of processing pipeline. From diffusion‐weighted imaging (DWI), we built the whole brain probabilistic tractography and then ran the convex optimization modeling for microstructure informed tractography (COMMIT). In parallel, we parcellated the T1 images with the Desikan–Killiany atlas. Combining the cortical parcellation with COMMIT weights, we built the quantitative connectomes used for the feature selection procedure. MS, multiple sclerosis; NS, node strength

Feature selection

The goal of the FS step was to identify the set of most important features within each network representation (WN, NS, LE) for the discrimination between healthy subjects and MS patients. To achieve this, we devised a robust procedure, based on a leave one out (LOO) cross‐validation scheme already applied in other scenarios when estimating feature stability. Briefly, given a population of N subjects, we applied N times a FS method, each time on a set in which a single element has been removed (ie, in the first round we applied the FS method to all subjects except the first, in the second round to all subjects except the second, and so on). Good features are those that are selected as important for discrimination in most of the rounds (ie, features which are extracted even if we have a small perturbation of the dataset). Given the nonidentical sets of features extracted in the different perturbations of the original dataset, only features that were selected in more than 50% of the rounds were considered as discriminative within our robust framework. To further increase the robustness of the results, we applied two different FS schemes:

Least absolute shrinkage and selection operator (LASSO): an effective and widely applied method for FS, based on a penalized regression model. To estimate the parameter of the model in each LOO run and increase the robustness against random fluctuations, we repeated 10 times a 10‐fold cross‐validation, keeping the value that minimized the mean square error on the left‐out fold (averaged over the 10 runs).

Infinite feature selection (IFS): a very recent graph‐based feature ranking approach that exploits properties of power series matrices and Markov chains derived from a graph representation of the FS problem. This method has shown to be very competitive and effective in many problems: for example, IFS has been successfully employed for cancer classification , that, similarly to our case, typically involves few subjects and many features: diagnosis of amyotrophic lateral sclerosis and exploring gender differences in cortical morphological networks, just to cite a few.

Only features selected by both methods were considered as robust for our results.

Classification

The classification analysis aimed at quantitatively measuring the discrimination capability of a given set of features. To compute the classification accuracy, we applied the LOO cross‐validation protocol to five different classifiers, ranging from the simple nearest neighbor technique up to the more complex support vector machines and random forests (nearest neighbor and Kth nearest neighbor with Euclidian distance, support vector machine with linear kernel, support vector machine with radial basis function kernel, and random forest with 100 trees). The parameters of the classifiers have been automatically set via cross‐validation on the training set.

In order to have a comparison with standard analyses, we also computed the accuracies of the different classifiers when applying the “Standard FS” procedure, which represents standard preprocessing done in this field. This scheme consists of postprocessing usually applied to the whole connectome (ie, on the WN representation), and it retains only those features (ie, connections) that are present in at least 50% of the subjects involved in the study. , Even if these results are not directly comparable with those obtained with the proposed approach (in the previous case, the LE/NS have been computed before applying the FS), with this comparison we can get an idea on how the proposed approach compares with established preprocessing pipelines (which implies filtering with Standard FS and computing quantities such LE/NS from the filtered net). To aid in this comparison, we also computed the No Information Rate (or Null Information Rate [NIR]). The NIR is the accuracy obtained with a “silly” classifier that assigns every object to the class that is most probable a priori, that is, to the class that is the most frequent in the training set. The NIR is typically used to have an idea of the minimum significant accuracy below which a classifier is considered useless.

Relationships between extracted features and clinical scores

To confirm the meaningfulness of the selected features, we tested their correlation with clinical scales. Bivariate correlations were run between clinical scores and selected features for each network representation (WN, NS, LE). Significance level was set at p < .002, that is, Bonferroni corrected for multiple comparisons (0.05/23 as the number of features tested).

RESULTS

Extracted features

We applied the proposed FS pipeline to the three representations (WN, NS, LE). Obtained features (ie, those features that were selected by both LASSO and IFS in at least 50% of the runs) are listed in Table 2 and graphically represented in Figure 2. For each representation and for each FS scheme, we reported the feature together with its selection percentage (ie, how many times such feature has been selected among the different LOO runs and repetitions—eg, 1.00 indicates that the feature has always been selected).

TABLE 2

Features selected for each representation (whole network, local efficiency, and node strength)

Whole net
Features	IFS	LASSO
L.superiorparietal → R.supramarginal	1.000	1.000
R.thalamus → R.paracentral	0.987	0.934
L.middletemporal → R.superiorparietal	1.000	0.914
L.supramarginal → R.postcentral	0.987	0.909
L.paracentral → R.caudalmiddlefrontal	1.000	0.757
R.caudate → R.lateraloccipital	0.975	0.701
L.lateraloccipital → L.thalamus	1.000	0.678
R.bankssts → R.precuneus	0.937	0.653
L.superiorparietal → L.caudate	1.000	0.605
R.thalamus → R.precuneus	0.987	0.590
L.rostralmiddlefrontal → R.superiorfrontal	1.000	0.524

Note: For each feature, we report how many times it has been selected among the different leave one out runs and repetitions for the two different feature‐selection schemes: infinite feature selection (IFS) and least absolute shrinkage and selection operator (LASSO). In the features name, L indicates left and R right.

FIGURE 2

Selected features for whole net (left), local efficiency (top right), and node strength (bottom right). For whole net, we graphically show the selected connection and their mean strengths in healthy controls and multiple sclerosis (MS) patients. For local efficiency and node strength, we show the selected nodes color‐coded as follows: right superior frontal gyrus in blue; right cuneus in green; right pericalcarine in yellow; right thalamus in red; right caudate in magenta; left superior parietal in orange; left parsopercularis in light blue; left thalamus in yellow

Classification results

In Table 3, we report the accuracy results for each classifier on WN, NS, and LE before and after applying the proposed FS methodology. In the former case, each subject is represented with all features derived from the quantitative connectomes; in the latter, each subject is represented only with the features surviving FS. Average and maximum accuracy computed across all classifiers for a given representation are also shown.

TABLE 3

Classification accuracies before and after the application of the proposed feature selection procedure

Representation	Max	Average	NN	KNN	SVM‐LIN	SVM‐RBF	RF
WN	0.8101	0.7443	0.6456	0.6709	0.7975	0.8101	0.7975
WN + FS	0.9114	0.8658	0.8481	0.8354	0.8608	0.9114	0.8734
LE	0.8354	0.7949	0.7468	0.7848	0.8354	0.7848	0.8228
LE + FS	0.8481	0.8329	0.8228	0.8481	0.8228	0.8481	0.8228
NS	0.7975	0.7139	0.6582	0.7342	0.6962	0.6835	0.7975
NS + FS	0.8608	0.8278	0.8354	0.7722	0.8608	0.8481	0.8228

Note: For each representation: whole net (WN), local efficiency (LE), and node strength (NS), we report the accuracy results before and after the application of the feature selection (FS) procedure for each classifier (NN, KNN, SVM‐LIN, SVM‐RBF, and RF) as well as the maximum (Max) and the average scores among them (Max and Average).

Abbreviations: KNN, Kth nearest neighbor; NN, nearest neighbor; RF, random forest; SVM‐LIN, Support Vector Machine linear kernel; SVM‐RBF, Support Vector Machine radial basis function kernel.

In Table 4, we report the accuracies of the different classifiers applied to WN processed with Standard FS as well as NS/LE representations derived from the network preprocessed with Standard FS.

TABLE 4

Classification accuracies obtained with the standard feature selection approach

Representation	Max	Average	NN	KNN	SVM‐LIN	SVM‐RBF	RF
WN + StandardFS	0.7722	0.7038	0.6582	0.7089	0.7089	0.6709	0.7722
StandardFS + LE	0.8354	0.7989	0.7884	0.7884	0.8101	0.7722	0.8354
StandardFS + NS	0.8228	0.7519	0.6582	0.7342	0.8228	0.7468	0.7975

Note: Classification accuracies of the different classifiers (NN, KNN, SVM‐LIN, SVM‐RBF, and RF) and maximum (Max) and average obtained on whole net (WN), local efficiency (LE), and node strength (NS) preprocessed with the standard thresholding procedure (StandardFS).

Abbreviations: KNN, Kth nearest neighbor; NN, nearest neighbor; RF, random forest; SVM‐LIN, Support Vector Machine linear kernel; SVM‐RBF, Support Vector Machine radial basis function kernel.

The NIR was 0.6962.

Clinical impact of selected features

Correlations between selected features and clinical scores are shown in Table 5.

TABLE 5

Correlations with clinical scores

		Whole net											Node efficiency					Node strength
		L.SPG	R.Th	L.MT	L.SM	L.PC	R.Ca	L.LO	R.BT	L.SPG	R.Th	L.RMF
		R.SM	R.PC	R.SP	R.PoC	R.CMF	R.LO	L.Th	R.PCu	L.Ca	R.PCu	R.SF	R_Th	R_Ca	R_Cu	R_PCAL	R_SFG	L_POP	L_SPG	L_Th	R_Th	R_Ca	R_Pu	R_SFG
EDSS	Spearman rho	–.524***	–.275	–.212	–.354*	.064	–.109	–.385*	–.112	–.226	–.346*	–.305*	–.012	–.150	.391*	.357*	.018	.175	.129	–.211	–.353*	–.198	–.405*	.039
	p	<.001	.051	.135	.011	.656	.447	.005	.434	.110	.013	.030	.935	.295	.005	.010	.898	.220	.367	.137	.010	.163	.003	.786
9HPT	Pearson r	–.369**	–.352*	–.129	–.190	–.021	–.009	–.303*	–.200	–.053	–.345*	–.316*	–.387**	–.329*	.172	.142	.007	–.056	.04	–.236	–.578***	–.403**	–.471***	.054
	p	.010	.014	.381	.195	.885	.952	.036	.172	.719	.016	.029	.007	.022	.242	.335	.964	.705	.789	.107	<.001	.005	.001	.717
25FWT	Pearson r	–.071	–.182	–.096	–.072	–.057	–.087	–.163	–.123	–.186	–.233	–.246	–.265	–.271	.014	.002	–.025	–.039	–.199	–.279	–.303*	–.244	–.290*	.218
	p	.627	.210	.512	.622	.698	.554	.263	.399	.202	.107	.088	.066	.06	.926	.99	.865	.792	.170	.052	.034	.091	.043	.133
SDMT	Pearson r	.213	.239	.141	.216	.159	.047	.425**	.215	.409**	.287*	.347*	.321*	.445***	–.259	–.226	.058	.237	.128	.380**	.515***	.460***	.363**	.263
	p	.133	.091	.325	.128	.266	.745	.002	.129	.003	.041	.013	.021	.001	.066	.111	.684	.094	.372	.006	<.001	.001	.009	.062
BVMT	Pearson r	.291*	.249	.224	.397**	.067	.144	.446§	.116	.380**	.208	.199	.305*	.415**	–.044	–.339*	.151	.232	–.047	.360**	.416**	.515***	.432***	.318*
	p	.036	.075	.111	.004	.637	.307	.001	.414	.006	.139	.157	.028	.002	.756	.014	.286	.098	.742	.009	.002	<.001	.001	.021
CVLT	Pearson r	.271	.278*	–.041	.209	–.003	.054	.350*	.204	.331*	.144	.09	.144	.171	–.167	–.090	.181	.104	–.055	.259	.280*	.304*	.226	0
	p	.052	.046	.771	.137	.984	.703	.011	.148	.017	.307	.527	.309	.225	.238	.525	.199	.464	.699	.064	.044	.029	.107	.999
COWAT	Pearson r	–.055	.422**	–.054	.053	.117	–.017	.224	–.030	.157	.106	.074	.494***	.436***	–.110	–.178	.242	.116	–.067	.224	.418**	.427***	.310*	.265
	p	.695	.002	.704	.704	.403	.902	.107	.829	.260	.452	.599	<.001	.001	.431	.202	.081	.408	.633	.106	.002	.001	.024	.055
CARD	Pearson r	.007	–.156	.354*	.033	.105	–.029	–.054	–.066	.544***	.016	.046	.207	.354*	–.076	–.076	.149	–.030	.145	.241	.229	.196	.157	.079
SORTING	p	.962	.271	.010	.815	.460	.837	.705	.641	<.001	.912	.748	.141	.010	.59	.593	.293	.831	.304	.085	.103	.164	.267	.579
MFIS	Pearson r	–.268	–.118	–.269	–.115	–.188	–.079	–.348*	–.150	–.164	–.038	.103	–.158	–.217	–.125	–.041	–.033	.071	–.033	–.314*	–.288*	–.359**	–.512***	.138
	p	.057	.410	.056	.420	.186	.581	.012	.293	.250	.790	.471	.268	.126	.383	.773	.817	.623	.816	.025	.040	.010	<.001	.336
BDI	Pearson r	–.260	–.036	–.285*	–.155	–.171	–.005	–.271	–.091	–.392**	–.194	–.063	–.059	–.275	.065	.016	.099	–.081	–.137	–.285*	–.235	–.272	–.402**	.113
	p	.066	.801	.042	.277	.23	.973	.054	.525	.004	.172	.658	.682	.051	.651	.909	.488	.571	.339	.043	.097	.054	.003	.429

*Correlation is significant at the 0.05 level (two‐tailed).

**Correlation is significant at the 0.01 level (two‐tailed).

***Correlations surviving Bonferroni correction (p < .002, two‐tailed).

Abbreviation: 25FWT, 25‐foot walk test; 9HPT, 9‐hole peg test; BDI, Beck depression inventory; BT, bankssts; BVMT, brief visuospatial memory test; Ca, caudate; CMF, caudalmiddlefrontal; COWAT, controlled oral word association test; Cu, cuneus; CVLT, California verbal learning test; EDSS, Expanded Disability Status Scale; L, left; LO, lateraloccipital; MFIS, modified fatigue impact scale; MT, middletemporal; PC, paracentral; PCAL, pericalcarine cortex; PCu, precuneus; PoC, postcentral; POP, parsopercularis; Pu, putamen; R, right; RMF, rostralmiddlefrontal; SDMT, symbol digit modalities test; SF, superiorfrontal; SFG, superiorfrontal; SM, supramarginal; SPG, superiorparietal; Th, thalamus.

DISCUSSION

In this study, we introduced a robust FS procedure that allowed us to improve the discrimination of MS patients from HC based on structural connectivity while identifying clinically meaningful features. We devised a robust procedure that combines the results of two classical methods of FS and tested their classification accuracy with five different classifiers. Finally, the different classifiers were also used to investigate the accuracy obtained with the standard thresholding procedure.

The investigation of structural connectivity through diffusion tractography and the extraction of meaningful properties to distinguish patients from HC is a well‐established procedure. However, this methodology selects features in the whole network via a thresholding procedure solely based on empirical assumptions, with local and global graph‐based metrics being extracted from the thresholded network. ,

On the contrary, our approach was based on the selection of connections and nodes’ local properties that were consistently extracted by two different FS procedures (LASSO and IFS) in more than 50% of the runs. Indeed, in the majority of cases LASSO and IFS agreement in the selection of individual features reached 70%. This approach offers two main advantages: (i) it reduces the dimensionality of the space of features (in our case from 3570 connections—deriving from 85 × 85 cells of the connectivity matrix minus 85 diagonal entries divided by 2—to 11 connections and from 85 nodes’ local properties to 5 or 7 [for node strength and local efficiency, respectively]); (ii) it offers insights on which WM tracts/GM nodes are more affected by the disease and, likely, represents the substrates of the observed clinical disability.

We achieved a higher accuracy when employing the proposed FS methodology than when applying the standard FS approach, thus confirming the importance of reducing the dimensionality to improve classification. Indeed, with our approach the mean accuracy was higher than 82% for all tested classifiers, which, compared to the NIR of 69.62%, confirms the potentiality of the proposed method. Here, a critical observation is that our analysis is limited due to the small sample size. Indeed, these accuracies might become even higher by employing a dataset with a larger number of both MS and HC. Nevertheless, the newly proposed procedure outperformed classification accuracies obtained with the standard FS procedure in all the investigated cases.

Perhaps our most relevant finding is the fact that the largest accuracy increments as well as the highest overall accuracy were obtained by using the FS procedure on the whole net. This suggests that the procedure of FS proposed here might be used as an alternative of the classical network analysis procedure, which usually reduces the dimensionality of the space of features by deriving global metrics like efficiency and strength. Indeed, analyzing the outcomes of the standard FS procedure, an accuracy improvement over classification attempted with the entire set of features was obtained only for NS.

The major potential of our approach is that it can directly select connections more affected by the disease instead of combining and merging the information of all the connections ending in the same node. In particular, the best performance was obtained by selecting only 11 connections. This model reached a classification accuracy of 91%, a percentage considerably higher than the ones previously reported by classification works based on structural/functional connectivity analysis. ,

Regarding the biological meaningfulness of the extracted features, those with higher selection percentage at both LASSO and IFS were more consistently associated with clinical scores, confirming that our approach was able to identify relevant indicators of the patients’ clinical status. Among the selected features, NS of deep GM (thalamus, putamen, caudate) showed the strongest correlations with a variety of clinical disability measures (eg, manual dexterity, fatigue, cognition—measured in terms of processing speed, visual memory, and verbal fluency). This finding is in line with the known relevance of deep GM damage for clinical disability in MS, , , and once again underscores how disconnection of this relay from cortical areas and wider brain circuits, here expressed by the reduction in NS, is a fundamental contributor to the clinical expression of disability. Of note, while the highest discriminative accuracy was reached with the set of features derived from whole net connectivity, in agreement with the prominent role of WM already reported with deep brain algorithms, NS constitutes the best correlate of clinical disability. This is probably because it represents a better synthetic measure of disconnection for salient regions than the damage of single WM bundles.

Among the extracted features, the connection between the left middle temporal gyrus and the right superior parietal gyrus, identified in 91.4% of the runs, is anatomically implausible and represents a false positive. The identification of false positive connections represents an unsolved issue in the tractography field, and it is known to exert a nonnegligible impact on connectivity analyses. The selected false positive connection was consistently present in the HC, while it was never identified in the MS patients. This can be explained by the presence of lesions that change the directions that can be followed by the tracking in particular voxels by altering the aspect of the fiber orientation distribution functions. Although the presence of false positive is a limitation of our analyses, none of the techniques currently available in the field are spared from this issue. Recently, the COMMIT2 framework has been proposed with impressive results in terms of false positive removal. Once validated (not only in a few healthy subjects but also in pathological brains), COMMIT2 could be used instead of COMMIT to create quantitative structural connectomes and investigate the stability of both performances and features extracted. Beyond the false positive issues, the inclusion of a small number of relapsing patients might limit the generalizability of our findings.

Despite these limitations, applying a robust FS procedure to input data derived from quantitative tractography, we were able to classify MS patients with excellent accuracy, while providing information on the WM connections and GM regions more affected by MS pathology.