NeuroImage: Clinical 38 (2023) 103376 Contents lists available at ScienceDirect NeuroImage: Clinical journal homepage: www.elsevier.com/locate/ynicl Deciphering multiple sclerosis disability with deep learning attention maps on clinical MRI Llucia Coll a, Deborah Pareto b, Pere Carbonell-Mirabent a, Álvaro Cobo-Calvo a, Georgina Arrambide a, Ángela Vidal-Jordana a, Manuel Comabella a, Joaquín Castilló a, Breogán Rodríguez-Acevedo a, Ana Zabalza a, Ingrid Galán a, Luciana Midaglia a, Carlos Nos a, Annalaura Salerno b, Cristina Auger b, Manel Alberich b, Jordi Río a, Jaume Sastre-Garriga a, Arnau Oliver c, Xavier Montalban a, Àlex Rovira b, Mar Tintoré a, Xavier Lladó c, Carmen Tur a, * a Multiple Sclerosis Centre of Catalonia (Cemcat), Hospital Universitari Vall d’Hebron, Universitat Autònoma de Barcelona, Barcelona, Spain b Section of Neuroradiology, Department of Radiology (IDI), Vall d’Hebron University Hospital, Spain, Universitat Autònoma de Barcelona, Barcelona, Spain c Research institute of Computer Vision and Robotics, University of Girona, Girona, Spain A R T I C L E I N F O A B S T R A C T Keywords: The application of convolutional neural networks (CNNs) to MRI data has emerged as a promising approach to Multiple sclerosis achieving unprecedented levels of accuracy when predicting the course of neurological conditions, including Structural MRI multiple sclerosis, by means of extracting image features not detectable through conventional methods. Addi­ Deep learning tionally, the study of CNN-derived attention maps, which indicate the most relevant anatomical features for Attention maps Disability CNN-based decisions, has the potential to uncover key disease mechanisms leading to disability accumulation. From a cohort of patients prospectively followed up after a first demyelinating attack, we selected those with T1-weighted and T2-FLAIR brain MRI sequences available for image analysis and a clinical assessment performed within the following six months (N = 319). Patients were divided into two groups according to expanded disability status scale (EDSS) score: ≥3.0 and < 3.0. A 3D-CNN model predicted the class using whole-brain MRI scans as input. A comparison with a logistic regression (LR) model using volumetric measurements as explan­ atory variables and a validation of the CNN model on an independent dataset with similar characteristics (N = 440) were also performed. The layer-wise relevance propagation method was used to obtain individual attention maps. The CNN model achieved a mean accuracy of 79% and proved to be superior to the equivalent LR-model (77%). Additionally, the model was successfully validated in the independent external cohort without any re- training (accuracy = 71%). Attention-map analyses revealed the predominant role of frontotemporal cortex and cerebellum for CNN decisions, suggesting that the mechanisms leading to disability accrual exceed the mere presence of brain lesions or atrophy and probably involve how damage is distributed in the central nervous system. 1. Introduction neurodegeneration, being one of the main non-traumatic causes of irreversible disability in young adults. The exact cause of MS and the Multiple sclerosis (MS) is a chronic disease of the central nervous pathological mechanisms ultimately leading to an irreversible accumu­ system characterised by inflammation, demyelination and lation of disability are still unknown. Furthermore, its disease course can Abbreviations: CIS, Clinically Isolated Syndrome; CNN, Convolutional Neural Network; DL, Deep Learning; EDSS, Expanded Disability Status Scale; FLAIR, Fluid- Attenuated Inversion Recovery; GM, Grey Matter; LR, Logistic Regression; LRP, Layer-Wise Relevance Propagation; MNI, Montreal Neurological Institute; MPRAGE, Magnetization-Prepared Rapid Gradient-Echo; MS, Multiple Sclerosis; MS PATHS, Multiple Sclerosis Partners Advancing Technology and Health Solutions; PDDS, Patient Determined Disease Steps; VHUH, Vall d’Hebron University Hospital; WM, White Matter. * Corresponding author at: Multiple Sclerosis Centre of Catalonia (Cemcat), Hospital Universitari Vall d’Hebron, Universitat Autònoma de Barcelona, Barcelona, Spain. E-mail address: ctur@cem-cat.org (C. Tur). https://doi.org/10.1016/j.nicl.2023.103376 Received 14 December 2022; Received in revised form 9 March 2023; Accepted 9 March 2023 Available online 15 March 2023 2213-1582/© 2023 The Authors. Published by Elsevier Inc. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by- nc-nd/4.0/). L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 be highly variable among individuals (Thompson et al., 2018). For this cross-sectional study, we included all patients after their first Magnetic Resonance Imaging (MRI) is an essential tool to diagnose demyelinating attack who had at least one brain MRI scan available for and predict prognosis of MS and is routinely used to assess disease ac­ image analysis and a clinical examination within six months after the tivity and treatment effectiveness through serial MRI analyses over time. scan. Clinical examinations included the assessment on the Expanded MRI biomarkers such as the number of brain lesions and their evolution Disability Status Scale (EDSS) score (Kurtzke, 1983). EDSS ranges from over time (Thompson et al., 2018; Tintore et al., 2015), and the quan­ 0 to 10. An EDSS of 3.0 is commonly recognised as the boundary be­ tification of brain volume (i.e., grey matter [GM] and white matter tween mild or no-disability vs moderate disability status (Amato et al., [WM] volume) (Calabrese & Castellaro, 2017; de Stefano et al., 2003), 2008; Tintoré et al., 2006; Tintore et al., 2015). are known to be good measurements to establish patients’ prognosis. We included 319 unique patients in the study, making a total dataset The accumulation of patient disability has been associated with pro­ of 382 scans i.e., for 33 patients we used one or more MRI scans, per­ gression on some of these biomarkers (Bonacchi et al., 2022; Sastre- formed at different time-points after the first attack. The experiments Garriga et al., 2020). were carried out with a fold cross validation strategy. Not to bias the In recent years, the combination of MRI and deep learning (DL)- results of the models during the cross validation strategy, we ensured based models, especially convolutional neural networks (CNN), has that all MRI scans from a given patient were included in the same fold, gained popularity thanks to their ability to solve complex classification preventing their use for training and testing at the same time. Each scan tasks. DL techniques are not dependent on pre-defined features, instead was matched with the first EDSS score obtained within the following six are capable of automatically extracting relevant information from raw months after the MRI acquisition. or minimally processed data. The whole in-house dataset was acquired in the same centre with five In neuroimaging, the use of CNN approaches has been useful not only different Siemens scanner models at two different magnetic fields (1.5 T to extract biomarkers from MRI images, but also to address the classi­ – 3.0 T). For all scanners, the MRI protocol included sagittal 3D T1 fication and the prediction of different diseases (Lian et al., 2021; van magnetization prepared rapid gradient-echo (MPRAGE) and transverse der Burgh et al., 2017; Venugopalan et al., 2021). However, the use of T2 fluid-attenuated inversion recovery (T2-FLAIR). The acquisition DL techniques to predict disease progression is still largely unexplored. protocols of each scanner are summarised in Table 1. Some recent studies have presented promising results for future pre­ diction of disability (Roca et al., 2020; Storelli et al., 2022; Tousignant 2.1.2. External validation (MS PATHS) dataset et al., 2019) and cross-sectional patient stratification (Cruciani et al., MS PATHS (Multiple Sclerosis Partners Advancing Technology and 2021; Marzullo et al., 2019). Health Solutions)(Mowry et al., 2020) is a learning health system in MS, New tendencies in the field have also put emphasis on producing started in 2016, comprising a collaborative network of 10 healthcare explainable DL-based techniques, with the goal of disentangling the centres, providing standardised routinely-acquired clinical and MRI reason behind DL decisions, therefore providing additional information data. From this large database, we randomly selected a subset, with that could be extremely useful to the end users enhancing data insights representation of all grades of disability, to be used as an independent (Bach et al., 2015; Simonyan et al., 2013; Springenberg et al., 2014). test set from the models trained with the in-house dataset. The resultant This task of deciphering the “black box” of DL-based models in MS has set was composed of 440 patients with 3D T1 MPRAGE and 3D T2-FLAIR recently been studied for disease diagnosis (Eitel et al., 2019; Lopatina sequences from four different Siemens 3 T scanner models from six et al., 2020) and MS phenotype signature decrypting (Cruciani et al., different sites (excluding the provider of our in-house dataset). All im­ 2021). All these studies concluded that the use of the layer-wise rele­ ages acquired using standardised image acquisition protocols, which vance propagation (LRP)(Bach et al., 2015), among other (similar) correspond to the 3D acquisition protocol of the Tim Trio in Table 1. All methods, is the most promising tool for these analyses providing indi­ patients completed the (patient-reported outcome) Patient Determined vidual heatmaps for each subject, called attention maps, reflecting the Disease Steps (PDDS) score (Rizzo et al., 2004). The PDDS score ranges voxel-specific relevance to the classification output, according to the DL between 0 and 8 and has been proven to have a strong correlation with model, in an easy and intuitive way. the EDSS score (Learmonth et al., 2013). The main goal of this work was to build a DL-model able to accu­ rately classify MS patients based on their disability level, while ana­ 2.2. Image Pre-processing lysing the CNN-derived attention maps in order to understand the reasons behind the decisions taken by the DL-model. We believe these The same image pre-processing was applied to both datasets, in- analyses may get us closer to deciphering the physiopathological house dataset and MS PATHS (external validation dataset from now mechanisms responsible for clinical progression in MS. on). All T1-weighted (T1-w) and T2-FLAIR sequences were pre- processed with (i) bias correction (Tustison et al., 2010), (ii) skull- 2. Materials and methods stripping (Isensee et al., 2019), (iii) registration (Jenkinson & Smith, 2001) to MNI152 space (1x1x1mm3), as well as co-registration of T2- 2.1. Datasets FLAIR sequences to T1-w space, and (iv) min–max voxel intensity normalisation. 2.1.1. In-house dataset The data used for our experiments is part of a larger in-house cohort 2.3. Proposed CNN method of patients followed over time after their first demyelinating attack (or clinically isolated syndrome, CIS), the Barcelona CIS cohort, from the The proposed pipeline to stratify patients based on EDSS < 3.0 and Multiple Sclerosis Centre of Catalonia (Cemcat), Vall d’Hebron Uni­ EDSS ≥ 3.0 is summarised in Fig. 1. For the external validation experi­ versity Hospital (VHUH) (Tintore et al., 2015). This is a prospective ment, patients were classified into PDDS < 3.0 (no or mild disability) cohort started in 1995 and is still in course. The selected subjects’ and PDDS ≥ 3.0 (moderate disability), since PDDS scores ≥ 3.0 indicate acquisition dates run from 2010 to 2020. The inclusion criteria was not moderate-marked disability, equivalent to EDSS ≥ 3.0 (Learmonth et al., dependent on the clinical phenotype at the scan acquisition time point. 2013). Our only requirement was that patients had to have experienced a CIS. After fully pre-processing the whole-brain T2-FLAIR and T1-w scans, Approval was received from the Vall d’Hebron Institute of Research – the volumes were cropped to the brain region obtaining the input patch Research and Ethics Committee (XMG-INT-2014-01 and PR(AG)389/ with a fixed size of 144x184x152mm3. These patches were used as input 2021) and informed consent was obtained from each patient conforming to train and test the CNN. When testing on new data, the output pre­ the cohort. diction was propagated to provide the LRP heatmap, complementing the 2 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Table 1 MRI acquisition parameters for each scanner used in the in-house database. Avanto Avanto Fit Symphony Symphony Tim Tim Trioa Field Strength (T) 1.5 1.5 1.5 1.5 3.0 MPRAGE TR (ms) 1980 2300 2700 1980 2300 TE (ms) 3.1 3.05 4.8 3.08 2.98 TI (ms) 1100 900 850 1100 900 Voxel size (mm) 1x1x1 1x1x1 1x1x1.2 1x1x1 1x1x1 Plane Sagittal Sagittal Sagittal Sagittal Sagittal T2-FLAIR TR (ms) 8500 8500 9000 8500 9000//2300 TE (ms) 92 99 114 95 87//392 TI (ms) 2439 2440 2500 2440 2500//1800 Voxel size (mm) 0.97x0.97x3 0.97x0.97x3 0.49x0.49x3 0.97x0.97x3 0.49x0.49x3//1x1x1 Plane Axial Axial Axial Axial Axial//Sagittal In dataset, n (%) 64 (17) 64 (17) 10 (3) 51 (13) 193 (50) TR: repetition time, TE: echo time, TI: inversion time. a With Tim Trio there are T2-FLAIR scans acquired with 2D and 3D. In table are specified by 2D parameters//3D parameters. The 3D T2-FLAIR acquisition pa­ rameters are the same used in the external validation dataset. Fig. 1. Pipeline followed in this work. a For training and inference, whole brain input is evaluated by the classifier model to predict the probability belonging to each disability status. We set a threshold of 0.5 for this probability. After inference, the probability obtained from the model is fed into the LRP algorithm to backpropagate through the CNN and obtain the attention map. As shown in b, the training procedure is only performed with the in-house database in a cross-validation strategy, hence obtaining different models. Afterwards, these models are used to evaluate the external database, where the final classification decision is obtained with the majority voting of the different models. P: probability, f: fold, LRP: layer-wise relevance propagation. numerical classification decision. extracted was projected in a global adaptive max pooling (GAP) layer to reduce feature dimensionality. The produced vector was fed into three 2.3.1. Network architecture successive 1x1x1 3D convolutional layers, with k = 128, 64, 2, where the The proposed DL network was inspired by the ResNet CNN archi­ first two were activated with a ReLU and the last one with a Softmax, tecture (He et al., 2016), built with three-dimensional (3D) layers. Each providing the probability to belong to one class or another. residual block is based on 3D convolutional layers that produce 3x3x3 and 1x1x1 kernel convolution layers, normalised with batch normal­ 2.3.2. Training procedure isation and activated with a leaky rectified linear unit (LeakyReLU). As The model was trained using only the in-house dataset. A 7-fold shown in Fig. 2, the architecture was composed of four residual blocks patient cross-validation strategy was used to train and test the model. with an increasing number of kernels k (16, 32, 64 and 128), followed by We sampled the folds to keep the same class distribution in each one, a 2x2x2 downscale max pooling operation. Afterwards, the feature map while following the distribution present along the dataset. In each 3 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Figure 2. Network architecture. Residual convolutional neural network architecture. Res block: residual block, conv: convolutional layer, k: kernels, GAP: global adaptive max pooling, ReLU: rectified linear unit. iteration five folds (275 scans) were used for training, one fold for support the decisions taken by the CNN. The method used for that was validation (55 scans) and the last fold for inference (55 scans). the LRP (Bach et al., 2015). This technique decodes the resulting clas­ To mitigate class imbalance and relatively small cohort data, we sification output through the network, propagating the relevance layer augmented the available data for training applying an axial flip to all by layer, obtaining a heatmap on the input space with each voxel subjects with an EDSS ≥ 3.0 (i.e., doubling the samples) and to 75% from contribution. the other class, EDSS < 3.0, considering the difference in class-size. A For this study, we extracted the individual LRP heatmaps for each random Gaussian noise (σ = 0.02) was also used as data augmentation to patient of the in-house dataset. The implementation of the method was create intensity variation on both channels (T1-w and T2-FLAIR). inspired in Pytorch-LRP (Böhle et al., 2019), and adapted to our specific We trained the model for a maximum of 200 epochs, with a fixed CNN, with β = 0, only reflecting the positive contributions that led the batch size of two, and an early stopping strategy based on the validation classification to the winning class. The decision of choosing only the loss behaviour to prevent overfitting. Every epoch the whole training set positive contributions was made based on previous literature findings was evaluated. The model was optimised with Adam (Kingma & Ba, (Böhle et al., 2019) and corroborated with preliminary analyses on our 2015) with a learning decay strategy depending on the validation per­ own data that showed that the negative contributions did not add formance and trained by minimising a weighted cross-entropy loss as additional information. cost function. The resulting LRP heatmaps were evaluated (i) individually and (ii) as a class-average prediction. 2.3.3. Inference Following the same sampling procedure described for training, the 2.6.1. Individual attention map analyses whole 3D brain was used as input through the trained model providing Individually, the attention maps were assessed visually and quali­ the output probabilities to belong to one class or another. The final tatively showing which voxels contributed the most to the classification classification was determined by the maximum of both probabilities, given by the model inference. Additionally, semi-quantitative analyses with a threshold fixed at 0.5. of the individual attention maps were carried out, multiplying such maps by a parcellation map. This allowed us to classify the relevant 2.4. Logistic regression model voxels for the CNN decision into the different anatomical areas. For this purpose, we set a threshold at the 95% percentile of positive relevance, To assess whether our CNN method was superior to a more con­ to catch up the most relevant areas in each case, although other ventional statistical approach, we also built a logistic regression (LR) thresholds have been studied (see Supplementary Material). model. For the LR, the disability class (EDSS < 3.0/EDSS ≥ 3.0) was considered as the dependent variable (output) and the different volu­ 2.6.2. Class-average attention map analyses metric measures were considered as the explanatory variables (input). Class-average attention maps were also built for visual inspection, Volumetric measures included: WM, GM, total intracranial, and lateral through averaging the individual values of the attention maps across the ventricle volumes, closely related to atrophy measures. These were subjects of each one of the prediction results: true positives (TP), false calculated from the extracted brain parcellation atlas (Henschel et al., positives (FP), false negatives (FN), and true negatives (TN). In addition, 2020) using the T1-w lesion filled (Prados et al., 2016). Brain lesion each class-average map was multiplied by the parcellation map, and a volume was also computed from the automatically extracted masks mean value per anatomical area was obtained. This allowed us to obtained using LST (Schmidt et al., 2012) and included as explanatory identify the anatomical areas with greatest relevance. variable. 2.6.3. Voxel-wise regression to explain individual attention maps 2.5. External validation We then carried out a quantitative analysis based on the LRP heat­ maps aiming at investigating to what extent the variability within each Inference on the external validation dataset was computed using voxel (across independent subjects) of the attention map could be directly the trained models on the in-house dataset. As seen in Fig. 1b, explained by the presence of lesions and atrophy, the best-known con­ the reported results were obtained using a majority voting of the seven tributors to disability in MS, at the voxel level. different models. For that, we first grouped our subjects by the predicted vs real class Additionally, the LR-models built on the in-house dataset were (TP, FP, FN, TN). We then computed voxel-wise regression models, one applied to the external validation dataset in exactly the same way. per each prediction group (i.e., TP, FN, FP, TN), where the LRP value at each voxel (considering all individual LRP heatmaps of the same group) was the dependent variable. As explanatory variables of these voxel- 2.6. Attention Maps: LRP wise regression models we included: i) voxel-wise binary indicator of lesion (obtained from individual T2-FLAIR scans), and ii) voxel-wise To investigate the decisions made by the CNN, attention maps were value describing the deformation suffered by the T1-w scan when computed in the input image space to show which were the regions that 4 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 moving it to a common space (MNI), calculated with the Jacobian Table 2 determinant (Lungu et al., 2019). The latter variable (ii) was used as Demographic, clinical and brain MRI characteristics of MS patients from the in- proxy for atrophy, being aware that the common template is based on house database. healthy controls. After this step we obtained four voxel-wise maps of R- Full cohort EDSS < 3.0 EDSS ≥ 3.0 p-value squared values, where each voxel indicated the proportion of the vari­ NPAC/SCAN = N = 215/ N = 104/ ability of the LRP value that could be explained by the presence of le­ 319/382 215 167 sions and native-to-MNI deformations. Additionally, we estimated the Female, n (%) 207 (65) 147 (64) 60 (58) 0.08 standardised beta coefficient for each one of the explanatory variables, Age at CIS, years, mean 32.3 [14–50] 32.4 32.2 0.78 [range] [14–49] [14–50] which indicated the relative importance of each one of these for the Confirmed diagnosisa, 260 (82) 160 (74) 100 (96) <0.001 prediction of the dependent variable, being the two standardised beta n (%) maps comparable (since they were in the same scale). Age at diagnose, years, 33.2 [14–59] 33.5 32.7 0.43 mean [range] [16–59] [14–55] CIS topography, n (%) 2.7. Evaluation measures of the classification models <0.001 Brain stem 83 (26) 59 (27) 24 (23) Optic nerve 98 (31) 74 (35) 24 (23) To evaluate the performance of the DL-based and LR-models, we Spinal cord 98 (31) 60 (28) 38 (37) used the following metrics: Other 40 (12) 22 (10) 18 (17) MS topography, n (%) <0.001 CIS 123 (39) 106 (49.5) 17 (16) • Balanced accuracy in correctly classifying each patient by means of SP 41 (13) 1 (0.5) 40 (39) their EDSS, RR 155 (48) 108 (50) 47 (45) Presence of OB, n (%) <0.001 TP + TN Positive 194 (61) 115 (53) 79 (75) Accuracy = (1) TP + TN + FP + FN Negative 75 (23) 64 (30) 11 (12) Unknown 50 (16) 36 (17) 14 (13) DD, years, mean (SD) 10.4 (7.0) 7.6 (6.6) 14.0 (5.6) <0.001 EDSS, median [range] 2.0 [0.0–9.0] 1.5 5.0 <0.001 [0.0–2.5] [3.0–9.0] • Sensitivity of correctly classified subjects with EDSS ≥ 3.0, Lesion load, ml, mean 27.5 (39.9) 10.4 (13.1) 49.6 (50.7) <0.001 TP (SD) Sensitivity = (2) Ventricles vol, ml, 29.7 (19.8) 21.5 (9.6) 40.2 (25.3) <0.001 TP + FN mean (SD) GM vol, ml, mean (SD) 787.9 (64.0) 813.1 755.4 <0.001 (48.3) (67.2) WM vol, ml, mean (SD) 694.8 (68.0) 707.5 678.3 <0.001 • Specificity of correctly classified subjects with EDSS < 3.0, (54.3) (79.6) Brain vol, ml, mean 1542.7 1572.0 1504.8 <0.001 TN (SD) (103.7) (83.7) (114.4) Specificity = (3) TN + FP Scanner model <0.001 Avanto, n (%) 64 (17) 19 (9) 45 (27) where TP are the number of correctly classified patients with EDSS ≥ Avanto Fit, n (%) 64 (17) 43 (20) 21 (13) Symphony, n (%) 10 (3) 7 (3) 3 (2) 3.0, TN are the number of correctly classified patients with EDSS < 3.0, Symphony Tim, n (%) 51 (13) 13 (6) 38 (23) FP are the number of patients with EDSS < 3.0 classified as they have Tim Trio, n (%) 193 (50) 133 (62) 60 (36) EDSS ≥ 3.0, and FN are the number of patients with EDSS ≥ 3.0 classified as they have EDSS < 3.0. The results are reported in terms of mean and PAC: patients; EDSS: expanded disability status scale; CIS: clinically isolated syndrome; SP: secondary progressive; RR: relapsing remitting; OB: oligoclonal standard deviation for the 7-fold cross validation computed on the in- bands; DD: disease duration; GM: grey matter; WM: white matter. house dataset and as majority voting of the 7 models for the external a MS confirmed diagnosis by McDonald 2017 criteria. Two patients were validation dataset. confirmed before their first demyelinating attack. 2.8. Descriptive statistics and other general aspects < 3.0 and 4.0 (in range 3.0–9.0) for patients with EDSS ≥ 3.0. Patients with EDSS ≥ 3.0 had a longer disease duration than patients with EDSS Descriptive statistics included comparison between the two groups < 3.0, with a similar age at CIS with a mean age of approximately 30 (EDSS < 3.0/EDSS ≥ 3.0). Chi-square or mixed-effects linear regression years old (SD = 8). Compared with patients with EDSS < 3.0, patients models accounting for repeated measures were used, as appropriate. with EDSS ≥ 3.0 had lower tissue volumes, GM, WM and total intra­ The proposed method and analyses were entirely implemented in cranial volume, and higher ventricle volume and lesion load. Python (https://www.python.org/), using the Pytorch library (Paszke The average accuracy of the in-house dataset across the seven folds et al., 2019). We ran all the experiments on a GNU/Linux machine box with the whole brain input patch was 79% (SD = 4%), with a sensitivity running Ubuntu 20.04, with 128 GB RAM. For training the model, we of 77% (SD = 5%) identifying patients with EDSS ≥ 3.0, and a specificity used a single Quadro RTX 5000 GPU (NVIDIA Corp, USA) with 16 GB of 81% (SD = 9%) identifying patients with EDSS < 3.0. VRAM memory. 3. Results 3.2. Comparison with a LR-model 3.1. Evaluation of DL-models: in-house dataset The LR-model built with the six brain volumes corresponding to WM, GM, ventricles, brainstem and cerebellum, intracranial volume and Out of the 319 patients, 104 were assigned with EDSS ≥ 3.0 and 215 lesion load, achieved an accuracy of 77% (SD = 7%), with a sensitivity of with EDSS < 3.0. This meant a volume of 215 MRI scans with EDSS < 68% (SD = 10%) when classifying patients with EDSS ≥ 3.0 and a 3.0, and 167 MRI scans with EDSS ≥ 3.0 (non-confirmed). The two specificity of 86% (SD = 6%). Therefore, the LR-model trained with MRI stratification groups showed comparable demographic, clinical and MRI pre-extracted features showed a 10% lower sensitivity than DL-based characteristics, as summarised in Table 2. The full cohort EDSS median models. When considering other input combinations with fewer fea­ was 2.0, with a median of 1.5 (in range 0.0–2.5) for patients with EDSS tures, the LR-model obtained always inferior results than the CNN 5 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 model, and when the LR-models only had a single feature, the accuracies 3.4.3. Voxel-wise regression analyses dropped to 50–55%. After carrying out the voxel-wise regressions on the in-house dataset, we obtained R-squared maps for the four prediction groups (Fig. 5). In 3.3. External validation dataset general, the R-squared values were low in all four groups, and mostly ranged between 0 and 0.2. In the correctly classified groups (TP, TN), Table 3 summarises the demographic, clinical and MRI data char­ the R-squared values were even lower (0–0.1). That means that in these acteristics of the patients included in the external validation dataset. In groups, at most the 80% of the variability of the attention map could not this external validation dataset, 220 patients were categorised to the be explained by the presence of lesions (observed in T2-FLAIR scans) or group with PDDS ≥ 3.0 and 220 with PDDS < 3.0. The full median the native-to-MNI deformation (of the T1-w scans), calculated through cohort PDDS was 2.5, with a median of 0.5 (in range 0.0–2.0) for the the Jacobian determinant. group of no or mild disability patients (PDDS < 3.0) and 5.0 (in range When focused on the maps of partial R-squared, we found that values 3.0–7.0) for the group of more disabled patients (PDDS ≥ 3.0). As for the were always greater for the presence of lesions than for the Jacobian in-house dataset, the disease duration was longer in patients with more determinant (Fig. 5). disability, with a mean age at diagnosis of 36 years old, and MRI char­ acteristics kept the same relation between groups. 4. Discussion In this validation dataset, the majority voting system after applying the seven DL models trained on the in-house cohort showed an accuracy In this work, we investigated the use of DL-models to classify MS of 71%, a sensitivity of 68% and a specificity of 75%. In contrast, when patients according to their disability status while trying to provide an we applied the different LR-models to the validation dataset, we ach­ explanation of the decisions taken by the CNN. Our findings showed that ieved accuracies close to 50% only. a DL-based model using only brain MRI scans, without any guidance, was able to stratify MS patients with high accuracy (79%) when testing 3.4. DL-based attention maps analyses data from the same cohort, and 71% when testing images from an un­ seen database (external validation dataset). Furthermore, the compari­ 3.4.1. Individual attention maps son with the LR-model brought to light the superiority of the CNN Individual attention maps showed that, in both disability groups, the model. Finally, our attention-map analyses revealed that the most most relevant voxels that led the classification decisions were mainly relevant anatomical areas that the CNN model used to decide the level of located in the periventricular WM regions, which often contained patients’ disability were the frontotemporal cortex and the cerebellum demyelinating lesions, and frontal and temporal cortical areas (see and did not depend on the mere presence of lesions or atrophy in these example Fig. 3a). Fig. 3b shows a case-example of the distribution of locations. relevant voxels across the different anatomical regions, considering all Compared with previously reported studies that used DL-based voxels with a relevance above a 95% percentile. Other thresholds were models to perform predictions of MS progression, our work presents also explored (see Supplementary Fig. 1), however insignificant changes some similarities and some major differences. Despite the differences in in the distribution were observed. the target of the classification, i.e., cross-sectional vs future predictions, all published studies (Storelli et al., 2022; Tousignant et al., 2019), like 3.4.2. Class-average attention maps ours, used the whole brain as input to the network. However, in our Class-average maps revealed that the most relevant areas were the study, we only used routinely-acquired T1-w and T2-FLAIR MRI data, frontal cortex, cerebellar cortex, periventricular WM, temporal cortex whereas other studies required the use of several MRI modalities and and lateral ventricles, as shown in Fig. 4. The ranking of these regions additional masks (Tousignant et al., 2019) to achieve similar accuracies. slightly varied across the different groups (TP, FP, TN, FN). For the Of note, we validated our CNN model on an external, unseen cohort, patients classified with EDSS ≥ 3.0 (TP, FP), the relevance of periven­ where the performance of the classification task was more than tricular and frontal WM, where most lesions are located, was particu­ acceptable, without needing any recalibration or re-training of the larly high. Instead, for no or mild-disability statuses (TN, FN), the model. This utterly necessary step of validating the CNN model relevance of the cortex, especially frontal and temporal cortical areas, (although uncommon in the literature) strengthens our findings, laying and that of the cerebellum were the highest. the foundation for future longitudinal, comprehensive models of MS disease progression. Furthermore, having at hand an excellent tool to Table 3 discriminate patients according to their disability level may be Demographic, clinical and brain MRI characteristics of MS patients from the extremely useful in contexts where the EDSS score is not obtained as part external database. of routine practice (Bove et al., 2023). Such a discriminating tool may Full cohort PDDS < 3.0 PDDS ≥ 3.0 p-value allow the fast identification of patients who already present a certain N = 440 N = 220 N = 220 degree of disability and who, therefore, are at high risk of reaching unfavourable clinical outcomes (Signori et al., 2023). To further assess Female, n (%) 310 (70) 170 (77) 140 (64) 0.007 Age at diagnosis, years, 36.8 36.1 37.6 0.24 our model robustness to the actual disability threshold used to classify mean [range] [19–69] [19–62] [19–69] patients or the distribution of disability scores in the given study pop­ DD, years, mean (SD) 11.5 (9.1) 8.5 (7.6) 14.9 (9.5) <0.001 ulation, we carried out a post-hoc analysis on the in-house cohort where PDDS, median [range] 2.5 0.5 5.0 <0.001 we only considered those patients with an EDSS ≤ 2.0 or ≥ 4.0. We [0.0–7.0] [0.0–2.0] [3.0–7.0] Lesion load, ml, mean 13.2 (24.7) 7.9 (11.7) 18.5 (32.0) <0.001 observed that the performance of the model only improved by 2%, i.e., (SD) we obtained an accuracy, sensitivity, and specificity of 81%, 79%, and Ventricles vol, ml, mean 35.9 (24.0) 28.7 (17.4) 43.1 (27.4) <0.001 83%, respectively. This slight increase in model performance may (SD) highlight the robustness of our CNN model, whose accuracy did not GM vol, ml, mean (SD) 804.6 828.4 779.7 <0.001 seem too dependent on the threshold used to classify patients into (73.5) (64.6) (74.0) WM vol, ml, mean (SD) 706.1 718.6 693.6 <0.001 disability groups or how the disability scores were distributed in our (53.4) (47.3) (56.3) cohort. Other data-related aspects that may have affected model per­ Brain vol, ml, mean (SD) 1496.4 1526.2 1466.5 <0.001 formance include the imbalance between classes in terms of scanner (92.8) (79.7) (95.7) models. That is, while 62% of patients in the group with an EDSS < 3.0 PDDS: patient determined disease steps; DD: disease duration; GM: grey matter; were scanned in the Tim Trio scanner, only a 36% of patients with an WM: white matter. EDSS ≥ 3.0 did so (Table 2). For this reason, we performed a post-hoc 6 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Fig. 3. Example of an individual attention map analysis. This MS patient was correctly classified as moderate disability with an EDSS = 6.0. a Different T2-FLAIR axial slices with their corresponding computed attention map. b Relevant-voxel accumulation by anatomical area, obtained from the product of the binarised attention map by the brain anatomical parcellation. In this case, the frontal cortex and periventricular WM were the most relevant areas leading the decision. analysis where we assessed, separately, the results obtained with the best feature representation to solve the classification, giving the possi­ Trim Trio scans and those obtained with the rest of the scanner models. bility to extract this representation. Notice that the databases used were From this post-hoc analysis we obtained an accuracy of 79% for the Tim conformed of different scanner models and magnetic fields, to which the Trio scans and an accuracy of 76% for the rest. Thus, considering that DL-models were more robust, having<10% drop in performance when the new accuracies were indeed very similar to that obtained with the testing on the external database. On the other hand, LR-models are whole in-house cohort, it might be assumed that the potential impact of based on obtained volumetric measurements that need to be previously such imbalance on model performance was not major. However, we computed using different tools which can be affected by changes on the acknowledge that any imbalances between classes are not desirable, MRI scanners and image protocols used to acquire the data. Thus, since they may introduce biases which are difficult to assess and over­ considering that volumetric measures have been frequently used as come a posteriori. For that, this may be seen as a limitation of our work. outcome measures in clinical trials (Tur et al., 2018b), our findings Future studies will evaluate the actual impact that such forms of suggest that DL-based discriminative tools may also be considered as imbalance, so common in studies which use routinely acquired data, trial outcomes, given that seem to explain clinical outcome better than have on deep learning-based models. conventional volumetric measures feeding into a LR model. Further Another important remark is that we demonstrated the superiority of research in this regard, evaluating the CNN model’s sensitivity to clin­ our CNN model when compared with a (conventional) LR-model. DL- ical change, is therefore warranted. based models are end-to-end systems which enable the learning of the Finally, in this study we attempted to understand the reasons behind 7 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Fig. 4. Class-average attention map analyses. a Class-average attention maps (TP, FP, TN, FN) binarised at 95% percentile. b Brain parcellation with the mean attention value (normalised across groups) attributed to each anatomical region. c Mean voxel accumulation by each anatomical area depending on the class- average group. the decisions made by the CNN through qualitative and quantitative using a CNN when trying to classify patients into different disability analyses of the individual and class-average attention maps derived categories. Some other studies aiming at solving similar classification from the CNN model. To the best of our knowledge, this is the first large problems have also focused on the visualisation of individual attention study focusing on giving an explanation to the performance obtained maps, but without performing further analyses on them (Storelli et al., 8 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Fig. 5. Voxel-wise regression analyses. Average map for the TP, FP, TN and FN of the T1-w scan, lesion probability map, Jacobian determinant, LRP heatmap and the R-squared (R2) map obtained from the voxel-wise regression model built with the individual attention maps as dependent variable, and the individual lesion masks and individual Jacobian determinants as explanatory ones. The partial R2 on each separate variable is also represented. The Jacobian determinant is represented in terms of expansion (values > 1.0) and compression (values < 1.0). 2022). Moreover, other studies have mainly investigated the classifica­ genuine differences between subjects, possibly related to atrophy, tion between healthy controls and patients (Böhle et al., 2019; Eitel although not only: looking at the results of the voxel-wise regression et al., 2019; Lopatina et al., 2020), whereas the task of classifying pa­ (Fig. 5), it is possible that other pathological aspects, for example tients into different disability classes, as in our study, is possibly much changes in image texture denoting underlying non-obvious demyelin­ more challenging, given the continuous nature of the disease. Thus, we ation, may be playing a role too. On the other hand, it must be went a step further trying to apply such methods on patients only, i.e., acknowledged that some of the highlighted attentions in the cortex moderately disabled and mildly or non-disabled, with a priori no clear might also be caused by non-pathological aspects such as registration pathological boundaries between them, making it very difficult their inaccuracies. Although the high accuracy of the model might lead us to identification by providing only brain MRI data to the model. With these think that the potential impact of registration inaccuracies on model analyses we found that the voxels with the highest relevance for the CNN performance was not major, this would need to be assessed in further to make the final decisions were in the frontal and temporal cortex, studies. That is, we are unable to discern if all the attentions were related followed by brainstem/cerebellum, and periventricular WM. or not to pathological aspects, which is a limitation of the current study. The importance of cortical GM for disease progression has been It should also be highlighted the involvement of brainstem/cerebellum, suggested in several studies, although the proposed underlying patho­ whose role in disability progression in MS has been repeatedly seen in logical mechanisms differ across studies. For instance, loss of cortical the literature (Palesi et al., 2015; Savini et al., 2019; Tintore et al., GM volume (Cordano et al., 2022), microstructural cortical damage, 2010). Future studies should investigate whether those relevant imaging reflected by increased total sodium concentration (Collorone et al., aspects present in the brainstem/cerebellum have the same pathological 2021), or grey matter demyelination (Gilmore et al., 2009), including translation as those of the relevant cortical areas. An interesting note is the presence of visible lesions in the cortex (Madsen et al., 2022), may all that the frontal cortex and the cerebellum have a strong anatomical contribute to the accumulation of disability. Furthermore, the loss of the connection through two distinct fronto-cerebellar pathways, which are physiological balance between the cortical thickness of different areas known to be damaged in toxic conditions involving the CNS, such as has also been proposed as a process potentially leading to clinical alcohol (Rogers et al., 2012) or heroin abuse (Wang et al., 2013). Future dysfunction in MS (Tur et al., 2018a; Tur et al., 2019). However, no studies will elucidate whether these paths are also altered in MS and, proper head-to-head comparisons across mechanisms have been carried especially, whether such alteration may be a key step in the develop­ out, which makes it difficult to understand the hierarchy or dynamics of ment of irreversible disability. such pernicious pathological events. In any case, it is plausible that our The relevance of periventricular WM may be related to the presence CNN has captured at least some of the imaging features related to them. of demyelinating lesions, which tend to appear in this location The consistency across individuals in terms of these areas likely reveal (Thompson et al., 2018), but also to lateral and third ventricle 9 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 enlargement, which has been shown to play a role too (Brown et al., vendor. Even though the data from the different MR scanner models 2017). Agreeing with that, cases with moderate disability presented a used for training were not balanced, they were so in terms of strength higher mean lesion load in periventricular WM than in mildly or non- field. This was confirmed by the excellent reproducibility of the models disabled cases. when we applied them to the external (i.e., completely unseen) valida­ Interestingly, these structures are followed by the thalamus, which tion cohort, despite the fact the disability score was not the same than also appears to be very relevant for the CNNs, in line with previous the one used for training. Another remark relates to the fact that we did studies showing the importance of deep GM structures for disease pro­ not account for potential disease-modifying treatment effects. Future gression in MS (Eshaghi et al., 2018). All these relevant areas were studies accounting for these are therefore warranted. Concerning the brought to light through both the individual attention map and the class- attention maps, they are limited by the lack of a ground truth. Surely, average map analyses. To the best of our knowledge, this is the first time other methodological choices could have been made, providing (prob­ that such areas have been identified as relevant for clinical progression ably) slightly different results, which deserve further research. More­ in MS in a completely hypothesis-free (unguided) manner, that is, over, we tried to relate the attention maps with well-known biomarkers without forcing the model to pay attention to them. of the disease, whereas future studies should investigate the association Regarding the voxel-wise regression analysis, the most striking with new, possibly more promising, biomarkers. finding was that the attention within a voxel, which indicates the importance of a given voxel to decide the disability class of a given 5. Conclusions patient, was not (only) explained by the mere presence of lesions or native-to-MNI deformation (used as proxy for volume change or atro­ In conclusion, our CNN model was able to stratify MS patients based phy) in that particular voxel. This may suggest that DL-based methods on their disability score solely using a single time-point brain MRI (T1-w pay attention to more general aspects of the image, possibly focusing on and T2-FLAIR sequences) providing an excellent performance. complex spatial relationships between voxel-wise information, consid­ Furthermore, our findings bring to light the potential of DL-based ering that distributional features of brain lesions might impact on models to provide key information about the mechanisms responsible disability progression (Tur et al., 2022), or image texture-related in­ for the accumulation of disability in MS, suggesting the relevant role of formation, maybe denoting microscopic processes such as underlying frontotemporal cortex and cerebellum for the development of irrevers­ demyelination (Gilmore et al., 2009). Of note, these more general as­ ible disability. Importantly, these findings may have immediate and pects of the image deserve further research and, anyway cannot be especially long-term implications for clinical practice, laying the foun­ summarised as presence/absence of lesions and/or atrophy in a given dations for building powerful predictive models of future disability. point. In any case, our findings strongly support the use of DL-based CRediT authorship contribution statement methods to perform classification/prediction tasks where the input data is derived from images, in line with the superiority observed by our Llucia Coll: Conceptualization, Methodology, Software, Data cura­ CNN with respect to the LR-model. This confirms the potential of DL- tion, Formal analysis, Validation, Visualization, Writing – original draft, based models to unveil key aspects of the disease which are uncatch­ Writing – review & editing. Deborah Pareto: Supervision, Conceptu­ able by the human eye. Future studies focusing on unveiling these as­ alization, Writing – original draft, Writing – review & editing. Pere pects are therefore warranted. Carbonell-Mirabent: Writing – original draft. Álvaro Cobo-Calvo: Potential implications for clinical practice of our findings may Writing – original draft. Georgina Arrambide: Writing – original draft. include, in the short term, the application of CNN-based models to Ángela Vidal-Jordana: Writing – original draft. Manuel Comabella: automatically classify patients according to their disability status using Writing – original draft. Joaquín Castilló: Writing – original draft. only routinely acquired brain MRI scans. This may be extremely useful Breogán Rodríguez-Acevedo: Writing – original draft. Ana Zabalza: in situations where large-scale therapeutic interventions, which may Writing – original draft. Ingrid Galán: Writing – original draft. Luciana vary depending on the disability status of a given patient population, Midaglia: Writing – original draft. Carlos Nos: Writing – original draft. need to be planned in a relatively rapid manner. Other, mid-term im­ Annalaura Salerno: Writing – original draft. Cristina Auger: Writing – plications include those derived from the development of CNN-based original draft. Manel Alberich: Writing – original draft. Jaume Sastre- models to predict future disability, which may be extremely powerful Garriga: Writing – review & editing. Arnau Oliver: Conceptualization, to manage patient expectations and design tailored treatment strategies. Writing – review & editing. Xavier Montalban: Writing – review & That is, we should acknowledge that this is a cross-sectional study and, editing. Àlex Rovira: Writing – review & editing. Mar Tintoré: Writing therefore, no strong statements about prediction of disease prognosis – review & editing. Xavier Lladó: Supervision, Conceptualization, based on our specific CNN model can be made. However, we believe that Writing – original draft, Writing – review & editing. Carmen Tur: Su­ our findings will help build powerful predictive models, possibly pervision, Conceptualization, Writing – original draft, Writing – review focusing on those areas identified as highly relevant through the & editing. attention map analyses. Up to now, most of those CNN-based models for future prediction that have been published so far show a relatively Declaration of Competing Interest limited accuracy, implying that more research is sorely needed (Roca et al., 2020; Storelli et al., 2022; Tousignant et al., 2019). The authors declare the following financial interests/personal re­ Among the methodological considerations and possible limitations lationships which may be considered as potential competing interests: of our study, it is worth mentioning the relatively small sample size, Ll. Coll: nothing to disclose. D. Pareto: has received a research contract considering the data needs of DL-based models. For this reason and to with Biogen Idec, and a grant from Instituto Salud Carlos III (PI18/ make the most of the data available, we did not apply any restriction to 00823). P. Carbonell-Mirabent: his yearly salary is supported by a grant the disease duration of our patients and considered each individual scan from Biogen to Fundació privada Cemcat for statistical analysis. A. which could be matched to an EDSS score as an independent piece of Cobo-Calvo: has received grant from Instituto de Salud Carlos III, Spain; information. As future work, we would like to analyse the impact of JR19/00007. G. Arrambide: has received compensation for consulting adding clinical data, such as disease duration, age at CIS or sex, on model services, participation in advisory boards or speaking honoraria from performance. Additionally, we applied data augmentation strategies Merck, Roche, and Horizon Therapeutics; and travel support for scien­ when possible. As a result of all these strategies, the final (in-house) tific meetings from Novartis, Roche, and ECTRIMS. G. Arrambide is database used for training, overall, had enough variability, proving editor for Europe of the Multiple Sclerosis Journal – Experimental, robustness and generalisation of our models to scans from the same Translational and Clinical; a member of the executive committee of the 10 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 International Women in Multiple Sclerosis (iWiMS) network, and a award. member of the European Biomarkers in MS (BioMS-eu) consortium steering committee. She is a recipient of grants PI19/01590 and PI22/ Appendix A. Supplementary data 01570, awarded by the Instituto de Salud Carlos III (ISCIII), Ministerio de Ciencia e Innovación de España. Á. Vidal Jordana: has engaged in Supplementary data to this article can be found online at https://doi. consulting and/or participated as speaker in events organized by Roche, org/10.1016/j.nicl.2023.103376. Novartis, Merck, and Sanofi. M. Comabella: has received compensation for consulting services and speaking honoraria from Bayer Schering References Pharma, Merk Serono, Biogen-Idec, Teva Pharmaceuticals, Sanofi- Aventis, and Novartis. J. Castilló: nothing to disclose. B. Rodríguez- Amato, M.P., Portaccio, E., Stromillo, M.L., Goretti, B., Zipoli, V., Siracusa, G., Battaglini, M., Giorgio, A., Bartolozzi, M.L., Guidi, L., Sorbi, S., Federico, A., de Acevedo: has received honoraria for consulting services from Wellspect. Stefano, N., 2008. Cognitive assessment and quantitative magnetic resonance A. Zabalza: nothing to disclose. I. Galán: nothing to disclose. L. Midaglia: metrics can help to identify benign multiple sclerosis. Neurology 71 (9), 632–638. nothing to disclose. C. Nos: has received funding for travel from Biogen https://doi.org/10.1212/01.wnl.0000324621.58447.00. Bach, S., Binder, A., Montavon, G., Klauschen, F., Müller, K.-R., Samek, W., Suarez, O.D., Idec and F. Hoffmann-La Roche, Ltd. and speaker honoraria from 2015. On pixel-wise explanations for non-linear classifier decisions by layer-wise Novartis. A. Salerno: nothing to disclose. C. Auger: has received relevance propagation. PLoS One 10 (7), e0130140. speaking honoraria from Novartis, Biogen and Stendhal. M. Alberich: Böhle, M., Eitel, F., Weygandt, M., Ritter, K., 2019. Layer-wise relevance propagation for explaining deep neural network decisions in MRI-based Alzheimer’s disease nothing to disclose. J. Río: has received speaking honoraria and personal classification. Front. Aging Neurosci. 10 https://doi.org/10.3389/fnagi.2019.00194. compensation for participating on Advisory Boards from Biogen-Idec, Bonacchi, R., Meani, A., Pagani, E., Marchesi, O., Filippi, M., Rocca, M.A., 2022. The role Genzyme, Merck- Serono, Mylan, Novartis, Roche, Teva, and Sanofi- of cerebellar damage in explaining disability and cognition in multiple sclerosis phenotypes: a multiparametric MRI study. J. Neurol. 269 (7), 3841–3857. https:// Aventis. J. Sastre-Garriga: serves as co-Editor for Europe on the edito­ doi.org/10.1007/s00415-022-11021-1. rial board of Multiple Sclerosis Journal and as Editor-in-Chief in Revista Bove, R., Poole, S., Cuneo, R., Gupta, S., Sabatino, J., Harms, M., Cooper, T., Rowles, W., de Neurología, receives research support from Fondo de Investigaciones Miller, N., Gomez, R., Lincoln, R., McPolin, K., Powers, K., Santaniello, A., Sanitarias (19/950) and has served as a consultant/speaker for Biogen, Renschen, A., Bevan, C.J., Gelfand, J.M., Goodin, D.S., Guo, C.-Y., Romeo, A.R., Hauser, S.L., Campbell Cree, B.A., 2023. Remote observational research for multiple Celgene/Bristol Meyers Squibb, Genzyme, Novartis and Merck. A. sclerosis: a natural experiment. Neurology(R) Neuroimmunol. Neuroinflamm. 10 Oliver: nothing to disclose. X. Montalban: has received speaking hono­ (2), e200070. raria and travel expenses for participation in scientific meetings, has Brown, J.W.L., Pardini, M., Brownlee, W.J., Fernando, K., Samson, R.S., Prados Carrasco, F., Ourselin, S., Gandini Wheeler-Kingshott, C.A.M., Miller, D.H., Chard, D. been a steering committee member of clinical trials or participated in T., 2017. An abnormal periventricular magnetization transfer ratio gradient occurs advisory boards of clinical trials in the past years with Abbvie, Actelion, early in multiple sclerosis. Brain 140 (2), 387–398. https://doi.org/10.1093/brain/ Alexion, Biogen, Bristol-Myers Squibb/Celgene, EMD Serono, Genzyme, aww296. Calabrese, M., Castellaro, M., 2017. Cortical Gray Matter MR Imaging in Multiple Hoffmann-La Roche, Immunic, Janssen Pharmaceuticals, Medday, Sclerosis. Neuroimag. Clin. N. Am. 27 (2), 301–312. https://doi.org/10.1016/J. Merck, Mylan, Nervgen, Novartis, Sandoz, Sanofi-Genzyme, Teva NIC.2016.12.009. Pharmaceutical, TG Therapeutics, Excemed, MSIF and NMSS. A. Rovira: Collorone, S., Prados, F., Kanber, B., Cawley, N.M., Tur, C., Grussu, F., Solanky, B.S., Yiannakas, M., Davagnanam, I., Gandini Wheeler-Kingshott, C.A.M., Barkhof, F., serves on scientific advisory boards for Novartis, Sanofi-Genzyme, Ciccarelli, O., Toosy, A.T., 2021. Brain microstructural and metabolic alterations Synthetic MR, Roche, Biogen, and OLEA Medical; has received speaker detected in vivo at onset of the first demyelinating event. Brain 144, 1409–1421. honoraria from Bayer, SanofiGenzyme, Merck-Serono, Teva Pharma­ https://doi.org/10.1093/brain/awab043. Cordano, C., Nourbakhsh, B., Yiu, H.H., Papinutto, N., Caverzasi, E., Abdelhak, A., ceutical Industries Ltd, Novartis, Roche, and Biogen; and is CMO and co- Oertel, F.C., Beaudry-Richard, A., Santaniello, A., Sacco, S., Bennett, D.J., Gomez, A., founder of TensorMedical. M. Tintoré: has received compensation for Sigurdson, C.J., Hauser, S.L., Magliozzi, R., Cree, B.A.C., Henry, R.G., Green, A.J., consulting services and speaking honoraria from Almirall, Bayer 2022. Differences in age-related retinal and cortical atrophy rates in multiple Schering Pharma, Biogen-Idec, Genzyme, Merck-Serono, Novartis, sclerosis. Neurology 99 (15) e1685–e1693. Cruciani, F., Brusini, L., Zucchelli, M., Retuci Pinheiro, G., Setti, F., Boscolo Galazzo, I., Roche, Sanofi-Aventis, and Teva Pharmaceuticals. MT is former co- Deriche, R., Rittner, L., Calabrese, M., Menegaz, G., 2021. Interpretable deep editor of Multiple Sclerosis Journal. X. Lladó: is currently being sup­ learning as a means for decrypting disease signature in multiple sclerosis. J. Neural ported by the ICREA Academia Program. He has also received support Eng. 18 (4), 0460a6. de Stefano, N., Matthews, P. M., Filippi, M., Agosta, F., de Luca, M., Bartolozzi, M. L., from the DPI2020-114769RBI00 project funded by the Ministerio de Guidi, L., Ghezzi, A., Montanari, E., Cifelli, A., Federico, A., & Smith, S. M. (2003). Ciencia, Innovación y Universidades. C. Tur: is currently being funded Evidence of early cortical atrophy in MS. Neurology, 60(7), 1157 LP – 1162. https:// by a Junior Leader La Caixa Fellowship (fellowship code is LCF/BQ/ doi.org/10.1212/01.WNL.0000055926.69643.03. Eitel, F., Soehler, E., Bellmann-Strobl, J., Brandt, A.U., Ruprecht, K., Giess, R.M., PI20/11760008), awarded by “la Caixa” Foundation (ID 100010434). Kuchling, J., Asseyer, S., Weygandt, M., Haynes, J.D., Scheel, M., Paul, F., Ritter, K., She has also received the 2021 Merck’s Award for the Investigation in 2019. Uncovering convolutional neural network decisions for diagnosing multiple MS, awarded by Fundación Merck Salud (Spain) and a grant awarded by sclerosis on conventional MRI using layer-wise relevance propagation. NeuroImage: Clin. 24, 102003 https://doi.org/10.1016/j.nicl.2019.102003. the Instituto de Salud Carlos III (ISCIII), Ministerio de Ciencia e Eshaghi, A., Prados, F., Brownlee, W.J., Altmann, D.R., Tur, C., Cardoso, M.J., De Innovación de España (PI21/01860). In 2015, she received an ECTRIMS Angelis, F., van de Pavert, S.H., Cawley, N., De Stefano, N., Stromillo, M.L., Post-doctoral Research Fellowship and has received funding from the Battaglini, M., Ruggieri, S., Gasperini, C., Filippi, M., Rocca, M.A., Rovira, A., Sastre- Garriga, J., Vrenken, H., Leurs, C.E., Killestein, J., Pirpamer, L., Enzinger, C., UK MS Society. She is a member of the Editorial Board of Neurology. She Ourselin, S., Wheeler-Kingshott, C.A.M.G., Chard, D., Thompson, A.J., Alexander, D. has also received honoraria from Roche and Novartis and is a steering C., Barkhof, F., Ciccarelli, O., 2018. Deep gray matter volume loss drives disability committee member of the O’HAND trial and of the Consensus group on worsening in multiple sclerosis. Ann. Neurol. 83 (2), 210–222. Follow-on DMTs. Gilmore, C.P., Donaldson, I., Bo, L., Owens, T., Lowe, J., Evangelou, N., 2009. Regional variations in the extent and pattern of grey matter demyelination in multiple sclerosis: a comparison between the cerebral cortex, cerebellar cortex, deep grey Data availability matter nuclei and the spinal cord. J. Neurol. Neurosurg. Psychiatry 80 (2), 182–187. He, K., Zhang, X., Ren, S., Sun, J., 2016. Deep residual learning for image recognition. In: Proceedings of the IEEE Computer Society Conference on Computer Vision and Data will be made available on request. Pattern Recognition, pp. 770–778. https://doi.org/10.1109/CVPR.2016.90. Henschel, L., Conjeti, S., Estrada, S., Diers, K., Fischl, B., Reuter, M., 2020. FastSurfer – A Acknowledgments fast and accurate deep learning based neuroimaging pipeline. Neuroimage 219, 117012. https://doi.org/10.1016/j.neuroimage.2020.117012. Isensee, F., Schell, M., Pflueger, I., Brugnara, G., Bonekamp, D., Neuberger, U., Wick, A., MS PATHS is funded by Biogen. This study has been possible thanks Schlemmer, H.P., Heiland, S., Wick, W., Bendszus, M., Maier-Hein, K.H., to a Junior Leader La Caixa Fellowship awarded to C. Tur (fellowship Kickingereder, P., 2019. Automated brain extraction of multisequence MRI using artificial neural networks. Hum. Brain Mapp. 40 (17), 4952–4964. https://doi.org/ code is LCF/BQ/PI20/11760008) by “la Caixa” Foundation (ID 10.1002/hbm.24750. 100010434). The salaries of C. Tur and Ll. Coll are covered by this 11 L. Coll et al. NeuroImage: Clinical 38 (2023) 103376 Jenkinson, M., Smith, S., 2001. A global optimisation method for robust affine Montalban, X., Rovira, À., 2020. MAGNIMS consensus recommendations on the use registration of brain images. Med. Image Anal. 5 (2), 143–156. https://doi.org/ of brain and spinal cord atrophy measures in clinical practice. Nat. Rev. Neurol. 16 10.1016/S1361-8415(01)00036-6. (3), 171–182. Kingma, D.P., Ba, J.L., 2015. Adam: A method for stochastic optimization. 3rd G. Savini M. Pardini G. Castellazzi A. Lascialfari D. Chard E. D’Angelo C.A.M. Gandini International Conference on Learning Representations, ICLR 2015 - Conference Wheeler-Kingshott Default Mode Network Structural Integrity and Cerebellar Track Proceedings. Connectivity Predict Information Processing Speed Deficit in Multiple Sclerosis Kurtzke, J.F., 1983. Rating neurologic impairment in multiple sclerosis. Neurology 33 Front. Cell. Neurosci. 13. (11), 1444–1452. https://doi.org/10.1212/WNL.33.11.1444. Schmidt, P., Gaser, C., Arsic, M., Buck, D., Förschler, A., Berthele, A., Hoshi, M., Ilg, R., Learmonth, Y.C., Motl, R.W., Sandroff, B.M., Pula, J.H., Cadavid, D., 2013. Validation of Schmid, V.J., Zimmer, C., Hemmer, B., Mühlau, M., 2012. An automated tool for patient determined disease steps (PDDS) scale scores in persons with multiple detection of FLAIR-hyperintense white-matter lesions in Multiple Sclerosis. sclerosis. BMC Neurol. 13, 37. https://doi.org/10.1186/1471-2377-13-37. Neuroimage 59 (4), 3774–3783. https://doi.org/10.1016/J. Lian, C., Liu, M., Wang, L., Shen, D., 2021. Multi-task weakly-supervised attention NEUROIMAGE.2011.11.032. network for dementia status estimation with structural MRI. IEEE Trans. Neural Signori, A., Lorscheider, J., Vukusic, S., Trojano, M., Iaffaldano, P., Hillert, J., Hyde, R., Networks Learn. Syst. 33 (8), 4056–4068. https://doi.org/10.1109/ Pellegrini, F., Magyari, M., Koch-Henriksen, N., Sørensen, P.S., Spelman, T., van der TNNLS.2021.3055772. Walt, A., Horakova, D., Havrdova, E., Girard, M., Eichau, S., Grand’Maison, F., Lopatina, A., Ropele, S., Sibgatulin, R., Reichenbach, J.R., Güllmar, D., 2020. Gerlach, O., Terzi, M., Ozakbas, S., Skibina, O., Van Pesch, V., Sa, M.J., Prevost, J., Investigation of deep-learning-driven identification of multiple sclerosis patients Alroughani, R., McCombe, P.A., Gouider, R., Mrabet, S., Castillo-Trivino, T., Zhu, C., based on susceptibility-weighted images using relevance analysis. Front. Neurosci. de Gans, K., Sánchez-Menoyo, J.L., Yamout, B., Khoury, S., Sormani, M.P., 14 https://doi.org/10.3389/fnins.2020.609468. Kalincik, T., Butzkueven, H., 2023. Heterogeneity on long-term disability trajectories O. Lungu P. Pantano F. Kumfor A. Gallo D. Joel Shaw M.R. Cek K. Czekóová K. Saxunová in patients with secondary progressive MS: a latent class analysis from Big MS Data M. Dufek Marě cek, R., VaníčekVaníˇVaníček, J., & Brázdil, M. Impaired Self-Other network. J Neurol Neurosurg Psychiatry 94 (1), 23–30. Distinction and Subcortical Gray-Matter Alterations Characterize Socio-Cognitive Simonyan, K., Vedaldi, A., & Zisserman, A. (2013). Deep Inside Convolutional Networks: Disturbances in Multiple Sclerosis Front. Neurol 10 2019 525 10.3389/ Visualising Image Classification Models and Saliency Maps. ICLR 2014. https://doi. fneur.2019.00525. org/10.48550/arxiv.1312.6034. Madsen, M.A.J., Wiggermann, V., Marques, M.F.M., Lundell, H., Cerri, S., Puonti, O., Springenberg, J. T., Dosovitskiy, A., Brox, T., & Riedmiller, M. (2014). Striving for Blinkenberg, M., Romme Christensen, J., Sellebjerg, F., Siebner, H.R., 2022. Linking Simplicity: The All Convolutional Net. ICLR 2015. https://doi.org/10.48550/ lesions in sensorimotor cortex to contralateral hand function in multiple sclerosis: a arxiv.1412.6806. 7 T MRI study. Brain 10 (145), 3522–3535. https://doi.org/10.1093/brain/ Storelli, L., Azzimonti, M., Gueye, M., Vizzino, C., Preziosa, P., Tedeschi, G., de awac203. Stefano, N., Pantano, P., Filippi, M., Rocca, M.A., 2022. A Deep Learning Approach Marzullo, A., Kocevar, G., Stamile, C., Durand-Dubief, F., Terracina, G., Calimeri, F., to Predicting Disease Progression in Multiple Sclerosis Using Magnetic Resonance Sappey-Marinier, D., 2019. Classification of multiple sclerosis clinical profiles via Imaging. Invest. Radiol. 57 (7), 423–432. https://doi.org/10.1097/ graph convolutional neural networks. Front. Neurosci. 13 (JUN), 1–15. https://doi. RLI.0000000000000854. org/10.3389/fnins.2019.00594. Thompson, A.J., Baranzini, S.E., Geurts, J., Hemmer, B., Ciccarelli, O., 2018. Multiple Mowry, E.M., Bermel, R.A., Williams, J.R., Benzinger, T.L.S., de Moor, C., Fisher, E., sclerosis. Lancet 391 (10130), 1622–1636. https://doi.org/10.1016/S0140-6736 Hersh, C.M., Hyland, M.H., Izbudak, I., Jones, S.E., Kieseier, B.C., Kitzler, H.H., (18)30481-1. Krupp, L., Lui, Y.W., Montalban, X., Naismith, R.T., Nicholas, J.A., Pellegrini, F., Tintoré, M., Rovira, A., Río, J., Nos, C., Grivé, E., Téllez, N., Pelayo, R., Comabella, M., Rovira, A., Schulze, M., Tackenberg, B., Tintore, M., Tivarus, M.E., Ziemssen, T., Sastre-Garriga, J., Montalban, X., 2006. Baseline MRI predicts future attacks and Rudick, R.A., 2020. Harnessing Real-World Data to Inform Decision-Making: disability in clinically isolated syndromes. Neurology 67 (6), 968–972. https://doi. Multiple Sclerosis Partners Advancing Technology and Health Solutions (MS org/10.1212/01.wnl.0000237354.10144.ec. PATHS). Front. Neurol. 11 https://doi.org/10.3389/fneur.2020.00632. Tintore, M., Rovira, A., Arrambide, G., Mitjana, R., Río, J., Auger, C., Nos, C., Edo, M.C., Palesi, F., Tournier, J.-D., Calamante, F., Muhlert, N., Castellazzi, G., Chard, D., Castilló, J., Horga, A., Perez-Miralles, F., Huerga, E., Comabella, M., Sastre- D’Angelo, E., Wheeler-Kingshott, C.A.M., 2015. Contralateral cerebello-thalamo- Garriga, J., Montalban, X., 2010. Brainstem lesions in clinically isolated syndromes. cortical pathways with prominent involvement of associative areas in humans in Neurology 75 (21), 1933–1938. https://doi.org/10.1212/WNL.0b013e3181feb26f. vivo. Brain Struct. Funct. 220 (6), 3369–3384. Tintore, M., Rovira, L., Río, J., Otero-Romero, S., Arrambide, G., Tur, C., Comabella, M., Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Nos, C., Jesú Arévalo, M., Negrotto, L., Galán, I., Vidal-Jordana, A., Castilló, J., Gimelshein, N., Antiga, L., Desmaison, A., Köpf, A., Yang, E., DeVito, Z., Raison, M., Palavra, F., Simon, E., Mitjana, R., Auger, C., Sastre-Garriga, J., Montalban, X., 2015. Tejani, A., Chilamkurthy, S., Steiner, B., Fang, L., … Chintala, S. (2019). PyTorch: An Defining high, medium and low impact prognostic factors for developing multiple Imperative Style, High-Performance Deep Learning Library. In Proceedings of the 33rd sclerosis. Brain 138 (Pt 7), 1863–1874. https://doi.org/10.1093/brain/awv105. International Conference on Neural Information Processing Systems (pp. 8026–8037). Tousignant, A., Lemaître, P., Doina, P., Arnold, D.L., Arbel, T., 2019. Prediction of Curran Associates Inc. https://doi.org/10.5555/3454287.3455008. Disease Progression in Multiple Sclerosis Patients using Deep Learning Analysis of Prados, F., Cardoso, M.J., Kanber, B., Ciccarelli, O., Kapoor, R., Gandini Wheeler- MRI Data. Proc. Mach. Learn. Res. 102, 483–492. In: http://proceedings.mlr.pre Kingshott, C.A.M., Ourselin, S., 2016. A multi-time-point modality-agnostic patch- ss/v102/tousignant19a.html. based method for lesion filling in multiple sclerosis. Neuroimage 139, 376–384. Tur, C., Eshaghi, A., Altmann, D.R., Jenkins, T.M., Prados, F., Grussu, F., https://doi.org/10.1016/j.neuroimage.2016.06.053. Charalambous, T., Schmidt, A., Ourselin, S., Clayden, J.D., Wheeler-Kingshott, C.A. Rizzo, M.A., Hadjimichael, O.C., Preiningerova, J., Vollmer, T.L., 2004. Prevalence and M.G., Thompson, A.J., Ciccarelli, O., Toosy, A.T., 2018a. Structural cortical network treatment of spasticity reported by multiple sclerosis patients. Mult. Scler. J. 10 (5), reorganization associated with early conversion to multiple sclerosis. Sci. Rep. 8 (1), 589–595. https://doi.org/10.1191/1352458504ms1085oa. 10715. https://doi.org/10.1038/s41598-018-29017-1. Roca, P., Attye, A., Colas, L., Tucholka, A., Rubini, P., Cackowski, S., Ding, J., Budzik, J.- Tur, C., Moccia, M., Barkhof, F., Chataway, J., Sastre-Garriga, J., Thompson, A.J., F., Renard, F., Doyle, S., Barbier, E.L., Bousaid, I., Casey, R., Vukusic, S., Lassau, N., Ciccarelli, O., 2018b. Assessing treatment outcomes in multiple sclerosis trials and in Verclytte, S., Cotton, F., Brochet, B., Casey, R., Cotton, F., De Sèze, J., Douek, P., the clinical setting. Nat. Rev. Neurol. 14 (2), 75–93. https://doi.org/10.1038/ Guillemin, F., Laplaud, D., Lebrun-Frenay, C., Mansuy, L., Moreau, T., Olaiz, J., nrneurol.2017.171. Pelletier, J., Rigaud-Bully, C., Stankoff, B., Vukusic, S., Marignier, R., Tur, C., Kanber, B., Eshaghi, A., Altmann, D.R., Khaleeli, Z., Prados, F., Ourselin, S., Debouverie, M., Edan, G., Ciron, J., Ruet, A., Collongues, N., Lubetzki, C., Thompson, A.J., Gandini Wheeler-Kingshott, C.A.M., Toosy, A.T., Ciccarelli, O., Vermersch, P., Labauge, P., Defer, G., Cohen, M., Fromont, A., Wiertlewsky, S., 2019. Clinical relevance of cortical network dynamics in early primary progressive Berger, E., Clavelou, P., Audoin, B., Giannesini, C., Gout, O., Thouvenot, E., MS. Mult. Scler. J. 26 (4), 442–456. https://doi.org/10.1177/1352458519831400. Heinzlef, O., Al-Khedr, A., Bourre, B., Casez, O., Cabre, P., Montcuquet, A., Tur, C., Grussu, F., de Angelis, F., Prados, F., Kanber, B., Calvi, A., Eshaghi, A., Créange, A., Camdessanché, J.-P., Faure, J., Maurousset, A., Patry, I., Charalambous, T., Cortese, R., Chard, D.T., Chataway, J., Thompson, A.J., Hankiewicz, K., Pottier, C., Maubeuge, N., Labeyrie, C., Nifle, C., Ameli, R., Ciccarelli, O., Gandini Wheeler-Kingshott, C.A.M., 2022. Spatial patterns of brain Anxionnat, R., Attye, A., Bannier, E., Barillot, C., Ben Salem, D., Boncoeur- lesions assessed through covariance estimations of lesional voxels in multiple Martel, M.-P., Bonneville, F., Boutet, C., Brisset, J.-C., Cervenanski, F., Claise, B., Sclerosis: the SPACE-MS technique. NeuroImage: Clinical 33, 102904. https://doi. Commowick, O., Constans, J.-M., Dardel, P., Desal, H., Dousset, V., Durand- org/10.1016/J.NICL.2021.102904. Dubief, F., Ferre, J.-C., Gerardin, E., Glattard, T., Grand, S., Grenier, T., Guillevin, R., Tustison, N.J., Avants, B.B., Cook, P.A., Zheng, Y., Egan, A., Yushkevich, P.A., Gee, J.C., Guttmann, C., Krainik, A., Kremer, S., Lion, S., Menjot de Champfleur, N., 2010. N4ITK: Improved N3 Bias Correction. IEEE Trans. Med. Imaging 29 (6), Mondot, L., Outteryck, O., Pyatigorskaya, N., Pruvo, J.-P., Rabaste, S., Ranjeva, J.-P., 1310–1320. https://doi.org/10.1109/TMI.2010.2046908. Roch, J.-A., Sadik, J.C., Sappey-Marinier, D., Savatovsky, J., Tanguy, J.-Y., van der Burgh, H.K., Schmidt, R., Westeneng, H.J., de Reus, M.A., van den Berg, L.H., van Tourbah, A., Tourdias, T., 2020. Artificial intelligence to predict clinical disability in den Heuvel, M.P., 2017. Deep learning predictions of survival based on MRI in patients with multiple sclerosis using FLAIR MRI. Diagn. Interv. Imaging 101 (12), amyotrophic lateral sclerosis. NeuroImage: Clinical 13, 361–369. https://doi.org/ 795–802. 10.1016/j.nicl.2016.10.008. Rogers, B.P., Parks, M.H., Nickel, M.K., Katwal, S.B., Martin, P.R., 2012. Reduced fronto- Venugopalan, J., Tong, L., Hassanzadeh, H.R., Wang, M.D., 2021. Multimodal deep cerebellar functional connectivity in chronic alcoholic patients. Alcohol Clin Exp Res learning models for early detection of Alzheimer’s disease stage. Sci. Rep. 11 (1), 36 (2), 294–301. https://doi.org/10.1111/j.1530-0277.2011.01614.x. 1–13. https://doi.org/10.1038/s41598-020-74399-w. Sastre-Garriga, J., Pareto, D., Battaglini, M., Rocca, M.A., Ciccarelli, O., Enzinger, C., Wang, Y., Zhu, J., Li, Q., Li, W., Wu, N., Zheng, Y., Chang, H., Chen, J., Wang, W., Wuerfel, J., Sormani, M.P., Barkhof, F., Yousry, T.A., De Stefano, N., Tintoré, M., Soriano-Mas, C., 2013. Altered fronto-striatal and fronto-cerebellar circuits in Filippi, M., Gasperini, C., Kappos, L., Río, J., Frederiksen, J., Palace, J., Vrenken, H., heroin-dependent individuals: a resting-state fMRI study. PLoS One 8 (3), e58098. 12