Predicting disease activity in patients with multiple sclerosis: An explainable machine-learning approach in the Mavenclad trials.

Multiple sclerosis (MS) is among the most common autoimmune disabling neurological conditions of young adults and affects more than 2.3 million people worldwide. Predicting future disease activity in patients with MS based on their pathophysiology and current treatment is pivotal to orientate future treatment. In this respect, we used machine learning to predict disease activity status in patients with MS and identify the most predictive covariates of this activity. The analysis is conducted on a pooled population of 1935 patients enrolled in three cladribine tablets clinical trials with different outcomes: relapsing-remitting MS (from CLARITY and CLARITY-Extension trials) and patients experiencing a first demyelinating event (from the ORACLE-MS trial). We applied gradient-boosting (from XgBoost library) and Shapley Additive Explanations (SHAP) methods to identify patients' covariates that predict disease activity 3 and 6 months before their clinical observation, including patient baseline characteristics, longitudinal magnetic resonance imaging readouts, and neurological and laboratory measures. The most predictive covariates for early identification of disease activity in patients were found to be treatment duration, higher number of new combined unique active lesion count, higher number of new T1 hypointense black holes, and higher age-related MS severity score. The outcome of this analysis improves our understanding of the mechanism of onset of disease activity in patients with MS by allowing their early identification in clinical settings and prompting preventive measures, therapeutic interventions, or more frequent patient monitoring.

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

Classical pharmacometric methods are computationally intensive and have a limited ability to exploit high‐dimensional covariate data. Hence, incorporating machine‐learning (ML) approaches that are better suited for big data in clinical analysis has the potential to improve prediction in model‐informed drug development.

WHAT QUESTION DID THIS STUDY ADDRESS?

This study presents a predictive modeling approach for identifying patients with MS who will have an onset of disease activity in 3 and 6 months along with identification of clinical covariates that drive the model prediction using explainable ML methods.

WHAT DOES THIS STUDY ADD TO OUR KNOWLEDGE?

This study illustrates how incorporation of interpretable ML methods such as Shapley Additive Explanations (SHAP) with complex, nonlinear “black box” ML models, such as XgBoost into Drug Discovery and Development questions can lead to efficient exploration of high‐dimensional patient covariates and assess their contribution to composite clinical end points in MS. The interpretability methods such as SHAP make the decision process of ML models transparent, increasing trust in the models, whereas ensemble methods such as XgBoost enable improved prediction by capturing nonlinear effects and covariate interactions in a data set with a high number of diverse covariates.

HOW MIGHT THIS CHANGE DRUG DISCOVERY, DEVELOPMENT, AND/OR THERAPEUTICS?

The study demonstrates a framework to use complex “black box” ML models and explainability methods to analyze high‐dimensional data and present clinical insights and interpretable covariate importance to further enrich the currently available domain knowledge of MS disease progression and drug response mechanisms in clinical pharmacology.

INTRODUCTION

Multiple sclerosis (MS) is an inflammatory neurological condition that is the leading disabling disease in young adults , with compounded socioeconomic impact. The pathology and clinical presentation of MS have been widely studied. , , Cladribine tablets (Mavenclad®) represent an attractive treatment owing to its benefit–risk balance and convenience of use as an oral disease‐modifying drug (DMD) for use in MS treatment. , The recommended dose is 3.5 mg/kg, consisting of two annual courses, each comprising two treatment weeks 1 month apart and starting at the beginning of two consecutive years. Reviews of cladribine tablets for relapsing MS treatment can be found in Deeks and Rammohan et al.

Cladribine functions by preferentially reducing lymphocytes, key immune cells underlying MS pathogenesis. The clinical pharmacology of cladribine tablets have been investigated in detail including the pharmacokinetics, pharmacodynamics, and exposure–response relationships. Network meta‐analysis of 10,825 articles covering 44 studies assessing 12 disease‐modifying treatments (DMTs) showed that cladribine tablets are comparatively effective and a safe alternative to other DMTs in both active relapsing–remitting MS (RRMS) and high disease activity populations with respect to multiple clinical end points such as annualized relapse rate, confirmed disease progression, and no evidence of disease activity.

