In this study, we investigated gene expression features associated with the responses to chemotherapy. The gene expression profiles were extracted from the datasets summarized in Table 1. The clinical outcome information was related to response on different chemotherapy regimens, linked with high throughput gene expression profiles for the individual patients.

Each patient was primarily labeled as either responder or non-responder to a treatment. For all the datasets taken from the GEO repository, we used the response criteria formulated in the respective original papers first publishing these data. Namely, for two breast cancer datasets, GSE25066 (Hatzis et al., 2011; Itoh et al., 2014) and GSE41998 (Horak et al., 2013), we considered partial responders as responders. For the first multiple myeloma dataset, GSE9782 (Mulligan et al., 2007), we took the (non)responder classification used by the authors, where patents with complete and partial response were annotated as responders, and with no change and progressive disease – as non-responders. For three other multiple myeloma datasets, GSE39753 (Chauhan et al., 2012), GSE68871 (Terragna et al., 2016), and GSE55145 (Amin et al., 2014), we considered complete, near-complete and very good partial responders as responders, whereas partial, minor and worse responders – as non-responders. For the datasets of pediatric Wilms kidney tumor, ALL and AML, extracted from the TARGET gene expression repository of National Cancer Institute (Goldman et al., 2015), the cases was classified according the distribution of the event-free survival time, which appeared to have two modes with different slopes (Supplementary Figure S1).

To preclude any possible bias that may affect the performance of machine-learning classifiers due to unequal representation of samples in two different classes (clinical responders and non-responders), numbers of responding and non-responding cases were equalized within each dataset. Equalization was done by taking the full smaller subset of those for the two classes (responders/non-responders), and then by random selection of samples from the bigger subset. Thus, each resulting dataset contained equal numbers of cases classified as responders and non-responders.

To engineer a plausible feature space, where the SVM can be applied efficiently, we proposed to select from tens of thousands of individual gene expression features only few of them, which produce a good separation of clinical responders from non-responders. To do so, for every dataset under investigation we selected its particular top 30 genes, whose expression levels taken one by one had the highest ROC AUC values for distinguishing responder and non-responder profiles. We made a number of top informative features equal to 30 because the usual number of samples in considered datasets was not lower than 50 (a direct heuristic number for degree of freedom). These 30 top marker genes, and response statuses (100 for a responder, 0 for a non-responder) for all selected patients from all datasets are listed on Supplementary Table S1.

To produce more robust feature selection, for each dataset having, say, N samples, the leave-one-out procedure has been performed. Each individual sample was removed from the investigation one at a time, so N subdatasets each having N-1 individuals were generated. For each subdataset, the ROC AUC test was performed between responders and non-responders for each gene. The genes were next sorted according to their ROC AUC, and top 30 were selected for each subdataset. The final list of such core informative genes was generated as an intersection between top 30 selected genes for all N subdatasets. For every dataset under investigation, these final core sets are listed in Supplementary Table S2; the number of core marker genes is also shown on Table 1.

We developed a first of its class data trimming³ tool termed FloWPS that has a potential to improve the performance of machine learning methods. Since extrapolation is a widely recognized Achilles heel of SVM (Arimoto et al., 2005; Balabin and Lomakina, 2011; Balabin and Smirnov, 2012; Betrie et al., 2013), FloWPS avoids it by using the rectangular projections along all irrelevant expression features that cause extrapolation during the SVM-based predictions for every validation point.

In this section we describe and investigate our data trimming procedure (FloWPS) as a preprocessing for SVM application.

Since the number of samples in most of the datasets used here was relatively low, we tested our classifier using the leave-one-out cross-validation method, which introduces lesser errors than the standard five-bin cross-validation scheme generally applied for bigger datasets. According to the leave-one-out approach, for each sample i = 1, N serves as a validation case whose response to the treatment had to be predicted, whereas all remaining samples, j = 1,…(i-1),(i+1),…,N, collectively acts as a training dataset, and this procedure is repeated for all the samples. For machine leaning without data trimming, in a predefined feature space F = (f₁,…, f_s ) every sample i, given for the test, is assigned by a classifier, constructed to (N-1) samples used for training.

