Systematic review of computational methods for drug combination prediction.

Synergistic effects between drugs are rare and highly context-dependent and patient-specific. Hence, there is a need to develop novel approaches to stratify patients for optimal therapy regimens, especially in the context of personalized design of combinatorial treatments. Computational methods enable systematic in-silico screening of combination effects, and can thereby prioritize most potent combinations for further testing, among the massive number of potential combinations. To help researchers to choose a prediction method that best fits for various real-world applications, we carried out a systematic literature review of 117 computational methods developed to date for drug combination prediction, and classified the methods in terms of their combination prediction tasks and input data requirements. Most current methods focus on prediction or classification of combination synergy, and only a few methods consider the efficacy and potential toxicity of the combinations, which are the key determinants of therapeutic success of drug treatments. Furthermore, there is a need to further develop methods that enable dose-specific predictions of combination effects across multiple doses, which is important for clinical translation of the predictions, as well as model-based identification of biomarkers predictive of heterogeneous drug combination responses. Even if most of the computational methods reviewed focus on anticancer applications, many of the modelling approaches are also applicable to antiviral and other diseases or indications.

Introduction

Drug combinations have become the standard of care for various complex diseases, including tuberculosis, malaria, HIV and other viral infections, as well as most of the advanced cancers[1]. In treatment of cancers and other multi-factorial complex diseases, multi-targeted treatments can offer therapeutic benefits both by enhancing treatment efficacy and by avoiding monotherapy resistance through inhibiting of multiple disease- or resistance driving signalling pathways or malignant cell subpopulations in heterogeneous patient and cell populations[2]. Especially in cancer treatment, therapy resistance is a critical clinical challenge for chemo-, targeted and immuno-therapies alike, and combination therapies are necessary to improve clinical benefits for most patients with refractory malignancies[3]. However, the current success rates of combinatorial trials remain limited because most clinically-used combinations have been identified based on empirical testing, and because of a lack of biomarkers to predict which patients would benefit from a specific combination regimen. Recent biomarker-based clinical trials, such as I-PREDICT, have demonstrated that targeting a larger faction of patient-specific molecular alterations, compared to limited gene panels used in monotherapy patient matching, may improve disease control rates without increasing adverse effects[4], [5]. However, how to preselect a panel of drugs for such combination trials remains an open question.

Since the number of possible drug combinations vastly exceeds that what could be tested clinically, identification of candidate drug combinations is often based on high-throughput screening (HTS) of the phenotypic effects of combinations in pre-clinical cell models (e.g. cell lines or patient-derived samples). The aim of HTS combinatorial discovery is to pre-select those drug combinations that provide higher-than-additive effect when compared to that using individual single-agents as monotherapies (so-called synergistic effect), yet at the same time would show minimal toxicity in healthy (non-malignant) cells. However, even with the automated HTS instruments, systematic screening of all possible drug combinations becomes impractical, both in terms of time and patient cells required, as the number of potential drug and dose combinations increases exponentially with the number of tested drug components and dose levels. Especially in anticancer applications, the inherent genetic and molecular heterogeneity of cancer cells makes drug combination synergy an exceptionally rare and highly context-dependent event[6], hence requiring the testing of the panels of combinations in various cellular contexts and genomic backgrounds to identify context-specific, selective combinations, rather than broadly active combinations that may lead to toxic effects.

Computational methods, such as those based on machine learning (ML) and other classes of artificial intelligence (AI), have shown their great potential for the prediction of most potent combinations to be prioritized in HTS for further preclinical or clinical development[7], [8], [9]. Even if such predictive methods can identify predictive markers and personalized combinations, their accuracy depends on the information available from the samples, drugs and targets, which poses further challenges for the practical use of these methods. Compared to the previous reviews, which have described the methodological details and data resources available for computational approaches to drug combination prediction[10], [11], [12], [13], [14], [15], this mini-review aims to provide a more practical guidance on the existing methods for various real-world applications in terms of the specific prediction tasks and input data required by the computational methods. For those readers interested in the algorithmic details of the methods, we refer to the previous methodological reviews[10], [12], [13]. In contrast, we will focus here on computational methods that provide also estimates of possible toxic effects of combinations, not only synergy predictions, and those methods that make predictions of multi-dose combination-response patterns, rather than only binary synergy classification, which are often required for clinical translation.