The efficacy and safety of cladribine tablets in the treatment of patients with MS in different stages of the evolution of the disease have been studied by integrating data from multiple phase III clinical trials. , The CLARITY (NCT00213135) study in patients with RRMS showed that annualized relapse rates and worsening sustained disability were reduced in patients treated with cladribine tablets compared with patients on placebo. , The efficacy observed in CLARITY was maintained without further active treatment during CLARITY‐Extension (NCT00641537). , ORACLE‐MS (NCT00725985) clinical data showed the efficacy of cladribine in delaying the conversion of patients with a first clinically isolated event to clinically definite MS (CDMS). A significant treatment effect of cladribine tablets was observed in patients who later converted to CDMS and were switched to a different DMD.

The pooled data set from the three trials results in a richly informative population with a long observation period of more than 6 years, with heterogeneous baseline and time‐varying covariates that can be used to explore the relationship of multiple efficacy end points with patient (baseline) characteristic and treatment effects. Incorporation of machine‐learning (ML) methods to efficiently handle the analysis of high‐dimensional data for model‐informed drug development has been widely supported and evidenced in the literature, , , , including in predicting the MS disease course in patients and the conversion to secondary progressive MS. , , , This multivariate data set presents an opportunity to characterize patients who will have disease activity in the future with data‐driven ML modeling. Earlier studies have assessed the relation of various baseline and time‐varying patient characteristics for overall clinical outcome through time‐to‐event modeling. , , , The focus of our current analysis is on the early identification of patients who will experience the onset of disease activity within a 6‐month period and explore which covariates contribute to their early identification.

In this article, we present our findings by applying interpretable ML methods such as Shapley Additive Explanations (SHAP) in combination with XGBoost, an ensemble tree‐boosting method, a state‐of‐the‐art approach for many ML challenges. , , These methods generate insights into the contribution of the input covariates into model predictions, thus rendering the decision process of the so‐called “black‐box” models transparent and informative for clinical patient management. To this end, we studied the top predictive covariates in our models affecting the future disease activity status of patients with MS.

METHODS

In this study, combining the patient populations from ORACLE‐MS, CLARITY, and CLARITY‐Extension, a disease activity event, for a patient while on cladribine treatment or placebo or observational follow‐up, is defined as meeting any of the following five criteria:

C1: one qualified relapse and at least one new T1 gadolinium enhancing (Gd+) lesion during the previous 48 weeks.

C2: one qualified relapse and at least two new or enlarging T2 lesions during the previous 48 weeks.

C3: two or more qualified relapses in absence of any magnetic resonance imaging (MRI) finding during the previous 48 weeks.

C4: 3‐month sustained Expanded Disability Status Scale (EDSS) progression.

C5: Required switching to an alternative DMT.

Sustained EDSS progression was defined as an increase in the EDSS score of ≥1 point if baseline EDSS was between ≥1.0 and ≤4.5, ≥1.5 points if baseline EDSS was 0, or ≥0.5 if baseline EDSS was ≥5.0 during a period of at least 3 months.

Although the first four criteria measure various dimensions of disease activity, the last criterion, C5, is a marker of treatment persistence and can be driven by lack of efficacy or a tolerability issue. Because of the complex and multidimensional aspect of disease progression in MS, there is a strong interest in studying composite end points that combine several clinical end points to address better treatment planning and MS patient monitoring. To note, our analysis involves all patients from the three cladribine trials, including those in the placebo arm not receiving cladribine.

Materials

In our analysis, we adopted a supervised ML framework which algorithmically learns the complex relationship between the high‐dimensional input covariate set with the dependent or output covariate. We trained and validated four models that make predictions of disease activity for patients approximately 3 months/12 weeks (T‐12) and 6 months/24 weeks (T‐24) in advance, as this matches the frequency of patient clinical visits and sampling schedules of various covariates in the trials. Our phase III T‐12 (P3‐T‐12) and T‐24 (P3‐T‐24) models are based on a set of 57 independent covariates, including various patient baseline and time‐varying covariates such as patient characteristics, neurological and MRI‐based assessments, laboratory measurements (refer to Table 1). The phase IV T‐12 (P4‐T‐12) and T‐24 (P4‐T‐24) models in contrast to the P3 models are based on 25 independent covariates, excluding laboratory covariates, to mimic the scenario of routine clinical practice or phase IV trials where laboratory covariates are not usually collected.

TABLE 1

Input covariates for the P3‐T‐24 and P3‐T‐12 models