According to the current data trimming approach, instead a fixed space F for all N testing samples, we propose using an individual space F_i, which contains individually adapted training data (of N-1 samples) for the testing sample i. It can be implemented using the following heuristics (Figure 1).

(1) From the whole predefined feature space F = (f₁,…, f_s ) we extract a subset F_i (m), where m is a parameter. A feature f_j is kept in F_i(m) if on its axis there are at least m projections from training samples, which are larger than f_j (i), and, at the same time, at least m, which are smaller than f_j (i). The procedure for extraction of a subset F_i(m) is illustrated in Figure 1A for a two-dimensional space F = (f₁, f₂). A violet point stands for the validation sample in the feature space. Turquoise dots represent scattering of the training points. For example, the m-condition for the feature f₂ is satisfied when m = 0,1,2,3,4,5 (projection of the training set on f₂ axis has five points both below and above the validation point), whereas for the feature f₁ it is satisfied only for m = 0 (projection of the validation point on axis f₁ lies outside of the cloud of training points).

(2) In F_i (m) we keep for training only k closest samples (from given (N-1) samples); k is also a parameter (Figure 1B; note that although for the sake of simplicity k = 4 in the picture, in the computational trials we varied k from 20 to N-1).

Hence, for every individual i = 1, N, and m and k parameter values, the predicted classification values are obtained [i.e., predictions P_i (m,k), i = 1, N]. Considering known response status for each sample i, it is possible to calculate AUC values for a whole set of samples as a function over whole range of the parameters m and k (Figure 2B). Since the predicted classification efficiencies depend upon the chosen values for m and k, it is possible to interrogate the AUC values over the full lattice of all possible (m, k) pairs.

We propose an algorithm of achieving the optimal (m,k)-settings for a final classifier (Figure 2A). The AUC threshold (θ) is set to θ = p ⋅ max(AUC), where max(AUC) is the maximal value of AUC, taken over the set of all possible (m, k) pairs, and the parameter p equals to a user-defined confidence threshold. To illustrate performance of this approach, we took two alternative values of p = 0.95 or 0.90, and then considered all the (m,k) pair positions on the AUC(m,k) topogram. We next screened for the positions where AUC exceeded the threshold θ, and the total combination of these positions was taken as the prediction-accountable set S (Figure 2B; prediction-accountable positions are shown in yellow). The final prediction of FloWPS (P_F) for a certain validation case should be calculated by averaging the SVM predictions, P(m,k), over the whole set of positions belonging to the prediction-accountable set S, according to the formula: P_F = mean_S(P(m,k)).

The usual SVM method, i.e., without FloWPS data trimming, corresponds to a very right and bottom corner of the AUC(m,k) topogram (Figure 2B), with the parameter settings m = 0, k = N - 1. On the example shown in Figure 2B, the classical SVM, without any doubt, provides essentially lower accuracy than FloWPS.

At first, we investigated performance of FloWPS on ten cancer gene expression datasets (Table 1) with the default SVM settings (linear kernel and cost/penalty parameter C = 1). During our calculations, the FloWPS classifier was first fitted for the training dataset without a sample (say, i) to be classified. For these all (N-1) samples AUC_i(m,k) was calculated as a function of data trimming parameters m and k (see Figure 2A). This enabled finding the prediction-accountable set S_i in the AUC_i(m,k) topogram (on Figure 2B, the set was marked with bright yellow). The m and k values from the set S_i were then used for data trimming and classifying of a single sample i. In parallel, we applied the standard SVM algorithm for leave-one-out cross-validation without data trimming, i.e., m = 0, k = N-1 for each training sub-dataset. The comparison is shown on Table 2, Supplementary Table S3, and Figures 3, 4.