Methods

We carried out a systematic literature review of computational methods published so far for drug combination prediction (see Supplementary Methods for details of the PubMed keyword search). A total of 117 original research papers satisfied our selection criteria; we excluded review articles, studies that did not introduce a new computational method for combination prediction, and those articles not available as full text (see Supplementary Table 1 for detailed annotations of all the 138 papers identified). In the next step, we classified the computational methods based on the following four aspects related to practical application of the methods:

What input data are required by the method for the predictions (e.g., single-drug and/or drug combination responses, transcriptomics and other omics profiles, or structural information of the drug molecules; see Supplementary Table 2 for detailed categories).

What is the prediction task: (i) synergy classification, (ii) synergy score prediction (either regression or ranking of combinations in terms of synergy), or (iii) prediction of dose–response combination matrices (or tensors for higher-order combinations).

Whether the methods provide dose-specific information of drug combinations, or use as their input data either multi-dose single-drug or combination response data, which can also include missing responses (see Fig. 1), which are common in HTS data.

Fig. 1

Examples of experimental designs for pairwise multi-dose combination assays[16]. The first column on left and the bottom row indicate the monotherapy responses of the single-drugs d1 and d2 in the pairwise combination as a function of increasing dose (arrows). Grey cells in the dose-combination matrices and single-drug dose–response vectors indicate non-measured dose-responses. (a) Fully-measured design, (b) single fixed-dose design, (c) double fixed-dose design, (d) diagonal design, and (e) shifted design used in a recent drug combination screen[17].

Whether the methods make disease-selective predictions of combination effects (i.e., with limited toxic effects on healthy cells or individuals), or whether they make disease, sample or subpopulation-specific combination predictions that are tailored for a particular disease, patient, cell line or cell subpopulation (e.g. malignant cells only).

Results

Below, we highlight some of our key observations made based on the systematic literature review and analysis of the methods, especially related to the practical use of the combination prediction algorithms. We highlight select example methods that address important real-world application cases for which we expect further developments; see Supplementary Table 1 for more detailed annotations, including classification of the computational methods in terms of the above four practical aspects, along with a listing of their algorithmic classes, and whether there are open-access data, open-source code or online web-applications available.

DL-based methods are gaining popularity in the field of combination prediction

Over the last 15 years, a wide range of both ML and non-ML-based computational methods have been developed to identify potential drug combinations for further experimental testing (Fig. 2, Supplementary Table 1 for the method classes). Notably, methods that do not use training data to learn classification or prediction model (i.e., non-ML based approaches), have remained popular, due to limited availability of pharmacogenomic data for model training. The recent AstraZeneca-Sanger Drug Combination DREAM Challenge[8] benchmarked 160 prediction methods, both for synergy prediction and synergy classification, using blinded experimental data from 910 combinations across 85 molecularly characterized cancer cell lines. The top-performing methods, including the winning method that was based on Random Forest (RF), incorporated prior knowledge of drug-target interactions, and managed to predict the combination synergy with an accuracy matching that of biological replicates for greater than 60% of combinations. However, 20% of drug combinations remained poorly predicted by all the methods, regardless of the method category. This is consistent with other DREAM Challenges, where it has been observed that the specific prediction algorithm per se is not the critical component, rather how the algorithm is used in practice and based on what information often defines the best-performing computational methods for many prediction tasks in biomedical research[18], [19].

Fig. 2

Number of publications per year that introduce new computational methods for drug combination prediction. DL, deep learning; ML, machine learning. Note: the time axis is not linear; 2022 bar, as of mid-May 2022. The first two methods, introduced in years 2007 and 2009, were based on regression model with misclassification-penalized posterior (MiPP) and Medicinal Algorithmic Combinatorial Screen (MACS) with genetic algorithm to iteratively optimize potent combinations, respectively, and therefore classified as non-DL ML-based. The total number of publications is 117.

