Machine learning for classification and quantification of monoclonal antibody preparations for cancer therapy

TL;DR

This study develops a rapid, noninvasive Raman spectroscopy-based machine learning method for quality control of monoclonal antibody preparations in cancer therapy, achieving high accuracy in classification and quantification.

Contribution

It introduces a collaborative data analytics workflow combining Raman spectroscopy and machine learning for antibody quality assessment, significantly reducing prediction errors.

Findings

Misclassification rate of 0.8%

Mean error rate of 4%

Effective collaboration platform used for model optimization

Abstract

Monoclonal antibodies constitute one of the most important strategies to treat patients suffering from cancers such as hematological malignancies and solid tumors. In order to guarantee the quality of those preparations prepared at hospital, quality control has to be developed. The aim of this study was to explore a noninvasive, nondestructive, and rapid analytical method to ensure the quality of the final preparation without causing any delay in the process. We analyzed four mAbs (Inlfiximab, Bevacizumab, Ramucirumab and Rituximab) diluted at therapeutic concentration in chloride sodium 0.9% using Raman spectroscopy. To reduce the prediction errors obtained with traditional chemometric data analysis, we explored a data-driven approach using statistical machine learning methods where preprocessing and predictive models are jointly optimized. We prepared a data analytics workflow and submitted the problem to a collaborative data challenge platform called Rapid Analytics and Model Prototyping (RAMP). This allowed to use solutions from about 300 data scientists during five days of collaborative work. The prediction of the four mAbs samples was considerably improved with a misclassification rate and the mean error rate of 0.8% and 4%, respectively.

Tables

Table 1. Table 1: Misclassification rates and concentration MAREs obtained by the traditional chemometric (linear) baseline method (Section  2.3.1 ) and by the crowdsourced machine learning ensemble (Section  2.3.2 ) on 429 test points.

molecule	miscl. rate (linear)	MARE (linear)	miscl. rate (ML)	MARE (ML)
Infliximab	$13.7\%$	$12.3\%$	$0.9\%$	$8.4\%$
Bevacizumab	$19.8\%$	$14.0\%$	$1.0\%$	$4.3\%$
Ramucirumab	$9.0\%$	$7.3\%$	$0.0\%$	$3.5\%$
Rituximab	$16.0\%$	$26.7\%$	$1.0\%$	$6.9\%$
Overall	$14.5\%$	$14.7\%$	$0.7\%$	$5.8\%$

Equations

$$\hat{\ell }={\mathrm{argmax}}_{j\in \{\text{A, B, Q, R}\}}{p}^j.$$

$$\widehat{y}=f^{\text{reg}}\big{(}{\bf x}^{\text{reg}}\big{)}\enspace.$$

$$R^{\text{clf}}(\mathcal{M},\mathcal{D})=\frac{1}{n}\sum _{i=1}^n\mathbb{I}\left\{{\ell }_i\ne {\hat{\ell }}_i\right\},$$

$$R^{\text{reg}}(\mathcal{M},\mathcal{D})=\frac{1}{n}\sum _{i=1}^n\frac{|y_i-{\hat{y}}_i|}{y_i}.$$

$$R^{\text{comb}}(\mathcal{M},\mathcal{D})=\frac{2}{3}{R}^{\text{clf}}(\mathcal{M},\mathcal{D})+\frac{1}{3}{R}^{\text{reg}}(\mathcal{M},\mathcal{D}).$$

Full text

1. 1 Introduction

2. 2 Material and methods

3. 2.1 mAbs Samples

4. 2.2 Raman spectroscopy

5. 2.3 Data analysis

6. 2.3.1 The traditional linear approach

7. 2.3.2 A crowdsourced machine learning approach

8. 3 Results

9. 4 Discussion

10. 5 Conclusion

1 Introduction

Cancer can be treated using numerous strategies, such as surgery, radiation therapy, immunotherapy, hormone therapy, stem cell transplant, chemotherapy, and, recently, targeted therapy. Targeted therapy, for example, using monoclonal antibodies (mAbs), is the foundation of precision medicine. At present, it constitutes one of the most important strategies to treat patients suffering from cancers such as hematological malignancies and solid tumors. Monoclonal antibodies are proteins which bind to specific substance on cancer cells and act by immune-mediated cell killing mechanisms. They may therefore help the immune system to destroy cancer cells, stop cancer cells from growing, stop signals that help to form blood vessels, deliver cell-killing substances to cancer cells, and cause cancer cell death. As classic chemotherapy drugs, mAbs are commonly used alone or with other cytotoxic drugs or radioactive substances to kill cancer cells.