Patient characteristics + baselines	Age, sex, race, dose (number of weeks of treatment), weight, age of onset of disease, time since first attack, lymphocytes_baseline, EDSS_baseline
Neurological assessment	Global Age‐Related Multiple Sclerosis Severity Score, KFSS1–Bowel and Bladder Functions, KFSS1–Brain Stem Functions, KFSS1–Cerebellar Functions, KFSS1–Cerebral or Mental Functions, KFSS1–Pyramidal Functions, KFSS1–Sensory Functions, KFSS1–Visual or Optic Functions
MRI assessment	Total number of T1 Gd+ lesions, total T1 hypointense (black holes), total number of T2/flair lesions, T1 Gd+ (volume in mm3), T1 hypointense lesions (volume in mm3), T2 lesions (volume in mm3), combined unique lesion count, new T1 hypointense (black holes)
Laboratory	Biochemistry: alanine aminotransferase, albumin, alkaline phosphatase, aspartate aminotransferase, bilirubin, blood urea nitrogen, calcium, creatine kinase, creatinine, sodium, potassium, urate, serum protein	Hematology: basophils, basophils/leukocytes, eosinophils, eosinophils/leukocytes, erythrocytes, hematocrit, hemoglobin, leukocytes, lymphocytes, lymphocytes/leukocytes, monocytes, monocytes/leukocytes, neutrophils, neutrophils/leukocytes, platelets	Urinalysis: urine pH, glucose

Note: The laboratory covariates are not collected in routine clinical practice. Hence, the input to P4‐T‐24 and P4‐T‐12 models have the same set of input covariates as P3 models except for the laboratory covariates.

Abbreviations: EDSS, Expanded Disability Status Scale; Gd+, gadolinium enhancing; KFSS, Kurtzke Functional Systems Scores; MRI, magnetic resonance imaging; P3‐T‐12, phase III 12 weeks; P3‐T‐24, phase III 24 weeks; P4‐T‐12, phase IV 12 weeks; P4‐T‐12, phase IV 24 weeks.

This covariate list includes those that were collected and can be matched across the combined population of 1935 individuals in the three trials. Covariates missing for greater than 20% of the patients are not included (e.g., lymphocyte subset CD4, CD8). This is because such covariates increase the dimension of the feature space without contributing to the predictive power of the model. These can further adversely affect generalizability of model performance. Last observation carried forward imputation was adopted to have a fixed number of input variables for the XGBoost algorithm. The input to the T‐12 models are baseline covariates and time‐varying covariates from 12 to 20 weeks before the first observation of disease activity for the patient or from the time of their last available observation record. Similarly, input to the T‐24 models are baseline covariates and time‐varying covariates from 21 to 30 weeks before the first observation of disease activity or from the time of their last available observation record. The representation of the cladribine dose to the model was chosen as the number of weeks of treatment received instead of cumulative dose because the former had fewer confounders such as patient body weight.

The dependent variable is a binary indicator of whether a patient will meet any of the five disease activity criteria in 3 or 6 months. The covariates used in the computation of the dependent variable are explicitly removed from the input set of model covariates to ensure that the retrospective analysis mimics how the information could be used prospectively, where the covariates for clinical determination of disease activity in patients in the future are not available. These covariates are qualified relapse count, new T1 Gd+ lesions count, new and enlarging T2 lesions count, EDSS, and whether the patient required switching to another DMT.

ML methods

The overview of our analysis framework is presented in Figure 1. First, the available patient data is split in an 80%–20% ratio using a stratified random sampling strategy for model training and testing, respectively. During the training process, GridSearch cross‐validation was used to select optimal model hyperparameters. We used 10 times repeated 10‐fold cross‐validation to optimize the model hyperparameters (fraction of covariates to subsample at each split = 0.8, tree depth = 6, minimum weight of leaf nodes for regularization = 15, learning rates = 0.05, and number of estimators = 150), which generalizes best the model performance on unseen data. In cross‐validation, the model is trained on a fraction of the training data, and its performance is estimated on the unseen “out of sample” (test) data, and this process is repeated by randomly shuffling the data splits to obtain an averaged performance estimate over several runs. In this way, the training of the classifier and the evaluation of its predictive performance are based on statistically unrelated training and test sets. Therefore, cross‐validation provides a more generalizable estimation of the performance of a model on new unseen data. Finally, the trained model is applied on the completely left out test data set to get an unbiased estimate of model performance on future unseen data. This training and evaluation protocol results in a good balance of preventing overfitting and ensuring generalizability of model performance on future unseen data, also termed the bias‐variance trade‐off in ML.