Like in many other fields of biomedicine, deep learning (DL) methods started to gain popularity in 2018 for drug combination predictions, and in 2020 these methods were already as common as the other computational prediction methods. The same trend seems to continue in 2022 in the proportion of methods classes (Fig. 2). However, comparative studies among various classes of drug combination prediction algorithms do not yet show consistent benefits from DL methods in terms of drug combination prediction accuracy[11], [13]. Especially in blinded test data, DL methods have shown limited accuracy for predicting the response of drugs that do not appear in the training dataset, indicating potential over-fitting to still limited amounts of the training data available[20]. In comparison to the non-DL and non-ML methods, the DL-based methods developed so far for drug combination prediction use more often molecular structures of drugs, as well as their physicochemical properties (Fig. 3, Supplementary Figure 1). One limitation of the complex DL algorithms is that their predictions are not always transparent, as they often rely on a large set of input data features, both structural and molecular, which may limit their practical implementation. Recently, more explainable DL methods have been developed for combination synergy prediction,[21] but there is still a need to make the learning algorithms and their outcomes more transparent before their widespread adoption[22].

Fig. 3

UpSet plot of input data used in the 22 DL methods for drug combination prediction. Left bars show the total number of methods that make use of an individual input data category. The filled dots and lines in the matrix indicate the common usage of specific input data combinations, and the top bars show the number of methods that use the particular input data intersection. Supplementary Fig. 1 shows the same analysis for the other method classes, and for all the 117 methods, and Supplementary Table 2 show the detailed annotations of the input data categories across the 117 prediction methods.

Methods that predict dose-dependent combination effects provide greater insights

Most of the current predictive modelling approaches treat combination prediction still as a binary classification or continuous regression problem, where the aim is to distinguish between synergistic combinations and those that show additive or antagonistic effects (Fig. 4, left). However, drug combination effects are not only context-dependent, but also highly dose-dependent[6], i.e., the same combination may show synergistic effects at one dose window and antagonistic effects at another dose window. Therefore, computational prediction methods should arguably provide dose-specific information of the combination responses, either the doses at which synergy exists, or ideally, predict combination effects across various dose combinations (see Fig. 4, right); such dose-specific combination predictions are also important for clinical applications, where lower dose combinations are often better tolerated by the patients. Furthermore, there exists a number of mathematical models for defining and quantifying synergy, each of which is formulated from different assumptions[1], [24], and therefore depending on the synergy model, one may end up making different interpretations of what is synergy, unless multi-dose responses are predicted. Finally, combinations can be highly effective in heterogeneous patient (or cell) populations also in the absence of any drug synergy or additivity[25], and such independent drug action of a combination can be sufficient to explain clinical benefit, especially in cases where the most effective single-agent varies across the patients[26], [27].

Fig. 4

Classification of the methods based on the combination prediction task. Drug combination synergy classification methods contain both binary and multi-class algorithms that aim to distinguish between synergistic and non-synergistic/antagonistic/additive combinations. Drug combination response prediction includes regression models for either synergy score or combination effect prediction, as well as ranking methods of combinations in terms of their combination effect. Dose-response combination effect prediction methods make predictions of multi-dose drug-response combination matrix (pairwise combinations) or tensor (higher-order combinations). Right panel: an example of pairwise combination dose–response matrix, where the highest synergistic effect is observed at lower doses of bryostatin 1 (dotted dose window), and the highest antagonistic effects at higher doses of nutlin-3; personalized and cancer-selective prediction for a leukaemia patient[23].

As an example of multi-dose combination response prediction methods (see Table 1), Zimmer et al. used a statistical regression model to predict cell line-specific effects of two or more antibiotics or anticancer drugs at multiple doses, based only on measurements of the drug pairs at a few doses, without requiring any mechanistic information or omics measurements[28]; they later extended their approach also to cases where measurements at only a single dose of drug pairs were used to predict high-order combination effects[29]. Similarly, Julkunen et al. used tensor learning to model the cell line-specific combination responses across various doses, through leveraging information from other combinations of similar drugs and using available omics data from the cancer cell lines, and achieved accurate predictions even in cases where no training data were available for the particular combination, as long as the single-drugs and cell lines had been tested in other combinations[30]. The same learning approach could also enable extensions to higher-order combinations, e.g., 3-drug cocktails, provided sufficient training data are available for the same drugs in similar cellular contexts (i.e., the dose–response tensors are densely populated). However, variability in combination screens, e.g., differences in the tested concentrations or frequent missing values, may pose challenges to learning approaches and lead to disagreements between synergy models[31]. As an alternative, prediction of dose-combination response surfaces, using either parametric or non-parametric models[16], [32], may lead to greater stability and insights into combination activity.