Those drugs designed for parenteral administration are aseptically prepared by pharmacists. Drugs are extemporarily reconstituted or diluted in 5% glucose or 0.9% sodium chloride to obtain a sterile final preparation at the dose prescribed by the physician. Even if pharmaceutical regulation does not require final characterization of each compounding drugs, numerous pharmacists have nevertheless implemented analytical control. Those controls have to be discriminant to ensure the nature of the drug in spite of similar physicochemical and spectral properties, sensitive to guarantee the dose even for low concentrations and fast to secure the medication process without delayed drug delivery. Consequently, different analytical strategies using direct flow injection analysis, high performance liquid chromatography with UV detection [Delmas et al., 2009], or vibrational molecular spectroscopies such as Raman and infrared caudron [Lê et al., 2016, 2014, Bourget et al., 2014] have been investigated to control cytotoxic compounding drugs. Despite their importance in the cancer therapy, only few studies reported on analytical methods to control mAbs preparations [Bazin et al., 2010, Jaccoulet et al., 2016b, a, 2015, Paul et al., 2012, Ashton et al., 2013]. All of these techniques were invasive and required a sampling of the preparation for analysis. Because of the inherent toxicity of cytotoxic drugs as cancerogenic, mutagenic and teratogenic properties, these drugs present a risk of occupational exposure for healthcare workers. Thus, we decided to explore the feasibility of a noninvasive, nondestructive, and rapid analytical method as Raman spectroscopy to ensure the quality of the mAbs preparations produced at hospital. This molecular vibrational spectroscopy based on inelastic scattering of monochromatic light, usually from a laser in the visible, near infrared, or near ultraviolet range, was successfully investigated for chemotherapy drug control [Lê et al., 2016, 2014, Amin et al., 2014, Bourget et al., 2014].