FIGURE 1

Overview of our analysis framework. The available data is split into a 80–20 fraction using stratified random sampling as training and testing data. The training data are used to select optimal XGBoost model parameters using repeated cross‐validation, and the final model performance is estimated on the completely unseen test data. In the final step, an explainable machine‐learning model SHAP is used to study the covariate contribution to the model predictions and assess covariate importance. gARMSSS, Global Age‐Related Multiple Sclerosis Severity Score; MRI, magnetic resonance imaging; SHAP, Shapley Additive Explanations

We employed the gradient‐boosting, tree‐based ensemble ML algorithm XGBoost that has been empirically demonstrated to be highly effective in a variety of problems. , Final model predictions are obtained by aggregating the predictions of individual decision trees by weighted voting. XGBoost can efficiently handle multicollinearity between covariates as a result of considering each individual covariate independently for the splits during the decision tree building. In fact, because of the ability of ensembles to handle correlated covariates and their interaction well, there is no requirement for covariate selection in the final model. The models have been built with the Python (Version 3.7) interface of the “xgb” library (Version 0.90).

Difficulties in interpreting complex ML models such as XGBoost and their predictions limit the practical applicability of and confidence in ML. In clinical settings, it is critical to understand not only the global aspects of covariate dynamics at a population level but also to understand it specifically for each individual patient. The goal of SHAP is to explain the model prediction for each patient by computing the contribution of each covariate to the prediction using cooperative game theory. The global ranking of important covariates is obtained by ordering the covariates by their mean absolute SHAP values for all patients (exemplified in Figure 3 and Figure 4). The results were obtained with the Python SHAP library (Version 0.35).

FIGURE 3

Covariates predictive of disease activity in patients 6 months in advance. The list of top predictive covariates sorted in decreasing order by their absolute mean SHAP values from the (a) P3‐T‐24 and (b) P4‐T‐24 models. We see that for both models there is a strong overlap of top predictive covariates, including the numbers of weeks of cladribine treatment received, the magnetic resonance imaging measures of new combined unique lesion count and new T1 hypointense lesion count, and other clinically well understood disability measures such as age‐related multiple sclerosis severity score. In absence of laboratory covariates in the P4 model, other well‐known predictive and prognostic covariates become more important such as T1 hypointense lesion volume, age of onset of disease, and time since first symptom. EDSS, Expanded Disability Status Scale; Gd+, gadolinium enhancing; KFSS, Kurtzke Functional Systems Scores; P3‐T‐24, phase III 24 weeks; P4‐T‐12, phase IV 24 weeks; SHAP, Shapley Additive Explanations

FIGURE 4

Dependency plots show the global relationship between top predictive covariates and the output variable for the P3‐T‐24 model. More positive SHAP values push model output toward a more confident prediction of disease activity in patients. For example, with the increasing number of weeks of cladribine treatment received, from 0 (placebo) up to 4 weeks, there is a decrease in model output toward a prediction of no disease activity for the patients. The missing values are imputed to population means only for visualization of these dependency plots and are highlighted with gray circles, noticeably for the new CUA lesion count and the new hypointense lesion counts. It is observed that patients who had a missing value for CUA lesion are most at risk for future disease activity events. It shows that missingness for CUA is not at random and in fact informative and related to the event of interest. CUA, combined unique active; MS, multiple sclerosis; SHAP, Shapley Additive Explanations

RESULTS

First, after combining the patient populations from the three trials, we examined the proportions of patients with disease activity in the three‐treatment arms: placebo, CT3.5 (cumulative cladribine dose of 3.5 mg/kg over 96 weeks), and CT5.25 (cumulative cladribine dose of 5.25 mg/kg over 96 weeks). Looking at Kaplan–Meier curves of time to first disease activity (Figure 2), there is a higher prevalence of disease activity among patients in the placebo arm compared with the two cladribine treatment arms. These results are in accordance with observations in studies demonstrating the efficacy of cladribine versus placebo in reducing clinical relapses, disability progression, and MRI‐assessed disease activity and some aspects of health‐related quality of life. , , , Of our combined population of 1935 patients with more than 6 years of observation, approximately 25% (497 patients) met one or more of the five criteria for disease activity. Second, most patients who had disease activity had it early during the observation period, with 80% of observed events in the first 96 weeks of our more than 6 years of observation.