Table 1

Examples of computational methods for predicting dose–response matrices (pairwise combinations) or dose–response tensors (higher-order combinations).

Method input data	Prediction task	Specific predictions	Algorithm class	Disease application	Experiment validation	Ref
3-drug combination effects with specific dose–response matrix design in leukaemia cell line	Dose-response matrix prediction	Cell line-specific, cancer-specific	Artificial neural network with one hidden layer (non-DL ML-method)	T-lymphoblastic leukaemia cell line	Yes	[33]
Signalling pathways, NMR imaging structures, 5 active compounds, 14 target proteins	Dose-response matrix prediction	General combinations	Pathway network algorithm (non-ML method)	Inflammation	No	[34]
Single-drugs and pairwise combinations at a few doses: greater than10 dose-combinations for each drug pair	Higher-order dose–response matrix prediction	Cancer-specific, disease-specific, cell line-specific	Higher order regression (non-ML method)	Lung cancer cell line, antibacterial infection model	No	[29]
6 drugs in 2 cancer cell lines; single and combination drug responses at multiple doses	Higher-order dose–response matrix prediction	Cell line-specific, cancer-specific	Generalization of the Bliss regression models (non-DL ML method)	Cancer cell lines (2 cancer types)	No	[28]
A total of 23,595 pairwise drug combinations, dose–response matrices of various dose dimensions	Pairwise dose–response matrix prediction	Cell line-specific, cancer specific, disease-specific	Composite non-negative matrix factorization (non-DL ML method)	Cancer cell lines (4 cancer types), malaria and Ebola infection models	Yes	[35]
NCI-ALMANAC drug combination data (50 FDA-approved drugs; 60 cancer cell lines; 333,180 combination dose-responses); 'estate' molecular fingerprint; gene expression data	Pairwise dose- response matrix prediction	Cell line-specific, cancer-specific	Higher-order factorization machines (non-DL ML method)	Cancer cell lines (9 cancer types)	Yes	[30]

NMR, nuclear magnetic resonance; NCI, National Cancer Institute (U.S.); FDA, Food and Drug Administration (U.S.); ML, machine learning; DL, deep learning.

Disease-selective combination predictions are important for clinical translation

Observations in cancer cell lines have suggested that both combinatorial inhibition of two compensatory signalling pathways and targeting cancer survival pathway crosstalk may lead to combination synergies[6], [36]. However, the same co-inhibition may also affect healthy cells, rendering the combinations non-selective against malignant cells, and causing them to elicit toxic side effects also on healthy cells. However, most of the current methods either do not use any control data when making predictions, i.e., non-selective combination predictions, or they make non-specific predictions not tailored for a cell line or patient sample (Fig. 5). Such computational methods may lead to predicting broadly toxic combinations that unselectively kill various cell types, with potential side effects to non-malignant cell. The situation is similar in experimental HTS efforts, where the current pre-clinical selection of optimal combinations relies merely on the observed synergy between drugs[26], even if efficacy and potential toxic effects are the key determinants for the therapeutic success and tolerability of drug treatments in clinical practice. The development of disease-selective methods is not only relevant for anticancer applications, but they are also applicable to other complex diseases, which require capturing heterogeneity of the disease progression between cell populations at various stages of pathogenesis to identify both safe and effective treatments; e.g., finding combinations that synergistically inhibit virus replication, with minimal effects on non-infected host cells[37].

Fig. 5