The discrimination threshold (τ), which is shown as a black horizontal line on Figure 3 (so that any sample with FloWPS prediction value above τ is classified as a responder, and below it – as a non-responder), was set to minimize the sum of FP and FN predictions.

For every dataset, confidence parameter p and scheme of gene selection, FloWPS classifier demonstrated the ROC AUC exceeding the corresponding value for the classical SVM (Table 2). For three datasets out of ten, AUC for classical SVM was between 0.64 and 0.68. For all these cases, application of FloWPS with confidence level p = 0.90 enabled obtaining essentially better AUC values ranging between 0.71 and 0.78.

The comparison of classifier’s quality by another metric, the FDR⁴, has demonstrated similar results: FDR was lower for FloWPS than for classical SVM for almost all the cases (Table 2, columns without boldface font). Other metrics, such as sensitivity (Sn), specificity (Sp), accuracy rate (ACC) and MCC⁵ also strongly tend to be higher for FloWPS than for classical SVM without data trimming (Supplementary Table S3).

Although the classifier quality tended to be higher for data trimming than for default SVM settings, the advantages were different in different cancer datasets. The FloWPS performance, therefore, was investigated for different SVM kernels (linear vs. polynomial) and different values for cost/penalty parameters C (ranged from 0.1 to 1000), Figure 5 and Supplementary Table S4. These calculations were done for the core marker gene datasets and FloWPS confidence parameter p = 0.90. The advantage of FloWPS over SVM is more essential in the conditions vulnerable to SVM overtraining, e.g., for linear kernel with high values of the cost/penalty parameter (C = 100 or 1000) or for polynomial kernel, where SVM may be easily overfitted. Fortunately, FloWPS precludes such overfitting, thus raising AUC and decreasing FDR. The same pattern was also seen for the Sn, Sp, ACC and MCC values (Supplementary Table S4).

Note that FloWPS is not the only possible data reduction/feature selection method, which may be used for preprocessing to improve the classifier’s quality. To try a simple alternative to FloWPS, which is, however, not specific to individual samples, we did calculations based on PCA mode rather than original features. The number of PCs taken for building the SVM model, may act as a parameter, which is optimized in a manner similar to optimization of m and k for FloWPS. Namely, a maximum for AUC as a function of PC number is found and then used as the optimal number of PCs for an SVM-based prediction.

Thus, we compared the classifier qualities for three methods, namely classical SVM without data reduction, PCA-assisted SVM with pre-trained PC number, and FloWPS with the confidence parameter p = 0.90 (Table 3; note that both classical SVM and FloWPS calculations were done using gene expression features rather than PCs). The calculations were done for core marker gene datasets and cost/penalty SVM parameters C = 1 and 100. For linear kernel, several datasets had comparable AUC for simple PCA-assisted data reduction and for FloWPS (Table 3). However, for polynomial kernel FloWPS essentially outperformed the PCA-assisted data reduction, most likely due to bigger risk of overtraining for SVM with nonlinear kernels.

For the datasets investigated using the Affymetrix microarray hybridization platforms, gene expression data were taken from the series matrices deposited in the GEO public repository and then quantile-normalized (Bolstad et al., 2003) using the R package preprocessCore (Bolstad, 2018). All pediatric datasets taken from the TARGET database (Goldman et al., 2015) contained results of NGS mRNA profiling at Illumina HiSeq 2000 platforms; they were normalized using R package DESeq2 (Love et al., 2014).

All the SVM calculations with linear and polynomial kernels were performed using the Python package sklearn (Pedregosa et al., 2012) that employs the C++ library ‘libsvm’ (Chang and Lin, 2011). The penalty parameter C varied from 0.1 to 1000 for different calculations. Other SVM parameters had the default settings for the sklearn package.

AUC(m,k) topograms, like Figure 2B, were plotted using mathplotlib Python library (Hunter, 2007). Violin plots for FloWPS predictions (see Figure 3) for responders and non-responders were plotted using the ggplot2 R package (Wilkinson, 2011).