FIGURE 2

Kaplan–Meier survival curves for disease activity in patients in the combined trial population from ORACLE‐MS, CLARITY, and CLARITY‐Extension. The survival curves are stratified by the treatment arm assignment at the start of the observation period for these three‐armed trials. We see that the disease activity free survival probability in the placebo arm (red) drops lower compared with the two treated arms (CT3.5 in blue and CT5.25 in green), showing that there is higher prevalence of disease activity in the placebo population. Vertical bars represent the time of censoring. CT3.5, cumulative cladribine dose of 3.5 mg/kg over 96 weeks; CT5.25, cumulative cladribine dose of 5.25 mg/kg over 96 weeks

In our analysis of disease activity, we treated all five criteria in the definition equally and assumed no temporal ordering among them. Hence, a patient may meet simultaneously or sequentially multiple criteria during the observation period to qualify as having disease activity, but the time of disease activity is defined as the one corresponding to the first occurrence of any criteria. To quantify the informativeness of each of these criteria, we performed a sensitivity analysis by dropping out one criterion at a time and calculating what percentage of our total number of patients with disease activity fail to be detected with the remaining four criteria.

Table 2 shows the contribution of each of the five criteria toward qualifying disease activity in patients. Dropping C1 criteria does not result in the loss of identification of any patients with disease activity, and in fact the remaining criteria are sufficient to identify this set of patients. This could be explained by the overlap of the relapse count and MRI lesion observations in criteria C2 and C3. On the other hand, dropping the C4 criterion, which is 3‐month sustained EDSS progression, results in missing the detection of approximately 42% of the patients with disease activity. This implies that C1 criteria do not uniquely identify disease activity in patients and that the C4 criterion is extremely important as patients meeting the disease activity criteria by C4 do not seem to exhibit the remaining four criteria at any timepoint. Simultaneously, it highlights the need to include the remaining criteria in our composite analysis objective given the multidimensional aspect of MS, which requires MRI, relapses, and DMT switch observations for qualifying disease activity.

TABLE 2

Percentage of patients with disease activity not detected because of dropping criterion X

No C1: one qualified relapse and one new T1 Gd+ in 48 weeks	0%
No C2: one qualified relapse and two NE T2 in 48 weeks	3.6%
No C3: two qualified relapses in 48 weeks	6.6%
No C4: 3‐month sustained EDSS progression	42.3%
No C5: switching DMT	16.9%

Abbreviations: DMT, disease‐modifying treatment; EDSS, Expanded Disability Status Scale; Gd+, gadolinium enhancing; NE, new and enlarging.

Next, we used a supervised multivariate ML method to train and validate models to predict which patients will have disease activity 3 and 6 months in advance of their clinical determination matching the frequency of patient visits according to the design of the clinical trials and current clinical practice. Specifically, we used XGBoost classifiers and several metrics suitable for studying model performance in cases of unbalanced data (Methods, Table 3, and Appendix S1: Table S1). Overall, we achieved comparable model performance for both 6‐month outcome prediction models P3‐T‐24 and P4‐T‐24 on test data at balanced accuracies of 80% (Table 3). Sensitivity and specificity are, respectively, 84% and 76% and 81% and 78%, whereas the areas under the receiver operating characteristic curve (ROC) for both models are also 80%. In case of the 3‐month outcome prediction models P3‐T‐12 and P4‐T‐12, P3‐T‐12 maintains a similar level of performance with balanced accuracy and area under ROC of 80%, but for P4‐T‐12, the performance reduces to around 75% for these metrics (refer to Appendix S1: Table S1). Balanced accuracy is a metric in binary classification problems to deal with imbalanced data sets and is defined as the average of sensitivity and specificity. Although intuitively it is expected that the 3‐month model should perform better than the 6‐month model, given we are closer to the time of event, there are other factors at play such as the frequency of covariate assessment that can impact the performance of the models.

TABLE 3

Performance estimation of models P3‐T‐24 and P4‐T‐24

		P3‐T‐24		P4‐T‐24
		Train (n = 1356)	Test (n = 340)	Train (n = 1356)	Test (n = 340)