Classification of the computational methods based on the selectivity and/or specificity of their combination predictions. Disease-selectivity means that the combination is predicted to co-inhibit mainly disease-related cells in patients, while having limited toxic effects on healthy cells or individuals (i.e., control data needed for toxicity evaluations). Disease, sample or subpopulation-specificity means that predictions are tailored either for a particular disease (e.g. AML), individual patient, cell line or cell subpopulations (e.g. blast cells in AML patients), instead of being general or non-specific combinations. Left bars show the number of methods in each category (color legend), and top bars their intersections.

As an example of cancer-selective drug combination prediction methods (see Table 2), He et al. identified patient-specific combination synergies in haematological cancer patients using genome-wide transcriptomic profiles, and estimated the toxic effects of drug combinations through differences in single-drug responses between cancer patients and healthy controls[23]; they also extended their approach to cell subpopulation-specific treatment predictions, using imaging-based drug testing assay to identify tumour cell-selective drug combinations for patients with solid tumours[40]. It was shown that simple RF algorithm provided relatively accurate predictions both for patient- and subpopulation-specific predictions. Single-cell profiling approaches, both for genotypic and molecular profiling and phenotypic drug responses, enable the identification of healthy cell populations and their co-inhibition effects directly from patient samples, respectively, and therefore avoid the need of having healthy control samples for toxicity estimation. To our knowledge, there is currently only one computational method that makes use of single-cell transcriptomic profiles for personalized combination predictions; Ianevski et al. used XGBoost algorithm to identify patient-specific and leukemic cell-selective drug combinations in primary samples from treatment-refractory leukaemia patients, each with different molecular backgrounds, by combining single-cell RNA-seq profiles with single-compound testing data through on– and off-targets of the compounds[41].

Table 2

Examples of computational methods for cancer-selective predictions of drug combination effects.

Method input data	Prediction task	Dose-specific	Subpopulation-specific	Control input data	Algorithm class	Cancer type	Experiment validation	Ref
Single-drug responses in 3 patients and 3 healthy donors for 218 drugs; bulk exome and RNA-sequencing	Drug combination synergy classification	No	No	Drug responses and transcriptomic profiles of healthy controls	Random forest (non-DL ML method)	T-cell prolymphocytic leukaemia (T-PLL)	Yes	[23]
114 single-drugs; 128 drug combinations at two doses; 155 combinations among nine drugs at three doses	Drug combination response prediction	Yes	No	Drug responses of control cell line	Quadratic phenotypic optimization (non-ML method)	Multiple myeloma	Yes	[38]
Single-drug and drug combination responses in pan-cancer cell lines	Drug combination response prediction	No	No	Average response over multiple cell lines	Multiobjective optimization (non-ML method)	NCI-ALMANAC cancer types	Yes	[36]
RNA-seq of prostate cancer patients and controls; gene expression responses; DTIs; PPIs, disease genes	Drug combination response prediction	No	No	Bulk transcriptomic profiles of healthy controls	Network integration and analysis (non-ML method)	Prostate cancer	Yes	[39]
Single-drug responses in 4 patient samples for 528 drugs; DTIs; Bulk gene expression and point mutations; scRNA-seq for one sample	Drug combination response prediction	No	Yes	Healthy subpopulation drug responses and transcriptomic profiles	XGBoost (non-DL ML method)	High-grade serous ovarian cancer	No	[40]
Single-drug responses in 4 patients for 456 drugs; DTIs, scRNA-seq for 4 samples	Drug combination response prediction	No	Yes	Healthy cell cell population transcriptomic profiles (scRNA-seq)	XGBoost (non-DL ML method)	Acute myeloid leukemia	Yes	[35]

DTI, drug-target interaction; PPI, protein–protein interaction; scRNA-seq, single-cell RNA-sequencing.

Conclusions and outlook

Given the increasing number of new computational methods introduced every year for drug combination prediction (Fig. 2), this is a timely review of their key features, before there appear too many articles for a systematic literature review. Rather than going into methods details, we focused here on their practical application: (i) identifying among the massive number of potential drug-dose combinations those that show maximal therapeutic potential and minimal toxic effects to be prioritized in the next phases of preclinical and clinical development; (ii) focusing on methods that enable learning from sparse and heterogeneous data sources (e.g., drug response data combined with genomic data) to provide patient-tailored solutions and response-predictive biomarkers for precision medicine; and (iii) highlighting the computational approaches that in the near future could bridge the gap to clinical translation and enable real-world applications, along with establishing their practical utility and impact in clinical decision making when optimizing combinations for patients.