Due to the complexity of Raman spectral data, multivariate analysis has to be used to extract pertinent information. We have tried traditional linear techniques (partial least squares discriminant analysis and regression) and we obtained unsatisfactory results at both the molecule classification and concentration regression tasks. In this paper we report results using machine learning methods. We have built an analytics workflow (Figure 1), and submitted the data to a collaborative data challenge platform called Rapid Analytics and Model Prototyping (RAMP; http://ramp.studio), developed by the Paris-Saclay Center of Data Science. A large number of feature extraction and prediction techniques were submitted by about 500 data scientists. The top results, reported in this paper, were obtained after about five days of collaborative work.

2 Material and methods

2.1 mAbs Samples

We evaluated four mAbs: Bevacizumab (Avastin${}^{\tiny{\text{\textregistered}}}$ 25 mg/mL from Roche), Infliximab (Remicade${}^{\tiny{\text{\textregistered}}}$ 100 mg from Schering-Plough), Ramucirumab (Cyramza${}^{\tiny{\text{\textregistered}}}$ 10 mg/mL from Lilly), and Rituximab (Mabthera${}^{\tiny{\text{\textregistered}}}$ 10 mg/mL from Roche). We prepared all drugs separately in aseptic condition and analyzed them after dilution in 0.9% chloride sodium at various concentrations that cover the therapeutic range currently prepared to treat patients (10 concentrations for each drug: Bevacizumab from 0.5 to 25 mg/mL, Infliximab from 0.3 to 10 mg/mL, Ramucirumab from 1 to 10 mg/mL and Rituximab from 0.4 to 10 mg/mL). We prepared 12 independent series (12 batch of 0.9% chloride sodium) for each drug. We conditioned all compounded solutions into 3 glass vials (Interchim${}^{\tiny{\text{\textregistered}}}$, Montluçon, France), conserved at +4 ∘C, and analyzed in accordance to laboratory requirements.

Throughout the paper, we will use the following labels for the four drugs: A= Infliximab; B=Bevacizumab; Q=Ramucirumab, R=Rituximab.

2.2 Raman spectroscopy

We performed Raman spectral acquisitions with a Labram HR evolution microspectrometer (Horiba Jobin Yvon, Lille, France). The excitation source was a 633 nm single-mode diode laser (Toptica Photonics, Germany) generating 35 mW on the sample. The microspectrometer was equipped with an Olympus microscope and measurements were recorded using 10 and 100X MPlan objectives (Olympus, Japan). We collected light scattered by the sample through the same objective. We intercepted Rayleigh elastic scattering by an edge filter. A Peltier cooled (-70 ∘C) multichannel CCD detector (Coupled Charge Device; 1024 $\times$ 256 pixels) detected the Raman Stokes signal dispersed with a 100 $\mu$m slit width and 600 grooves/mm holographic grating. The spectral resolution from the full width at half maximum of the silica wafer band at 521 cm*-1* was 2 cm*-1*. We studied the spectral region of 400 - 4000 cm*-1*. The acquisition time of each spectrum was 2 $\times$ 15 s per collected spectrum to have better sensitivity. We conducted spectral acquisition and data pre-processing with Labspec6 software (Horiba Jobin Yvon SAS, Lille, France). Because the sample compartment was not adapted to analyze vials, we used a homemade vial adaptor to centralize the vial and to secure the position of the sample on the base plate. In order to correct spectral variation due to focalization, we pre-treated all Raman spectra using LabSpec6 by normalization based on total area. The result of this collection step was a spectrum containing $1866$ values.

The raw spectra are shown in Figure 2 for all concentrations and all molecules. Raman spectra of mAbs, similarly to that of other proteins, are very difficult to interpret. Inter and intra molecular effects including peptide-bond angles and hydrogen-bonding patterns may influence Raman band positions. Structural information can nevertheless be deduced from Raman vibrational bands as amide I and amide III bands at 1650 and 1300 cm*-1*, respectively.

2.3 Data analysis

In this section we formalize the data analysis problem and report overall results using traditional linear techniques and statistical machine learning approaches. We will present a more detailed analysis using the optimal technique in Section 3.

The goal of the data analysis is to predict the drug type $\ell$ and its concentration $y$ based on the recorded Raman spectrum ${\bf s}\in{\mathbb{R}}^{1866}$. We used the following experimental protocol to develop and tune the analysis techniques and to test and compare the different methods. We collected a total of 360 spectra for each of the four mAbs agents, except for Ramucirumab (348 spectra). Following the standard machine learning protocol that consists in evaluating model performance on a held out set of data, the 1428 measurements were randomly split into two data sets. The 429 test instances were used to evaluate the developed techniques. All results reported in this paper are coming from this test set, but they were hidden from the analysts throughout development. The remaining 999 points were further split into eight random cross-validation folds. Each training and validation set contained 799 and 200 points, respectively. Each technique was trained eight times on the training sets. The eight trained models were then evaluated on each validation and test point. The final prediction is the average of these eight predictions on each validation and test point. This technique, commonly known as cross-validation bootstrap aggregation (CV bagging), takes advantage of the variance reduction property of averaging and allows us to develop robust predictors despite the small data size. To further improve the results, we blended the best models using the technique of Caruana et al. [2004]. In a nutshell, we take the pointwise mean prediction of the best models selected in a greedy loop until the validation result stops improving.

To allow the analysts to focus on various steps of the data analysis, we split the workflow into four modules (Figure 1). The first feature extractor $g^{\text{clf}}$ converts the raw spectra ${\bf s}$ into a fixed size feature vector ${\bf x}^{\text{clf}}=g^{\text{clf}}({\bf s})$. This is then followed by a classifier $f^{\text{clf}}$ that outputs a vector {\bf p}=(p^{\text{A}},p^{\text{B}},p^{\text{Q}},p^{\text{R}})=f^{\text{clf}}\big{(}{\bf x}^{\text{clf}}\big{)}, indexed by the four drug types, representing the estimated probability that the spectrum belongs to each of the four types. The predicted class is then

[TABLE]

The second feature extractor $g^{\text{reg}}$ also receives the raw spectra ${\bf s}$ but also the probability vector ${\bf c}=f^{\text{clf}}(g^{\text{clf}}({\bf s}))$, and converts them into a fixed size feature vector ${\bf x}^{\text{reg}}=g^{\text{reg}}({\bf s},{\bf c})$. The rationale is to let the regression algorithm know the (estimated) identity of the molecule whose concentration it has to predict. The final module is a regression model $f^{\text{reg}}$ that takes this second feature vector ${\bf x}^{\text{reg}}$ and estimates the drug concentration

[TABLE]

Given a validation or test set {\cal D}=\big{\{}({\bf s}_{i},\ell_{i},y_{i})\big{\}}_{i=1}^{n} containing triplets of spectra, molecule types, and concentrations, we evaluate the performance of the model {\cal M}=\big{(}g^{\text{clf}},f^{\text{clf}},g^{\text{reg}},f^{\text{reg}}\big{)} as follows. We first compute the classification error

[TABLE]

where the indicator function ${{\mathbb{I}}\left\{{X}\right\}}$ is one if its argument $X$ is true, and zero otherwise. To assess the concentration predictor, we use the mean absolute relative error (MARE)

[TABLE]

The analysts were instructed to minimize a combined error

[TABLE]

The coefficients are reflecting the more stringent requirements posed on the classification task.

2.3.1 The traditional linear approach

To establish baseline performance, we applied linear discrimination and regression methods using the Matlab software R2011a, common in traditional chemometry. That is to say, both the preprocessing $g^{\text{clf}}$, $g^{\text{reg}}$ and the predictor functions $f^{\text{clf}}$, $f^{\text{reg}}$ were linear. Classification was performed using partial least squares discriminant analysis (PLS-DA) and principal component analysis discriminant analysis (PCA-DA). We used PLS regression for estimating the concentration. Both for classification and regression, we used various spectral processing techniques to limit the non-informative spectral background: baseline correction with a $n$th degree polynomial curve, first and second Savitsky-Golay derivatives, standard normal variate (SNV), and combined preprocessing. For each model, the optimal number of latent variables was determined by 10-fold cross validation.

The overall misclassification rate of the baseline method was $14.5\%$ while its MARE was $14.7\%$. More detailed analysis of the performance can be found in Section 3.

2.3.2 A crowdsourced machine learning approach

Machine learning is essentially a trial-and-error-based science. The scikit-learn library [Pedregosa et al., 2011], which most of our analysts used, contains dozens of different regression and classification techniques. Selecting and tuning the best predictor takes experience and a lot of experiments. The choice of feature extractors and filters is even larger. The Paris-Saclay Center for Data Science has developed a unique collaborative studio which allows a large number of data scientists to collaborate on various scientific and industrial data analytics workflows. Participants have access to a brief description, the public training data set, and a first (non-optimized) solution in a starting kit. They study the problem, develop and tune solutions, and submit them through a web interface (Figure 3(a)). The models are then trained and evaluated, and the validation performance score is fed back to the participants through a public leaderboard (Figure 3(b)). The test score obtained on the held out test set, reported in this paper, remains hidden from the participants to avoid overfitting.

Unlike most of the data challenge sites that expect solutions to analytics problems, we ask the data scientists to submit code. This gives us more flexibility at the evaluation, generates a working prototype, and, most importantly, lets participants to collaborate by looking at each other’s code and combining the different ideas. The code is publicly available at the website (after free sign-up), making the results fully reproducible by anyone. Moreover, any new method can be easily tested by submitting it on the site, so our results are fully transparent.

The drug spectra challenge was submitted to more than three hundred students in three different classroom RAMPs at Université Paris-Saclay, Ecole Polytechnique, and Mines ParisTech. Each RAMP involved about five days of work. In each class, we started with a closed period when students did not see each other’s code, only their scores. This allowed us to grade the students individually but also to make them explore the space of solutions independently. This was followed by an open collaborative period. On the other hand, the classes could not see the solutions developed by the other classes. This, interestingly, led to different collective strategies in solving the problem. We blended the models in each RAMP separately which made it possible to compare the collaborative score of each group.

The overall winner of the contest was the Mines RAMP with a classification error of $0.7\%$ (only three low-concentration test samples were missed!) and a concentration MARE of $5.8\%$. The dominant solution that the group converged to was a log transformation of the spectra (but no other smoothing or preprocessing), followed by either factor analysis or principal component analysis to extract 10 features, followed by a small neural network. Interestingly, the same pipeline worked for both classification and regression, with, of course, different parameterizations of the neural net. In the regression step, all top models used the same strategy: to learn a different regression model (although parametrized the same way) for each molecule.

The Saclay RAMP (chronologically the first) achieved a classification error of $1.6\%$ and a concentration MARE of $4.9\%$. This group did more preprocessing. Smoothing of the spectra was done with linear filters such as the Savitzky–Golay filter [Schafer, 2011] that allows to preserve the amplitudes of the peaks in the signal. Certain models used simpler strategies based on convolution with Hanning windows. Following the smoothing step, the baseline correction was done by subtracting a polynomial least square fit of the data or simply by removing a constant. The polynomial order was 0 (constant) or 1 (linear). The order 0 corresponds to the subtraction of the mean of the spectra. Some of the solutions proposed to subtract the median of the spectra. As the machine learning models employed in the second step are sensitive to the scale of the data, the spectra, for some models, were normalized so that each spectrum after smoothing and baseline correction was of norm one. For prediction, the majority of winning solutions employed non-linear kernel-based techniques: kernel PCA for non-linear dimensionality reduction [Schölkopf et al., 1998] and support vector machines (SVM) for prediction [Cristianini and Shawe-Taylor, 2000]. The best proposed approaches used Gaussian kernels and polynomial kernels of order up to 4. The low MARE was achieved by exploiting the fact that concentrations were discretized (see Section 4).

Finally, the Polytechnique RAMP achieved a classification error of $2.1\%$ and a concentration MARE of $12.2\%$. The reason of the suboptimal performance was that the group was only exploring forest-based regression models (extra trees [Geurts et al., 2006], random forests [Breiman, 2001], gradient boosting [Friedman, 2000]). These otherwise popular and usually well-performing nonparametric classifiers seemed to be a suboptimal choice for the functional data we had in this problem.

Note again that all the models and code are available freely at the RAMP site, so anybody can consult their detailed parameterizations, reuse them, and reproduce the results. Furthermore, anybody can resubmit competing solutions which will be tested using the same protocol as during the RAMPs, on the same training and test sets.

3 Results

Table 1 shows that the large improvement obtained by machine learning techniques is uniform across all molecules. The three misclassified molecules (that correspond to the overall misclassification rate of $0.7\%$) all came from samples with low concentration ($y\leq 0.6\,\mathrm{mg/ml}$). Using the linear method, $63\%$ of the samples below concentration $y\leq 1\,\mathrm{mg/ml}$ had larger relative error than $15\%$ whereas the figure was only $5.6\%$ in the case of the ML method.

Figure 5 shows the estimation bias (the mean of the relative error $(y-\widehat{y})/y$) and the resolution (the standard deviation of the relative error) as a function of the concentration. Besides a downward slope explained by regression to the mean, the bias (Figure 5(a)) is negligible above 0.1 mg/ml. The relative resolution (Figure 5(b)) drops sharply with the concentration (not surprisingly) from 15% (below $y<1\,\mathrm{mg/ml}$) to 5% (above $y>3\,\mathrm{mg/ml}$).

4 Discussion

Over the last decade, the interest of statistical machine learning for biological applications has grown considerably. As part of both computer science and statistics, machine learning (ML), as a scientific discipline, aims to develop algorithms that can make sense of data. A typical outcome of machine learning is a model that can make prediction from new data, after having learnt from many training examples. Based on the hypothesis that spectra are influenced by the nature and the concentration of mAbs in solution, we decided to explore how machine learning can help to resolve the pharmaceutical challenge of classification and concentration estimation for pharmaceutical drugs particularly difficult to discriminate by traditional analytical methods.

We found that combining Raman spectroscopy with machine learning methods present an interesting potential to secure the medication process by the identification and the quantification of chemotherapy preparations.

As opposed to traditional analysis using CLHP/UV or LC/MS/MS, Raman spectroscopies, such as near infrared spectroscopy, allow direct measurement through the glass and plastic packaging Broad et al. [2000], Vergote et al. [2002]. Because of rapid analysis and non-destructive and non-invasive measurements, those spectroscopies are largely used for the Process Analytical Technology De et al. [2011] with at-line and on-line measurements to control primary and secondary manufacturing process in the pharmaceutical industry.

Despite promising results to ensure the nature and the dose of the drug in solution by classification and regression analysis respectively, machine learning has some particular limitations and pitfalls that should be avoided. Machine learning algorithms usually require more training data to train than linear models. In addition, one has to be careful to collect data under the same distribution as when the predictors are used in practice. First, if the training concentrations are over- or under-sampled in certain regions, the nonparametric concentration estimate may be biased towards the oversampled values. This is due to the fact that machine learning models try to minimize the errors on average. Second, if only certain concentration levels are used for certain molecules, the classifier can learn these values and use the concentration information for classification. Third, if the concentration levels are discretized, the regression model can learn those discretization levels, and quantize the continuous estimate, creating a bias in case true concentrations fall between those levels. In fact, one of the student groups boosted the prediction accuracy by forcing the model to predict only concentrations present in the training set. More concretely, they rounded the predicted concentration by a regression model to the closest discrete concentration. In our study, the predictors were so precise that these effects were negligible, but in more complicated measurements, with larger uncertainties, they have to be taken into consideration. One practical way to consider in future data collection protocols in order to prevent performance gains by rounding is to randomize concentration values using a sampling distribution that matches realistic scenarios.

5 Conclusion

Due to similar structures and low concentration values, discrimination and quantification of mAbs preparations posed difficult challenges. In addition to the type of analytical technique, this study have highlighted the importance of the data analysis method and more particularly the power of statistical machine learning associated to Raman spectroscopy to ensure the chemical quality of medications produced at hospital before patient administration.

Acknowledgement

We would like to thank the 500 students and researchers participating in the initial hackaton and in the three classroom RAMPs. Their individual performance and credits can be assessed through their combined scores and contributivities at the three leaderboards (here, here, and here). We are also grateful for the Paris-Saclay Center for Data Science for supporting this research.