Specificity	TN/(TN + FP)	0.76	0.76	0.77	0.78
Sensitivity	TP/(TP + FN)	0.81	0.84	0.78	0.81
Balanced accuracy	(Sensitivity + Specificity)/2	0.79	0.80	0.78	0.8
AUC‐ROC	Area under curve of ROC	0.79	0.80	0.78	0.8

Note: The table lists the model performance on training and test data with several metrics.

Abbreviations: FP, false positive; FN, false negative; P3‐T‐24, phase III 24 weeks; P4‐T‐12, phase IV 24 weeks; ROC, receiver operating characteristic curve; TP, true positive; TN, true negative.

The top 20 most predictive covariates and their relationship to the predicted output by the models are presented in Figure 3 (P3‐T‐24 and P4‐T‐24) and in Appendix S1: Figure S1 (P3‐T‐12 and P4‐T‐12) in descending order for the models, sorted by their mean absolute SHAP values. This represents a global view at the population level of the important covariates. Covariate importance is understood here as the contribution of covariates in the model to discriminate between the patients who in the future will have disease activity and those who do not. It is interesting to note that there is a strong overlap of the top predictive covariates for all the models including the number of weeks of cladribine treatment received, the new combined unique active (CUA) lesions count (of T1 Gd+ lesions and T2 lesions, without double counting), and the new T1 hypointense lesions count as well as disability measures such as the Age‐Related MS Severity Score (ARMSS). Figure 4 shows the global relationship of the top six predictive covariates with the P3‐T‐24 model output. The model‐derived probability for a patient to have disease activity in the future decreases with the increasing number of weeks of treatment up to 4 weeks, and there is no substantial reduction in model prediction for more than 4 weeks of treatment. Similarly, increasing counts of new CUA lesions and new T1 hypointense lesions raise the models' predicted probability of the patient having disease activity in the future. ARMSS, a clinically well‐understood representation measure of disability ranging from 0 to 10, also shows that a more advanced disease stage relative to a patient's age cohort (i.e., increasing ARMSS values) increases the model‐predicted probability of the patient to have disease activity. The other top predictive factors are laboratory variables such as platelets, creatine kinase, albumin, urate, and so on and prognostic factors such as time since first symptoms and age of onset of disease as well as neurological measures such as Kurtzke Functional Systems Scores Pyramidal Functions.

DISCUSSION

In this study, we developed multivariate ML models to identify predictors of future disease activity in patients with MS 3 and 6 months in advance of clinical determination. It is important to note our observations of disease activity in patients are both right censored as a result of patient dropout as well as interval censored because of missing covariates such as MRI lesions counts, EDSS, or the time of qualified relapses, which are unavailable in the interval between the end of the CLARITY trial and enrollment in the CLARITY‐Extension trial. The dropout of patients from the CLARITY and CLARITY‐Extension trials (approximately 90% of enrolled patients completed all three arms) were previously investigated, and the efficacy results with the treatment switching were found to be not substantially biased by informative dropout.

The ML models achieve comparable levels of performance with about 80% balanced accuracy on both the training and testing data sets, indicating no obvious model overfitting. XGBoost deals with a large number of input covariates while efficiently handling correlated covariates and missing values without an impact on performance. Given the imbalance in the data where disease activity was observed in only about 25% of the patient population, in addition to using a boosting algorithm, we also used majority class undersampling as well as a cost‐sensitive objective function. For example, mispredicting that a patient will have disease activity incurs a 10% additional penalty compared with mispredicting that a patient will not have disease activity, resulting in higher sensitivity than specificity of the models. Such ML models in conjunction with explainability methods (SHAP) allow us to gain insights on how the most predictive covariates can discriminate between patients who will have disease activity and those who will not.

Among the top predictive covariates, the number of weeks of cladribine treatment received has the largest impact. It is not surprising given the demonstrated effectiveness of the drug in patients experiencing a first demyelinating event and patients with RRMS with four weekly treatments for a sustained effect observed over 4 years. , , , We found that the other top predictive covariates for early detection of disease activity in patients are higher count of new CUA lesions, higher count of new T1 hypointense (black hole) lesions, and clinical disability measures such as higher ARMSSS in all our models. These results are explained by the clinical understanding of the onset of action of cladribine, which shows a reduction in CUA lesion count in treated patients. In the P4 models, which excluded the laboratory covariates not routinely collected in clinical practice, we observed other well‐known markers such as age of onset of disease, time since first symptom, volumes of T2 lesions, and T1 hypointense lesions become more predictive , in the absence of the laboratory covariates that have not been studied in the MS clinical literature. The role of laboratory covariates such as urate, a predictor of disease activity in our P3‐T‐12 model (Appendix S1: Figure S1), in MS disease has been previously investigated, , , , , but its precise understanding is currently lacking. Further investigations should focus on the assessment of the prognostic (treatment independent) versus predictive (treatment dependent) nature of these factors. It is interesting to note the strong overlap of predictive covariates among all our models, which gives further evidence to their robustness.