Most of the current pre-clinical screening efforts emphasize merely the combination synergy as the key determinant of drug combination performance[26], while neglecting potential toxicity or selective efficacy (difference between efficacy and toxicity between malignant and non-malignant cells), even though these are critical factors for the clinical success. Notably, around 20% of drugs fail in the early development phase because of safety concerns (e.g., non-tolerated toxicity), whereas more than 50% fail due to lack of sufficient efficacy. Since there is a fundamental trade-off between clinical efficacy and tolerable toxicity, prioritization of drug combinations both in terms of their therapeutic and toxic effects during the pre-clinical investigation is critical to speed-up and de-risk the drug combination discovery before entering into lengthy and costly animal or clinical studies. Approaches that can guarantee maximal cancer-selectivity should thereby significantly accelerate the future design and testing of combination therapies, as well as increase the likelihood of their success in clinical studies.

Most of the methods reviewed here predict pairwise drug combination effects only, and often neglect dose-dependent responses, and many of the approaches cannot be easily extended to higher-order combinations of more than 2 drugs that are often required for the treatment of patients with e.g. advanced cancers. In general, combinatorial therapy is witnessing a paradigm shift from the traditional ‘two drugs in combination’ to the more complex ‘multi-drug cocktails’[42]. Therefore, there is a timely need for upgraded computational and experimental approaches that can effectively reduce the massive sampling space of higher-order drug combinations, and to identify not only synergistic but also safe multi-drug combinatorial therapies. Since systematic testing of hundreds of drug-dose combinations is impossible in scarce patient cells, novel computational methods are required that can make use of partial measurements of the combinatorial drug-dose spaces, along the lines initiated by Zimmer et al.[28], with the aim of identifying the most potent combinations across various dose windows.

Another methodological limitation is that the combinations are typically predicted at a single time point only (e.g. at diagnosis or relapse), and the methods do not take into account the dynamic process of treatment resistance and disease progression. This is because most methods are developed and tested in established cancer cell line models. However, translating the combination prediction results (or predictive biomarkers) from cell lines to individual cancer patients is not at all straightforward, whereas experimental testing of multiple combinations in patient-derived cells is often impossible in practice. Therefore, dynamic models directly applicable to patient-derived samples are highly needed, along the lines of a recent experimental-computational study that made use of high-throughput drug screening of cancer patient biopsies using a microfluidic assay, combined with logic-based modelling of signalling pathways to generate patient-specific dynamic models for predicting personalized combinatorial treatments with a limited number of cells from pancreatic cancer patients[43].

For practical use, the computational methods should be implemented as easy-to-use web-applications that also explain to the users how the combination predictions were made. Transparency in the drug combination response prediction methods involve many aspects of the method development, such as clear description of the predictions objectives (synergy, efficacy and/or toxicity), as well as quantitative performance and confidence evaluation (cross-validation, wet-lab validations, and based on clinical data), which will help experimental and translational professionals to decide when and how to use the methods to obtain valid results and to improve either combinatorial designs or clinical decision making. For routine clinical implementation, there is also a need to develop computational methods that ensure cost-efficiency, explainability and interpretability of drug combination predictions, through scoring feature importance and sparse modelling approaches. For instance, feature selection procedures can facilitate pinpointing the molecular markers most predictive of combination responses, both for wet-lab validations and economical clinical implementation.

In conclusion, although a number of computational methods have been developed to address many important experimental and translational challenges, there is still a need to implement novel computational solutions and to demonstrate their feasibility and benefits in translational applications, e.g., in cancer and antiviral applications, where there is an urgent need to identify combinatorial therapies for each patient individually based on patient-specific biomarkers. This is expected to lead to novel combinatorial design methodology with decision support tools to optimally exploit patient-specific therapeutic vulnerabilities to identify drug combinations that selectively co-inhibit disease cells but avoid severe targeting of healthy cells.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.