A data-driven model of disability progression in progressive multiple sclerosis.

See P. Oxtoby and Barkhof Leckey and Zetterberg ( https://doi.org/10.1093/braincomms/fcae474 ) for a scientific commentary on this article.

Introduction

Despite significant advancements in understanding the underlying mechanisms of multiple sclerosis (MS) and the availability of treatments to prevent relapses, there remains an unmet need to halt and reverse disability accrual. Accumulated evidence suggests that disability accrual in MS is not driven by a single, uniform disease mechanism but rather by a combination of several mechanisms that vary among patients and even within individual patients over time. 1 Quantifying disability accrual presents a notable challenge, particularly in the context of progressive multiple sclerosis (PMS) patients. The heterogeneous nature of disability manifestations, 2-10 the prolonged progression spanning years and the limitations inherent in the measures employed to quantify progression collectively pose a significant challenge in the analysis of PMS data. In the face of the persistent unknowns surrounding disease progression, elucidating the timeline of disability accrual holds the promise of unveiling crucial insights into the underlying mechanisms of disease progression, thereby facilitating patient stratification and monitoring in clinical trials.

Data-driven disease progression models 11 are a set of computational tools specifically designed to quantify the long-term timeline of progression from short-term multimodal clinical data. These models can reveal disability patterns in a wide variety of sporadic and familial neurological diseases including Alzheimer’s disease, 12-21 Huntington’s disease, 22 , 23 Parkinson’s disease, 24 frontotemporal dementia 25 and progressive supranuclear palsy. 26 , 27 Two recent studies involving PMS patients 28 , 29 determined the most likely sequence in which brain regions become atrophic; however, they analysed only MRI structural metrics and did not include features derived from other modalities, such as disability or cognition measures, or lesion volumes. In a more recent study, 30 the authors considered a broader set of structural, functional and cognitive outcomes, on a limited dataset of relapsing-onset MS subjects.

This study aims to elucidate the typical pattern of disability accrual in PMS by analysing a large dataset of randomized clinical trial data encompassing longitudinal multimodal disability assessments (walking ability, upper limb function and cognitive function) and MRI measurements (lesion and brain volume) collected within an International Progressive MS Alliance project. We employ Gaussian process progression model (GPPM), a probabilistic data-driven disease progression model, previously used to uncover disease progression in Alzheimer’s Disease 17-19 as well as longitudinal evolution of ageing. 20 The aim of this analysis is to better characterize the typical progression of primary progressive MS (PPMS) and identify potential sub-groups of outliers with respect to such progression.

Materials and methods

Training dataset

Three randomized clinical trials of participants with PPMS collected within the International Progressive MS Alliance project (IPMSA, award reference number PA-1603–08175) were included in this analysis: ARPEGGIO, OLYMPUS and ORATORIO ( ClinicalTrials.gov numbers, NCT02284568 , NCT00087529 and NCT01194570 respectively). The study design and inclusion/exclusion criteria of the trials were previously described. 31-33 The ethics committees and institutional review boards of all participating centres approved the study protocols. All participants provided written informed consent.

Longitudinal metrics were analysed and included the expanded disability status scale (EDSS), the Timed 25-Foot Walk test (T25FW), the 9-Hole Peg Test (9HPT) (dominant and non-dominant hand), the Symbol Digit Modalities test (SDMT, a cognitive test of attention), and MRI markers including T1 and T2 lesion volume and normalized brain volume (NBV). Participants with missing disease duration were excluded from the analysis (22 participants from ORATORIO trial and one each from the ARPEGGIO and OLYMPUS trials). The variables that were not normally distributed (lesion volumes, T25FW and 9HPT measures) were log-transformed for the analysis. Finally, all the measurements were converted into quantile scores (0–1 for normal to abnormal values), according to their distribution. Disease duration was collected as well, and is considered from the onset of symptoms (as reported by the participants).

Test dataset

The PPMS participants ( N = 227) from the SPI2 study ( ClinicalTrials.gov number NCT02936037 ) were used as an independent external validation set. 34 The same variables as the training set were collected, except for 9HPT measurements that were not assessed in this trial.

Disease progression modelling using GPPM

GPPM, a data-driven Bayesian progression method, was used to model PPMS progression on the training dataset. GPPM is a mixed-effect model designed to handle repeated measures data by incorporating both fixed effects (e.g. the progression of the pathology common to all individuals) and random effects (e.g. individual-specific variability), as well as scores measuring the individual rates of progression, relative to the population. 17 The model does so by assuming that measurements are taken at an unknown ‘disease-time’, represented by a hidden variable to be estimated for each individual, while fitting the set of group-wise trajectories across subjects. The resulting disease-time axis parametrizes the natural evolution of the disability common to all individuals, while the individual disease-time scores (also called ‘time reparameterization parameters’) are measures of individual progression severity with respect to the common disease-time.

Formally y s ( t 1 ) , … , y s ( t T s ) represent the longitudinal multimodal measurements associated with each individual s at corresponding study time points ( t 1 , … , t T s ). It is then possible to consider all measurements obtained at a particular visit of one individual s to occur at an ideal disease time τ j s , where j s = 1 , … , T s and the mapping from study time t to disease time τ is modelled by a subject-specific time reparametrization function, assumed to be a time-shift: t j s = τ j s + d s , for some parameters d s . Under this assumption, the observations (for all clinical and MRI progression markers k = 1 , … , K ) for subject s at a given time point t (indices omitted but implied) can be modelled as a random sample from the model:

Here, f ( t ) = ( f 1 ( t ) , … , f K ( t ) ) is a Gaussian Process; η s ( t ) = ( η 1 s ( t ) , … , η K s ( t ) ) are individual Gaussian random effects and ε = ( ε 1 , … ε K ) is the measurement noise, assumed to be Gaussian heteroskedastic. The model was previously described in detail together with the optimization scheme to recover the probabilistic estimates of the parameters for the fixed effect, the random effect, and the individual time reparameterization parameters. 17 , 18 , 35 Identifiability of the model is ensured by enforcing global monotonicity on the population-level trajectories. Also, as the covariance of the Gaussian Process can be specified in order to account for incomplete data, the model naturally accounts for missing information and no data imputation is needed. Explicit estimations of the predictive performance of the model in assessing the individual parameters with respect to follow-up assessments and missing values can be found in. 17 Further, the Bayesian formalism naturally allows for probabilistic estimation of the disability trajectories and uncertainty quantification of the predicted individual time reparametrization parameters, without the need for cross-validation or bootstrapping. Finally, Gaussian processes naturally allow for probabilistic predictions of unseen individuals given the individual clinical and MRI progression marker measurements of the training sample. 17

In this study, GPPM was run on the training set and used the observed participants’ disease duration (ODD) at visit as reference for input study time frame t . Consequently, GPPM estimates the individual reparametrized disease durations (RDD) as measures of progression stages, representing a proxy for the disease time frame τ . Because the model requires monotonicity constraints to fit population-level trajectories, monotonic non-increasing constraints were enforced on SDMT and NBV, and non-decreasing constraints on the six remaining clinical and MRI progression markers.

Outliers identification using reparametrized disease duration

Following GPPM modelling on the training dataset, the individual time reparameterization parameters were examined, which measure the differences between the individual ODD and the model’s RDD. By doing so, one ‘normative’ sub-group was identified, in which the individuals had an estimated time reparameterization parameter, which was within 1 SD from the population level average parameter. Additionally, two ‘outlier groups’ with respect to the normative progression estimated by the model were identified, i.e. individuals with an estimated time reparameterization parameter, which was above or below 1 SD from the population level average parameter. Specifically, these three groups were denoted as:

Slow progressors (G1)—participants with an associated reparameterization parameter falling below 1 SD from the population-level average parameter.

Normative progressors (G2)—participants whose time reparameterization parameter is within ±1 SD of the population-level average parameter.

Fast progressors (G3)—participants with an associated reparameterization parameter exceeding 1 SD from the population-level average parameter.

Statistical analysis

The demographic and clinical characteristics of the three populations defined by the GPPM were described using appropriate measures such as mean (standard deviation [SD]), median (interquartile range [IQR]), or frequencies (percentage), depending on the data distribution. These characteristics were compared among sub-groups by employing a χ 2 test for categorical variables, and utilizing ANOVA or Kruskal–Wallis test, depending on the nature of the variables, for continuous measures.

To assess whether using the ODD alone was effective to identify the normative and outlier subpopulations, a ‘null model’ was implemented. This involved classifying three distinct subpopulations based solely on the ODD, where the average ODD for the entire population was computed along with its SD. Specifically, three ‘null’ sub-groups were defined (following what we did for G1, G2 and G3): a normative ‘null’ sub-group, in which the individuals had an ODD within 1 SD from the population-level average ODD; and two outliers ‘null’ groups, in which the individuals had an ODD above or below 1 SD from the population-level ODD. Demographic and clinical characteristics of the three ‘null’ populations and their 3-month confirmed disability progression (CDP) were then assessed and compared to what was obtained when looking at the three populations defined by the GPPM.

A Cox proportional hazards model was used to evaluate if there was a difference in the risk of 3-month CDP computed on EDSS, T25FW, SDMT and 9HPT (both dominant and non-dominant), among the three populations defined by the GPPM. CDP calculated on EDSS was defined as a 1-point increase in EDSS if the baseline score was 3.0–5.0 or a 0.5-point increase if the baseline score was 5.5–6.5, confirmed at a scheduled visit at least 3 months later. CDP computed on T25FW and 9HPT was defined as a worsening of at least 20% from baseline confirmed at a scheduled visit at least 3 months later. CDP computed on SDMT was defined as 6-month confirmed worsening of ≥3 points or 10% of the SDMT score. Age at symptoms onset, sex and treatment arm were included in the multivariable Cox regression model as covariates.

The median time and age at which an EDSS of 6 occurred were computed using the Kaplan–Meier approach. To compare the survival curves, the log-rank test was employed, and a two-sided P -value was reported without any correction for multiple comparisons.

Results

A total of 1521 PPMS participants (18 581 time-points) (373 from ARPEGGIO, 438 from OLYMPUS and 710 from ORATORIO) were analysed, with complete or partial availability of the 8 disease evolution measures. 998/1521 (65.61%) participants were treated (ARPEGGIO 233/373 (62.47%); OLYMPUS: 292/438 (66.67%); ORATORIO: 473/710 (66.2%)). Demographic and clinical characteristics of the training cohort are described in Table 1 . Mean follow-up duration of the pooled training dataset was 28.20 months (SD: 10.79 months), with a median frequency of EDSS assessments of 2.76 months (IQR: 2.66–2.79 months). Median frequency of 9HPT, T25FW and SDMT assessments was every 2.76 months (IQR: 2.70–2.79 months) for all three disease measures. MRI measurements have comparable median frequencies. SDMT values were collected solely in the ARPEGGIO trial, leading to a relatively small number of participants with this information available ( N = 363). Similarly, the volume of lesions in T1 was not available for the OLYMPUS trial and was available for the ORATORIO trial and for 58 participants in the ARPEGGIO trial.

Progression of primary progressive multiple sclerosis patients

Figure 1 top shows the estimated clinical and MRI progression markers trajectories, over an estimated disease-time axis of roughly 50 years (48.7 years). It also shows individual measurements, shifted according to the estimated RDD. Indeed, it is characterized at the initial stages by a sharp increase for EDSS, which then reaches a plateau after about 20–30 years of disease duration. The other disability markers (T25FW and both dominant and non-dominant 9HPT) get steadily worse along progression, similarly to lesion volumes. These latter measures are however heterogeneously distributed in the population and with very large variability. The evolution is further characterized by a steady worsening of SDMT, starting at roughly 25 years of disease duration, followed by NBV loss. These measures are again quite variable in the population, likely due to the large number of missing data (see Table 1 ).

Similar conclusions can be drawn by looking at Fig. 1 bottom-left, which shows the temporal ordering in which the clinical and MRI progression markers reach their highest rate of change along the estimated timeline. Except for lesion volumes, whose variability is spread across the entire timeline, the model predicts a sequence disease progression characterized by initial changes of EDSS, followed by T25FW and SDMT. Changes in 9HPT dominant and non-dominant appear at a later stage—the latter one with higher associated uncertainty—and finally NBV loss. Supplementary Fig. S1 supports this results using an alternative discrete disease progression model, 11 the z-score SuStain 13 , 36 constrained to a single progression type. The estimated sequence broadly aligns with that derived from GPPM, with EDSS and disability-related features appearing early, while the most significant changes in NBV, SDMT, and lesion markers occur later in the progression, although some initial decline is also observed.

Figure 1 bottom-right shows the global percentages of change in the estimated timeline of 48.7 years for six out of the eight markers: EDSS, NBV, SDMT, T25FW and 9HPT (dominant and non-dominant).

Clinical characteristics of the outlier groups with respect to the normative population

Figure 2 shows the distribution of the time reparametrization parameters (the difference between ODD and RDD) as estimated by the GPMM in the training set. Scatter plot of RDD versus ODD is shown in Supplementary Fig. S2 . We identify the presence of two outlier groups (in magenta and green in Fig. 2 ) with respect to the normative progression population (in blue in Fig. 2 ). Specifically, 167 participants (11%) were classified as slow progressors (G1), a sub-group of 193 participants (13%) were classified as fast progressors (G3), while the remaining 1161 participants (76%) were classified as normative progressors (G2).

Table 1 shows the heterogeneity in clinical presentation among the two outlier groups and the normative population. No differences were observed among the groups in terms of the rate of relapses or in the proportion of treated participants during the follow-up period. We note that fast progressors (G3) as compared to slow progressors (G1) are older at symptom onset (mean: 35.0 years versus 38.5 years, P < 0.001) and are males in a higher proportion (61.1% versus 48.5%, P = 0.013).

Prognostic differences between fast and slow progressors

Table 2 and Fig. 3 show a statistically significant difference in the rate of CDP computed on EDSS, T25FW, SDMT and 9HPT (both dominant and non-dominant), among the three populations defined by the GPPM (log-rank test for three-group comparison, P < 0.05 for all the disease activity measures). In particular, participants in the slow progressors group faced a 71% reduced risk of experiencing 3-month CDP as determined through EDSS when compared to those in the normative-progressor group (HR = 0.29, 95% CI: 0.17–0.47, P < 0.0001). Conversely, when we compared the fast progressors group to the normative-progressors group, the risk doubled (HR = 2.09, 95% CI: 1.66–2.62, P < 0.0001). Risk of CDP computed on T25FW and 9HPT is reduced by more than half when comparing slow progressors to normative-progressor individuals and, on the other hand, the risk increased by more than 3-fold comparing fast progressors to normative progressors (see Table 2 ). Only 42 CDP events were observed using SDMT as study outcome, with 3 occurring in the slow progressors and only 1 in the fast progressors. Such limited data precluded statistically robust comparisons between the groups.

To assess the potential for interactions between treatment and the sub-groups of participants defined by the GPPM on 3-month CDP, a sensitivity analysis was utilized with a treatment-by-sub-group interaction term as a covariate within the Cox regression model. Notably, this interaction term did not demonstrate statistical significance for any of the outcomes evaluated, even when the ARPEGGIO trial (the only trial without anti-CD20 drugs) was excluded (data not shown).

When classified in three populations using solely the ODD (the ‘null model’ described in ‘Statistical Analysis’), 263 (17.29%) of the participants were classified as slow progressors, 1021 (67.13%) as normative progressors and the remaining 237 (15.58%) as fast progressors. Differences in the rate of CDP calculated for the study outcomes among the three ‘null’ populations were not found ( Supplementary Table S2 ).

The median age at EDSS 6 was 53.5 years (IQR: 48.0–56.5 years) in G1 participants, 49.0 years (IQR: 43.0–53.0 years) in G2 participants, and 44.0 years (IQR: 38–51 years) in the G3 sub-group. Supplementary Fig. S3 shows that the median time from the onset of disease symptoms to reaching a confirmed EDSS 6 was 32.78 years (95% CI: 30.08–35.48 years), 11.62 years (95% CI: 10.95–12.29 years) and 6.75 years (95% CI: 5.83–7.68 years) for G1, G2 and G3, respectively (log-rank test for three-group comparison P < 0.0001).

Model testing in the external independent dataset

A total of 227 PPMS participants from the SPI2 trial were used as an external testing set and their demographic and clinical characteristics are detailed in Supplementary Table S1 . The relapse rate and T1 lesion volume were similar between the training and the external testing datasets, while the population in SPI2 trial appeared to be slightly older and with a longer disease duration than the training set.

The predicted RDD of the testing individuals was computed, based on the trained model, estimating a probabilistic model for the individual temporal staging and computing its expectation. 17 Then, following what was done on the training set, the normative and outlier sub-groups of progression were identified on the test set by examining the difference between individual ODD and model’s RDD. 25/227 (11.01%) late slow progressors (G1); 37/227 (16.30%) fast progressors (G3); and 165/227 (72.69%) normative progressors (G2) were identified.

Table 2 shows the difference in rate of 3-month CDP computed on EDSS, T25FW and SDMT, among the three test set sub-groups identified by the model. The results were all consistent with those found on the training set, even if the statistical significance was not reached, likely because of the small size of the test sample. Disease duration, EDSS, T25FW and SDMT (all at baseline) are significantly different in the sub-groups. Further, the slow progressors face a reduced risk of experiencing 3-month CDP determined through EDSS (HR = 0.21), while fast progressors have an increased risk (HR = 1.45). Similar effect sizes were obtained when determining the 3-month CDP through T25FW.

Discussion

In this study, we identified the group-wise, multimodal longitudinal evolution timeline of disability, cognition and neurodegeneration in PPMS along an estimated timeline of roughly 50 years. To do so, we employed GPPM, a probabilistic data-driven disease progression model, on a large dataset encompassing longitudinal multimodal disability assessments (walking ability, upper limb function and cognitive function) and MRI measurements (lesion and brain volume).

The outcomes of our investigation revealed compelling insights, compatible with previous findings in longitudinal studies in PPMS. First, the EDSS reaches its highest rate of change at around 15 years along the estimated progression, plateauing after ∼30 years of disease duration. This suggests that the utility of EDSS diminishes in monitoring disease progression beyond this temporal threshold, advocating for the research of potentially superior alternative measures in later stages. Second, the 9-Hole Peg Test performed with the dominant hand displayed less variability compared to the non-dominant hand in the population, underscoring the significance of considering hand dominance in evaluating disease-related functional changes. 37 NBV emerged as the metric with the slowest progression to its peak rate of change along the disease trajectory. To elucidate this phenomenon, it is crucial to examine NBV in relation to clinical deficits. Notably, tissue damage extend beyond primary lesions, encompassing injuries and damage in seemingly normal white matter (intralesional). 38 This observation underscores the intricate interplay between structural changes and functional impairments. Furthermore, NBV emerges as a shared final pathway implicated in diverse processes, including the secondary damage associated with lesions. 39 Finally, our model suggests that, although the increase in lesion volume over time is highly variable within the population, it remains relatively constant overall. This finding is consistent with studies on PPMS which also show a relatively constant increase in lesion volume over time, particularly during the early stages of the disease. 40 These considerations highlight the multifaceted nature of neurodegenerative mechanisms and underscores the necessity for further investigation into the temporal dynamics of brain volume alterations correlated with clinical manifestations.

Our model allowed the data-driven identification of two subsets of outliers, representing individuals with either fast or slow progressions relative to the normative trajectory. Slow progressors can be interpreted as those individuals who, despite having a long disease duration, exhibit a level of disability comparable to those with a shorter disease duration—in other words, they appear to be performing better than normative progressors with the same disease duration. Conversely, fast progressors are participants who, despite a short disease duration, display a level of disability akin to those with a longer disease duration, as the GPPM aligns them accordingly. This observation tentatively suggests the presence of disease processes in these sub-groups that evolve at rates differing from the normative population. However, more targeted studies are required to elucidate the underlying pathophysiological mechanisms in these groups.

Distinctive characteristics were observed between these sub-groups in terms of baseline (defined as the time of enrolment in the trial) clinical and demographic features. Fast progressors (which have a shorted disease duration) were characterized by an older age at symptom onset and a higher male sex as compared to slow progressors. This statement is consistent with a prior study on PPMS participants, 10 indicating a tendency for a higher prevalence of males and older individuals at symptom onset among those with more severe disability. Prognostically, the two sub-groups demonstrated divergent outcomes based on the 3-month CDP evaluated on EDSS, with fast progressors exhibiting a significantly worse prognosis.

Despite identifying two subsets of outliers based on their progression rate with respect to the normative population, Fig. 1 shows considerable variance in the estimated GPPM trajectories/sequence, suggesting potential existence of multiple sub-groups within the overall population, each potentially following distinct (sub)group-level trajectories. This observation suggests two key points: (i) the primary clinical utility of the current GPMM version lies in group-level applications such as clinical trials rather than individual monitoring and (ii) exploring these subtypes using a clustering-based version of GPPM could be a promising direction for future research 13 , 41 (for clustering-based extensions of different disease progression models).

Machine learning models were previously applied to MS clinical trial data, to identify sub-groups of patients with different prognosis or different propensity to respond to treatment. The first paper applied an unsupervised machine learning algorithm, called Subtype and Staging Inference (SuStaIn), 29 to uncover data-driven assessment of the pathological changes visible on a large pooled dataset of MRI scans. They showed that it is possible to identify three data-driven MS subtypes that do not coincide with traditional MS classifications (relapsing remitting, primary progressive, etc.) and that are much more strongly predictive of response to a particular treatment. While most of the clinical trials in PMS failed to detect any relevant treatment effect in the whole populations enrolled in these trials, the ‘lesion-led’ subtype showed a significant treatment response in the three trials in PMS ( n = 2900) included in the analysis as the test set. However, another study observed that although the lesion-led subtype was associated with increased inflammatory activity and correlated with disease severity and functional impairment at baseline, this subtype did not reliably predict disability progression nor effectively differentiate treatment response heterogeneity. 42 Conversely, the sub-groups delineated using GPPM displayed divergent outcomes based on 3-month CDP, yet they were unable to discern treatment effects. This finding aligns with individual clinical trial results, where the treatment effect estimated on clinical endpoints in the PPMS population was either absent or weak at best. These observations suggest that the molecular mechanisms targeted by the specific drugs in these clinical trials may have limited impact on disease progression rates, which is the case for Laquinimod (ARPEGGIO trial 31 ). When we removed the ARPEGGIO trial and repeated the analysis we still observed no significant interaction for Rituximab and Ocrelizumab (ORATORIO and OLYMPUS trials). This indicates that while GPPM may not predict treatment response when it does not impact disease progression rate, it remains a valuable tool for understanding disease progression itself.

A novel deep learning framework was recently used to estimate an individual’s treatment effect on disability progression using readily available clinical information (demographic characteristics and clinical disability scores) and scalar MRI metrics (lesional and volumetric) collected in clinical trials. 43 An intriguing result was that in this experiment a non-linear model (multi-layers perceptron) outperformed other linear (and log-linear) combinations of baseline variables, suggesting that complex relationships exist between the baseline features and treatment effect. Recently, an analysis using Bayesian deep learning was used for estimating the posterior distribution over factual and counterfactual outcomes on several treatments. 44 The model was trained and tested to predict future new and enlarging T2 lesion counts on the Progressive MS Alliance dataset.

The present study has important limitations, including the relatively small size of the test set that contributed to the lack of statistical significance in testing on the external dataset. Addressing this limitation requires examining more extensive clinical trial data. One limitation of our model is the enforced global monotonicity of population-level trajectories to ensure identifiability. While this allows for small fluctuations within the whole timeline, it may not fully capture the non-monotonic nature of disease progression, even in PPMS, where progression is generally considered monotonic. PPMS patients can still experience significant fluctuations in disease activity, and enforcing global monotonicity may oversimplify these dynamics, potentially affecting the model’s accuracy. Future work should explore more flexible trajectory assumptions while maintaining identifiability to better reflect the heterogeneity and variability of disease progression. 18 Another significant limitation is the use of the same markers (notably, the EDSS) for computing both the RDD and subsequent a posteriori statistical tests (CDP). This dual use introduces potential circularity, where EDSS scores used for initial group categorization may bias outcomes assessment. While ensuring internal consistency, this approach limits generalizability and robustness. Furthermore, the limitations of the EDSS are well known, particularly its reduced sensitivity to the clinical accrual of non-motor disability in subjects with progressive forms of the disease, compared to newer scales like the MS Functional Composite. 45 However, we focused on the EDSS here because it remains the most widely used metric of overall disability in clinical trials. Future research, including data from broader datasets, will explore alternative metrics for a posteriori statistical testing to mitigate these limitations and enhance the study’s reliability.

Another limitation of this study is the use of +/− 1 SD as a cutoff point to define fast and slow progressor groups. We note that GPPM models the distribution of the time reparametrization parameters as Gaussian, simplifying the identification of the tails of the distribution by means of SDs or quantiles and obviating the need for outlier detection techniques (e.g. clustering). Nonetheless, it is possible to seek an optimized threshold by, for instance, maximizing the separation between the ‘outlier’ populations and the normative one. This would, however, introduce a bias on the analysis of such differences, and thus, we have not explored this strategy. Other ‘hard’ cutoff values (1.5 and 2 SD) were explored, and the results remained consistent, although the sample size was too small to draw definitive conclusions (see Supplementary Table S3 ). It is crucial to note that further research with larger samples may be necessary to confirm these findings and establish more robust criteria for defining progression groups.

Another noteworthy aspect for consideration pertains to the absence of ‘advanced’ metrics, such as MRI metrics related to the smoldering disease, 46 , 47 chronic inflammation, 48 regional brain and spinal cord atrophy 49 and emerging fluid biomarkers. 50 The rationale behind this omission stems from the inherent challenges associated with quantifying these parameters in randomized clinical trials. Despite these limitations, it is fundamental to acknowledge that the exclusion of these advanced measures may impact the comprehensive understanding of disease dynamics. Future research may benefit from employing innovative methodologies or collaborative frameworks to overcome these challenges and incorporate these critical metrics, thus contributing to a more nuanced and comprehensive understanding of disease progression.

In summary, the group-wise, multimodal longitudinal evolution timeline of disability, cognition and neurodegeneration in PPMS was described along an estimated timeline of roughly 50 years leading to identification of two sets of outliers experiencing either hastened or slowed progressions with respect to the normative population.

Supplementary Material

Acknowledgements

The authors are grateful to the companies who generously provided the clinical trial data that made this work possible: Hoffmann-La Roche/Genentech, Biogen and Teva. MedDay was the sponsor for the SPI2 clinical trial whose dataset is now curated at the University of California San Francisco.

Contributor Information

Sara Garbarino, Life Science Computational laboratory, IRCCS Ospedale Policlinico San Martino, 16132 Genoa, Italy; MIDA, Dipartimento di Matematica, Università di Genova, 16146 Genoa, Italy.

Carmen Tur, Multiple Sclerosis Centre of Catalonia (Cemcat), Department of Neurology, Hospital Universitari Vall d’Hebron, Universitat Autònoma de Barcelona, 08035 Barcelona, Spain.

Marco Lorenzi, Universitè Côte d’Azur, Inria, Epione Research Project, 06902 Sophia Antipolis, France.

Matteo Pardini, Department of Neuroscience, Rehabilitation, Ophthalmology, Genetics, Maternal and Child Health, Università di Genova, 16132 Genoa, Italy; IRCCS Ospedale Policlinico San Martino, 16132 Genoa, Italy.

Michele Piana, Life Science Computational laboratory, IRCCS Ospedale Policlinico San Martino, 16132 Genoa, Italy; MIDA, Dipartimento di Matematica, Università di Genova, 16146 Genoa, Italy.

Antonio Uccelli, Department of Neuroscience, Rehabilitation, Ophthalmology, Genetics, Maternal and Child Health, Università di Genova, 16132 Genoa, Italy; IRCCS Ospedale Policlinico San Martino, 16132 Genoa, Italy.

Douglas L Arnold, McGill University, H3A0G4 Montréal, Canada.

Bruce A C Cree, UCSF Weill Institute for Neurosciences, Department of Neurology, University of California San Francisco, 94143 San Francisco, USA.

Maria Pia Sormani, IRCCS Ospedale Policlinico San Martino, 16132 Genoa, Italy; Department of Health Sciences (DISSAL), Università di Genova, 16132 Genoa, Italy.

Francesca Bovis, Department of Health Sciences (DISSAL), Università di Genova, 16132 Genoa, Italy.

Supplementary material

Supplementary material is available at Brain Communications online.

Funding

This investigation was supported by an award from the International Progressive Multiple Sclerosis Alliance (award reference number PA-1412-02420). This work was supported by #NEXTGENERATIONEU (NGEU) and funded by the Ministero dell’Università e delle Ricerca, Italy (MUR), National Recovery and Resilience Plan (NRRP), project MNESYS (PE0000006)—A Multiscale integrated approach to the study of the nervous system in health and disease (DN. 1553 11.10.2022). The authors acknowledge financial contribution from the Italian Ministry of Health, with the project NeuroArtP3 (NET-2018-12366666). S.G. acknowledges financial contributions from the Istituto Nazionale di Alta Matematica—Gruppo Nazionale di Calcolo Scientifico (INdAM–GNCS) Project (CUP_E53C22001930001). C.T. is currently being funded by a Miguel Servet contract, awarded by the Instituto de Salud Carlos III (ISCIII), Ministerio de Ciencia e Innovación de España (CP23/00117). She has also received a 2023 FORTALECE grant, awarded by the ISCIII (FORT23/00034) and a 2020 Junior Leader La Caixa Fellowship (fellowship code: LCF/BQ/PI20/11760008), awarded by ‘la Caixa’ Foundation (ID 100010434), a 2021 Merck’s Award for the Investigation in multiple sclerosis, awarded by Fundación Merck Salud (Spain), and a Research Grant awarded by the ISCIII, Ministerio de Ciencia e Innovación de España (PI21/01860). M.L. acknowledges funding from the Michael J. Fox Foundation for Parkinson’s Research (grant ID MJFF-021683). F.B. received the 2022 Biostatistic/Informatics Junior Faculty Award (grant code BI-2107-38160) awarded by the National Multiple Sclerosis Society.

Competing interests

S.G., F.B., ML, MPi and AU and have nothing to disclose. CT has received honoraria from Roche, Novartis, Merck, and Bristol Myers Squibb, and is a steering committee member of the O’HAND trial and of the Consensus group on Follow-on DMTs. MPa received research support from Novartis and Nutricia and honoraria from Biogen, Eisai, Merck and Novartis. DLA reports consulting fees from Biogen, Celgene, Frequency Therapeutics, Genentech, Merck, Novartis, Race to Erase MS, Roche, and Sanofi-Aventis, Shionogi, Xfacto Communications, grants from Immunotec and Novartis, and an equity interest in NeuroRx. BACC has received consulting fees from Alexion, Atara, Autobahn, Avotres, Biogen, Boston Pharma, EMD Serono, Gossamer Bio, Hexal/Sandoz, Horizon, Immunic AG, Kyverna, Neuron23, Novartis, Sanofi and TG Therapeutics and research support from Genentech. MPS has received consulting fees from Biogen, Genzyme, GeNeuro, MedDay, Merck, Novartis, Roche and Teva.

Data availability

Data used in this study were collected within the International Progressive MS Alliance project (IPMSA, award reference number PA-1603-08175). Access requests should be forwarded to the relevant data controllers. SPI2 study data were made available by the SPI-2 investigators and are available upon request to qualified researchers. Source code for the GPPM algorithm is freely available at https://gitlab.inria.fr/epione/GP_progression_model_V2.git . This repository does not contain dataset-specific code, since the data we used are not publicly available.

References

1. Kuhlmann T, Moccia M, Coetzee T, et al. Time for a new mechanism-driven framework to define multiple sclerosis progression. Lancet Neurol. 2023;22(1):78–88.

2. Kremenchutzky M, Rice GPA, Baskerville J, Wingerchuk DM, Ebers GC. The natural history of multiple sclerosis: A geographically based study 9: Observations on the progressive phase of the disease. Brain. 2006;129(3):584–594.

3. Tremlett H, Paty D, Devonshire V. The natural history of primary progressive MS in British Columbia, Canada. Neurology. 2005;65(12):1919–1923.

4. Soldán MMP, Novotna M, Zeid NA, et al. Relapses and disability accumulation in progressive multiple sclerosis. Neurology. 2015;84(1):81–88.

5. Confavreux C, Vukusic S. Age at disability milestones in multiple sclerosis. Brain. 2006;129(Pt 3):595–605.

6. Koch M, Kingwell E, Rieckmann P, Tremlett H. The natural history of primary progressive multiple sclerosis. Neurology. 2009;73(23):1996–2002.

7. Stellmann JP, Neuhaus A, Lederer C, Daumer M, Heesen C. Validating predictors of disease progression in a large cohort of primary-progressive multiple sclerosis based on a systematic literature review. PLoS One. 2014;9(3):e92761.

8. Koch MW, Greenfield J, Javizian O, Deighton S, Wall W, Metz LM. The natural history of early versus late disability accumulation in primary progressive MS. J Neurol Neurosurg Psychiatry. 2015;86(6):615–621.

9. Harding KE, Wardle M, Moore P, et al. Modelling the natural history of primary progressive multiple sclerosis. J Neurol Neurosurg Psychiatry. 2015;86(1):13–19.

10. Signori A, Izquierdo G, Lugaresi A, et al. Long-term disability trajectories in primary progressive MS patients: A latent class growth analysis. Mult Scler. 2018;24(5):642–652.

11. Young AL, Oxtoby NP, Garbarino S, et al. Data-driven modelling of neurodegenerative disease progression: Thinking outside the black box. Nat Rev Neurosci. 2024;25(2):111–130.

12. Young AL, Oxtoby NP, Daga P, et al. A data-driven model of biomarker changes in sporadic Alzheimer’s disease. Brain. 2014;137(9):2564–2577.

13. Young AL, Marinescu RV, Oxtoby NP, et al. Uncovering the heterogeneity and temporal complexity of neurodegenerative diseases with subtype and stage inference. Nat Commun. 2018;9(1):4273.

14. Jedynak BM, Lang A, Liu B, et al. A computational neurodegenerative disease progression score: Method and results with the Alzheimer’s disease neuroimaging initiative cohort. Neuroimage. 2012;63(3):1478–1486.

15. Donohue MC, Jacqmin-Gadda H, Le Goff M, et al. Estimating long-term multivariate progression from short-term data. Alzheimers Dement. 2014;10(5):S400–S410.

16. Oxtoby NP, Garbarino S, Firth NC, Warren JD, Schott JM, Alexander DC. Data-driven sequence of changes to anatomical brain connectivity in sporadic Alzheimer’s disease. Front Neurol. 2017;8(NOV):580.

17. Lorenzi M, Filippone M, Frisoni GB, Alexander DC, Ourselin S. Alzheimer’s disease neuroimaging initiative. Probabilistic disease progression modeling to characterize diagnostic uncertainty: Application to staging and prediction in Alzheimer’s disease. Neuroimage. 2019;190:56–68.

18. Garbarino S, Lorenzi M. Investigating hypotheses of neurodegeneration by learning dynamical systems of protein propagation in the brain. Neuroimage. 2021;235:117980.

19. Garbarino S, Lorenzi M. The Alzheimer’s disease neuroimaging initiative for. Modeling and inference of spatio-temporal protein dynamics across brain networks. Springer; 2019:57–69. doi: 10.1007/978-3-030-20351-1_5

20. Garbarino S, Lorenzi M, Oxtoby NP, et al. Differences in topological progression profile among neurodegenerative diseases from imaging data. Elife. 2019;8:e49298.

21. Koval I, Bône A, Louis M, et al. AD course map charts Alzheimer’s disease progression. Sci Rep. 2021;11(1):1–16.

22. Wijeratne PA, Garbarino S, Gregory S, et al. Revealing the timeline of structural MRI changes in premanifest to manifest Huntington disease. Neurol Genet. 2021;7(5):e617.

23. Wijeratne PA, Young AL, Oxtoby NP, et al. An image-based model of brain volume biomarker changes in Huntington’s disease. Ann Clin Transl Neurol. 2018;5(5):570–582.

24. Oxtoby NP, Leyland LA, Aksman LM, et al. Sequence of clinical and neurodegeneration events in Parkinson’s disease progression. Brain. 2021;144(3):975–988.

25. Young AL, Bocchetta M, Russell LL, et al. Characterizing the clinical features and atrophy patterns of MAPT-related frontotemporal dementia with disease progression modeling. Neurology. 2021;97(9):e941–e952.

26. Saito Y, Kamagata K, Wijeratne PA, et al. Temporal progression patterns of brain atrophy in corticobasal syndrome and progressive supranuclear palsy revealed by subtype and stage inference (SuStaIn). Front Neurol. 2022;13:814768.

27. Scotton WJ, Shand C, Todd E, et al. Uncovering spatiotemporal patterns of atrophy in progressive supranuclear palsy using unsupervised machine learning. Brain Commun. 2023;5(2):fcad048.

28. Eshaghi A, Marinescu RV, Young AL, et al. Progression of regional grey matter atrophy in multiple sclerosis. Brain. 2018;141(6):1665–1677.

29. Eshaghi A, Young AL, Wijeratne PA, et al. Identifying multiple sclerosis subtypes using unsupervised machine learning and MRI data. Nat Commun. 2021;12(1):1–12.

30. Dekker I, Schoonheim MM, Venkatraghavan V, et al. The sequence of structural, functional and cognitive changes in multiple sclerosis. Neuroimage Clin. 2021;29:102550.

31. Giovannoni G, Knappertz V, Steinerman JR, et al. A randomized, placebo-controlled, phase 2 trial of laquinimod in primary progressive multiple sclerosis. Neurology. 2020;95(8):E1027–E1040.

32. Hawker K, O’Connor P, Freedman MS, et al. Rituximab in patients with primary progressive multiple sclerosis: Results of a randomized double-blind placebo-controlled multicenter trial. Ann Neurol. 2009;66(4):460–471.

33. Montalban X, Hauser SL, Kappos L, et al. Ocrelizumab versus placebo in primary progressive multiple sclerosis. N Engl J Med. 2017;376(3):209–220.

34. Cree BAC, Cutter G, Wolinsky JS, et al. Safety and efficacy of MD1003 (high-dose biotin) in patients with progressive multiple sclerosis (SPI2): A randomised, double-blind, placebo-controlled, phase 3 trial. Lancet Neurol. 2020;19(12):988–997.

35. Lorenzi M, Filippone M.. Constraining the dynamics of deep probabilistic models. In: Proceedings of the 35th International Conference on Machine Learning. 2018:3227–3236.

36. Aksman LM, Wijeratne PA, Oxtoby NP, et al. Pysustain: A python implementation of the subtype and stage inference algorithm. SoftwareX. 2021;16:100811.

37. Lamers I, Timmermans AAA, Kerkhofs L, Severijns D, Van Wijmeersch B, Feys P. Self-reported use of the upper limbs related to clinical tests in persons with multiple sclerosis. Disabil Rehabil. 2013;35(23):2016–2020.

38. Kutzelnigg A, Lucchinetti CF, Stadelmann C, et al. Cortical demyelination and diffuse white matter injury in multiple sclerosis. Brain. 2005;128(11):2705–2712.

39. Lie IA, Weeda MM, Mattiesing RM, et al. Relationship between white matter lesions and gray matter atrophy in multiple sclerosis. Neurology. 2022;98(15):E1562–E1573.

40. Tur C, Khaleeli Z, Ciccarelli O, et al. Complementary roles of grey matter MTR and T2 lesions in predicting progression in early PPMS. J Neurol Neurosurg Psychiatry. 2016;82(4):423–428.

41. Poulet PE, Durrleman S. Mixture modeling for identifying subtypes in disease course mapping. In: Lecture notes in computer science (including subseries lecture notes in artificial intelligence and lecture notes in bioinformatics). Springer International Publishing; 2021:571–582. doi: 10.1007/978-3-030-78191-0_44/TABLES/2

42. Jiang X, Shen C, Caba B, et al. Assessing the utility of magnetic resonance imaging-based “SuStaIn” disease subtyping for precision medicine in relapsing-remitting and secondary progressive multiple sclerosis. Mult Scler Relat Disord. 2023;77:104869.

43. Falet JPR, Durso-Finley J, Nichyporuk B, et al. Estimating individual treatment effect on disability progression in multiple sclerosis using deep learning. Nat Commun. 2022;13(1):1–12.

44. Durso-Finley J, Arnold DL, Arbel T. Saliency based deep neural network for automatic detection of gadolinium-enhancing multiple sclerosis lesions in brain MRI. In: Lecture notes in computer science (including subseries lecture notes in artificial intelligence and lecture notes in bioinformatics). Springer International Publishing; 2020:108–118. doi: 10.1007/978-3-030-46640-4_11/FIGURES/6

45. Goldman MD, Larocca NG, Rudick RA, et al. Evaluation of multiple sclerosis disability outcome measures using pooled clinical trial data. Neurology. 2019;93(21):e1921–e1931.

46. Arnold DL, Belachew S, Gafson AR, Gaetano L, Bernasconi C, Elliott C. Slowly expanding lesions are a marker of progressive MS—No. Mult Scler. 2021;27(11):1681–1683.

47. Elliott C, Belachew S, Wolinsky JS, et al. Chronic white matter lesion activity predicts clinical progression in primary progressive multiple sclerosis. Brain. 2019;142(9):2787–2799.

48. Absinta M, Maric D, Gharagozloo M, et al. A lymphocyte–microglia–astrocyte axis in chronic active multiple sclerosis. Nature. 2021;597(7878):709–714.

49. Cagol A, Schaedelin S, Barakovic M, et al. Association of brain atrophy with disease progression independent of relapse activity in patients with relapsing multiple sclerosis. JAMA Neurol. 2022;79(7):682–692.

50. Jiang X, Shen C, Teunissen CE, et al. Glial fibrillary acidic protein and multiple sclerosis progression independent of acute inflammation. Mult Scler. 2023;29(9):1070–1079.

Associated Data

Supplementary Materials

Data Availability Statement

Data used in this study were collected within the International Progressive MS Alliance project (IPMSA, award reference number PA-1603-08175). Access requests should be forwarded to the relevant data controllers. SPI2 study data were made available by the SPI-2 investigators and are available upon request to qualified researchers. Source code for the GPPM algorithm is freely available at https://gitlab.inria.fr/epione/GP_progression_model_V2.git . This repository does not contain dataset-specific code, since the data we used are not publicly available.