Such a predictive model will aid in treatment planning, for example, in deciding the frequency of patient visits for at‐risk patients, or prompt other clinical interventions. Although our analysis population is different from the label population for cladribine in several countries where it has been approved, this high‐dimensional data set enables us to present a feasibility study on how data‐driven ML models could help in the future as a clinical decision support tool by presenting an interpretable prediction about how clinically available covariates drive the probability of future disease status in patients with MS. In this regard, it is also interesting to note the updates to clinically relevant covariates attributed to the evolving understanding of MS. For example, recent trials include newer neurological measures of disability progression such as a timed 25‐m walk and the nine hole‐peg test that have been found to be robust and stable compared with traditional clinical trial end points, such as 3‐month EDSS progression.

Prior work has investigated various clinical end points in patient populations with MS using ML techniques to access covariate importance. However, these studies used the in‐built covariate importance of various tree‐based ensemble methods that have been known in the literature to carry various biases. Our work differs in two key aspects. First, the analysis objective that serves as the dependent output variable for our ML models encompasses multiple clinical end points, which is important as demonstrated by the obtained results to have a more complete picture of disease activation and/or progression. Second, we use the state‐of‐the‐art explainable ML technique SHAP, , , which overcomes the known limitations of the in‐built covariate importance assessment methods of the ensemble decision tree methods.

Our findings further highlight the importance of quality MRI evidence in the management of MS disease progression and patient monitoring as highlighted in the MRI in MS guidelines for MS. Although such quality and detailed information regarding counts of various new lesion types are difficult to collect in routine clinical practice and are usually only available in research or during drug development, they are important predictive covariates for MS disease activity in patients. The current model treats the task of predicting the future probability of disease activity as a binary classification. Dynamic prediction models such as recurrent neural networks taking as input longitudinal covariates and predicting an updated probability of the risk of disease activity could be explored in future work.

CONCLUSION

In summary, we conducted an analysis on integrated data from multiple clinical trials at various stages of MS (patients experiencing a first demyelinating event and patients with RRMS) to investigate the predictive covariates of onset of disease activity. Disease activity, as defined by the five criteria, was found to be more frequent among the placebo population, and 3‐month sustained EDSS progression was the most informative among the five criteria. The novelty of our work lies in the use of explainable ML methods for training multivariate predictive models combining patient baseline characteristics as well as longitudinal MRI readouts and neurological and laboratory measures to identify patients 3 and 6 months before their clinical observation of disease activity. We trained and validated multivariate nonlinear models with the XgBoost algorithm at 80% balanced accuracy and area under ROC, with different subsets of input covariates available during phase III and phase IV trials to predict patient outcomes of disease activity 3 and 6 months in advance. Predicted disease activity probability in patients was found to decrease with cladribine treatment duration increase, the most important predictor. Additional top predictive covariates from explainable ML models were the count of new CUA lesions, new T1 hypointense lesions, and ARMSS score. Such multivariate predictive models to identify patients early who in the future will have disease activity can improve the available knowledge of the underlying mechanism of disease activity in patients with MS to enable better patient monitoring and treatment planning.

CONFLICT OF INTEREST

S.B., A.M., P.G., and N.T. are employees of Merck Serono S.A., Lausanne, Switzerland, an affiliate of Merck KGaA, Darmstadt, Germany. A.‐F.B.‐A. and S.R. are employees of Ares Trading SA, Eysins, Switzerland, an affiliate of Merck KGaA, Darmstadt, Germany.

AUTHOR CONTRIBUTIONS

N.T., A.M., and A.‐F.B.‐A. designed the research. N.T., S.R., P.G., and S.B. performed the research. S.B. analyzed the data. S.B., A.M., A.‐F.B.‐A., S.R., P.G., and N.T. wrote the manuscript.

Appendix S1

Click here for additional data file.

Appendix S2

Click here for additional data file.