SYNERGxDB: an integrative pharmacogenomic portal to identify synergistic drug combinations for precision oncology.

Drug-combination data portals have recently been introduced to mine huge amounts of pharmacological data with the aim of improving current chemotherapy strategies. However, these portals have only been investigated for isolated datasets, and molecular profiles of cancer cell lines are lacking. Here we developed a cloud-based pharmacogenomics portal called SYNERGxDB (http://SYNERGxDB.ca/) that integrates multiple high-throughput drug-combination studies with molecular and pharmacological profiles of a large panel of cancer cell lines. This portal enables the identification of synergistic drug combinations through harmonization and unified computational analysis. We integrated nine of the largest drug combination datasets from both academic groups and pharmaceutical companies, resulting in 22 507 unique drug combinations (1977 unique compounds) screened against 151 cancer cell lines. This data compendium includes metabolomics, gene expression, copy number and mutation profiles of the cancer cell lines. In addition, SYNERGxDB provides analytical tools to discover effective therapeutic combinations and predictive biomarkers across cancer, including specific types. Combining molecular and pharmacological profiles, we systematically explored the large space of univariate predictors of drug synergism. SYNERGxDB constitutes a comprehensive resource that opens new avenues of research for exploring the mechanism of action for drug synergy with the potential of identifying new treatment strategies for cancer patients.

INTRODUCTION

Despite tremendous investments by pharmaceutical companies into anticancer drugs, patients often still either fail to respond to first-line chemotherapy or acquire resistance after initially responding to monotherapy in clinical settings (1,2). Treatment options are chosen on the basis of cancer type, TNM staging, and/or the physician's experience regardless of the patient's characteristics. This one-size-fits-all approach of chemotherapy is challenging because it assumes all drugs elicit the same response, regardless of patient characteristics. Recommendations, including the first-line treatments, generally lack specificity and may result in adverse drug reactions for a given patient (3). On the other hand, molecularly targeted agents that exert more-specific and less-toxic effects toward patients often elicit promising initial pathological responses (4). Nevertheless, targeted agents are only beneficial in the subset of patients who possess the targetable mutations (5). These problems have led to the use of combinations of approved drugs and/or investigational compounds for the rapid development of new therapeutics for cancers for which effective chemotherapy interventions are not available (6).

Many studies have demonstrated that patients who undergo combination therapy show favorable survival outcomes compared with monotherapy in the treatment of tumors, which gives this approach huge potential for overcoming cancer treatment failures (7–9). Although combination therapies represent promising treatment options in clinical settings, experimentally testing all possible combinations is not practical because of the large space of possible drug combinations and the high cost and resources required for testing these therapeutic strategies in clinical trials. Thus, there is an urgent need for an efficient approach to prioritize combination therapies with therapeutic potential.

Recent advances in high-throughput drug screening have provided an unprecedented opportunity to mine huge numbers of pharmacological profiles for predicting the synergistic effects in combining approved drugs and investigational chemical compounds in preclinical model systems (10). However, simultaneously analyzing pharmacogenomics datasets obtained in independent studies is made difficult by the lack of standardization of drug and cell line identifiers, and the diversity of computational methods used to quantify drug synergistic effects. To alleviate these limitations, the DrugComb (11) and DrugCombDB (12) databases have recently been released, which include pharmacological data for 437 932 and 448 555 drug combinations tested on 93 and 124 cell lines, respectively. However, these databases suffer from multiple limitations. First, neither the DrugComb nor DrugCombDB database includes molecular profiles of the preclinical model systems screened against the reported drug combinations. Second, the databases mainly consider the prediction of synergism/antagonism for a given combination, and they do not provide analysis modules to merge or compare multiple datasets for identifying robust biomarkers or drug combinations. Third, they do not allow users to directly compare the effect of drug combinations due to the lack of standardization of drug identifiers, which comprises a mixture of compound names (e.g. gemcitabine), trade names (e.g. zolinza [vorinostat]), investigational compound names (e.g. AZD1775 [adavosertib]), and abbreviations of the drug names (e.g. 5-FU [fluorouracil]). A comprehensive computational platform therefore needs to be developed for the systematic investigation of the impact of genetics on variability in combination responses—this represents a key challenge to achieving precision oncology.

Here we developed SYNERGxDB (http://SYNERGxDB.ca/), which is a Web application that includes the largest database of nine collections of pharmacological (477 839 drug combinations tested on 151 cell lines) and molecular profiles of preclinical model systems. This application enables clinicians and researchers to explore and predict the synergy of drug combinations in preclinical models (Figure 1). To overcome the lack of standardization across datasets, we leveraged Cellosaurus (13), which is the most-comprehensive catalog of cell lines, and PubChem and DrugBank (Supplementary Figure S1). (14,15) to uniquely identify and comprehensively annotate the names of cell lines and drugs, respectively. The SYNERGxDB web server allows users to (i) query drug combinations in multiple datasets, (ii) visualize sensitivity data of the drug combinations across datasets, (iii) identify drug combinations with therapeutic potential, (iv) discover candidate predictive biomarkers for a given drug combination and (v) mine existing data to optimize and design future drug screening studies. Moreover, SYNERGxDB provides analytic tools to leverage molecular and pharmacological profiles of a large panel of cancer cell lines on top of a transparent architecture with an optimal framework for predicting the synergistic combinations through repurposing and re-evaluating existing drug combinations.

Figure 1.

Schematic overview of SYNERGxDB. Pharmacological and molecular profiles were integrated into the database using standardized cell-line and compound names. SYNERGxDB provides a query interface where users can select cell line(s), compound(s), or dataset(s), along with analysis modules including for biomarker discovery, and effective visualization functions so as to produce a prioritized list of drug combinations.

MATERIALS AND IMPLEMENTATION

Compendium of drug combination pharmacogenomic datasets

Both academic groups and pharmaceutical companies who generated pharmacological and molecular profiles of large panels of immortalized cancer cell lines have reported the synergism of cancer cell lines to drug combinations. Despite their considerable therapeutic potential, the lack of standardized cancer cell line and drug/compound annotations and quantification methods hinders the systematic calculation and prediction of drug synergism or antagonism (hereafter referred to as synergy scores). The ability to simultaneously analyze multiple studies can increase the statistical power and improve the robustness of predictors of drug combinations, which would constitute a major step toward designing new therapeutic strategies in cancer.

To address these issues, we integrated nine of the largest drug combination datasets and then applied a semiautomated curation process (Table 1) that maximized the overlap among datasets. In addition to the standardization, unified methods to quantify synergy scores are required to compare and combine across datasets obtained in studies that have applied different protocols for high-throughput drug combination screening. First, the number of doses in the combinations differs across datasets regardless of the drug dose used in each experiment. Lanevski et al. (DECREASE dataset) (16) and Langdon et al. (YALE-PDAC dataset) (17) performed 8-by-8 and 3-by-3 combinations of the drug concentrations, respectively, while Forcina et al. (STANFORD dataset) (18) and Licciardello et al. (CLOUD dataset) (19) performed 1-by-1 combinations for synthetic lethality screening, which resulted in the inability to apply ZIP (zero interaction potency) and Bliss methods since they require a fitted curve (hill slope) and IC50, respectively. Friedman et al. (MIT-MELANOMA dataset) (20) introduced a unique design that involved performing two pairs of experiments with drugs in low and high concentrations for a given combination. Wali et al. (YALE-TNBC dataset) (21) used a single fixed dose of FDA-approved drugs in combination with the experimental drugs at five doses. Second, two datasets included replicated experiments (MERCK and VISAGE datasets) (22,23), while the other datasets did not include replicates. For MERCK and VISAGE datasets, the median viability from the experiments was taken as the observed value, resulting in lower numbers of experiments and measurements. The median viability from the experiments was taken as the observed value and average and standard deviation of the synergy scores in 34 preclinical model systems are provided in the Supplementary Figure S2 to assess its quality. Third, the drug doses do not match between mono- and dual-therapy screening in the MERCK dataset (22), requiring the imputation of the viability for the corresponding doses in dual therapy from the drug dose–response curve of the monotherapy using logLogisticRegression and getHill functions implemented in the R PharmacoGx package (version 1.6.1) (25).

Table 1.

Datasets and statistics for screening drug combinations

Dataset	No. of cell lines	No. of compounds	No. of combinations	No. of experiments	No. of measurements	Experimental design	Reference
NCI-ALMANAC	60	101*	5354	311 407	4 567 145	3-by-3 or 5-by-3	(22)
MERCK	39	38	583	22 737	570 645	4-by-4	(23)
MIT-MELANOMA	36	8	5778	201 254	1 407 123	2-by-2	(20)
VISAGE	34	2	1	34	2040	10-by-6	(24)
DECREASE	13	33	36	210	13 440	8-by-8	(16)
YALE-TNBC	6	130	768	4576	54 912	1-by-5	(21)
YALE-PDAC	4	41	861	3326	50 707	3-by-3	(17)
STANFORD	1	1818	1818	1818	7272	1-by-1	(18)
CLOUD	1	55	1327	1327	5308	1-by-1	(19)
TOTAL (unique)	151	1977	22 507	536 596	6 678 592

*Screening data not available for three compounds (4′-Epiadriamycin, Eribulin mesylate, and Idarubicin hydrochloride).

To handle the variability arising from experimental design issues, we constructed a unified database of drug screening in combination with the curated datasets, resulting in the inclusion of a considerable number of experiments that were common among the studies as well as a substantial number of preclinical models. All of the datasets that we have integrated into the SYNERGxDB are accessible on the Datasets page of the application, or as individual articles in the Supplementary Data.

Integration of the pharmacological and molecular profiles and data statistics

The cancer cell lines used in the major drug-combination studies included in SYNERGxDB have been extensively profiled for their molecular features at the DNA, RNA, methylation, metabolomic, and proteomic level (26). We therefore integrated the pharmacological profiles of the cell lines with 136 metabolomic activities, 19 068 protein-coding gene expressions, copy number in 18 281 genes, and 114 568 variants harbored in 17 894 genes for the corresponding cancer cell lines from Cell Model Passport (27) and Cancer Dependency Map (28). In order to standardize the cell line and compound names, we developed a semiautomated system in which they were uniquely mapped to Cellosaurus (13), while drug names were mapped to PubChem (14) and DrugBank (15).

Our curation process yielded 151 cell lines, 15 tissue types, and 1977 unique drugs/compounds. From the nine experimental datasets integrated into SYNERGxDB, 22 507 unique drug combinations were screened (536 596 experiments involving 6 678 592 measurements), and synergy scores were calculated to evaluate combination synergy using the following four statistical models: ZIP (29), Bliss independence (30), Loewe additivity (31), and HSA (highest single agent) (32). Synergy scores were calculated using the synergyfinder package (version 2.0.12) of R software (version 3.6.2) (33).

Web implementation

The web interface was implemented using Node (version 10.16.0) and Express (version 4.17.1) on the back end, with ReactJS (version 16.12.0) used on the front end to ensure rapid rendering and high performance. Visualizations were plotted using the JavaScript libraries d3.js (version 5) and PlotlyJS, which are both interactive and allow dynamic visualizations; they were constructed in HTML, CSS, and SVG. The database is hosted on a MySQL server to support large batch queries against multiple tables using the InnoDB storage engine to ensure ACID compliance and transactional support. In addition, we built a responsive website that will display components of the portal with ease and provide a better visual and interactive experience for the user.

The MySQL dump (version 0.2.1) is available on Zenodo at http://doi.org/10.5281/zenodo.3780920, and the relational schema is provided in Supplementary Figure S3. All data in the tables in SYNERGxDB are downloadable in CSV format. The Web application is hosted on Azure cloud services, and it leverages all of the benefits provided by PaaS (platform as a service) solutions specifically designed to guarantee high levels of security, performance, and flexibility for web resources. The back-end server runs under Azure Web App Service, while the database layer of the application is hosted using Azure Database for MySQL Server Service. In addition, some of the real-time computation for the application is performed by the OpenCPU server that is available on a separate Azure Virtual Machine.

Code and documentation

The SYNERGxDB code is open source and publicly available on GitHub (https://GitHub.com/bhklab/SYNERGxDB). Detailed documentation is available in the SYNERGxDB Web application (http://SYNERGxDB.ca/documentation/), which includes examples of use cases and URLs for searching, compounds, samples, datasets and synergy score summary pages. The documentation also details the analysis pages and visualizations, such as biomarker discovery, cell line sensitivity analysis, tissue-specific enrichment analysis and further synergy score analysis.

WEB-INTERFACE AND ANALYSIS FUNCTIONS

Search synergy scores across datasets and retrieve detailed information

The SYNERGxDB search engine allows users to explore potential biomarkers and drug combinations according to synergy scores, by querying cell line(s), two compounds in a combination, and/or the dataset of choice. Compounds can be identified by name, ATC classification, PubChem CID and DrugBank ID, while cell lines can be identified by name and Cellosaurus ID. For example, users can retrieve the synergy scores for the drug combination of bortezomib and topotecan (BT-combination) from the NCI-ALMANAC and MERCK datasets. A query from the front page would be displayed by a query panel at the top of the page and the list of respective combination drugs with synergy scores at the bottom of the page.

Users can further sort the table according to the synergy scores, cell line(s), compounds, or tissues, to produce a page that displays detailed metadata of each experiment, which can be downloaded as a table in CSV format. Users can retrieve detailed information by clicking each row in the table and the names of the cell lines, compound names, SMILES, InChl Key, and synergy scores, along with the source of the datasets, which is displayed at the top of the page. A horizontal bar plot in the middle of the page displays the cell-line rank based on synergy scores. Synergy matrices are provided for each score that is calculated for an experiment, where concentrations of the two analyzed compounds and corresponding inhibition percentage values are presented.

Users are able to further investigate the relationship between drug combinations and synergy scores through 3D surface plots, on which the degree of synergism is depicted using color intensity. The bottom of the page displays the dose–response curve for the two given drugs, along with the inhibition percentage values in a heat map. Additionally, SYNERGxDB supports the RESTful (Representational State Transfer) API (Application Programming Interface), which allows users to directly query the database without having to use a Web application interface (Supplementary Data).

Analysis modules to evaluate combination drug screening

SYNERGxDB provides the following four analysis modules when users select more than ten cell lines in a query: (i) biomarker discovery, (ii) cell-line sensitivity analysis, (iii) tissue-specific enrichment analysis and (iv) consistency in synergy scores.

Biomarker discovery

This module displays the associations between gene expression (FPKM [fragments per kilobase of exon model per million reads mapped]) and synergy score in each dataset and provides correlation metrics, including the concordance index (C-index) (34), in order to identify potential predictive biomarkers in combination therapies (Figure 2). Association testing between gene expression and each synergy score metric in each dataset was carried out, where users can retrieve the strength and significance of a correlation in each tab for four different scoring methods (ZIP, Bliss, Loewe and HSA) to facilitate the identification of potential predictive biomarkers. Genes that are significant in multiple datasets are displayed on the top of the table regardless of their P-values. The association between gene expression and synergy scores is displayed on a scatter plot, where users can define two groups of cell lines, such as those with high and low synergy scores. Differential gene expression between the two groups can then be analyzed, with the expression distributions being displayed interactively.

Figure 2.

BT-combination associations for biomarker discovery. The association between gene expression and Bliss synergy score was analyzed for each gene in each dataset, with the resulting table sorted by P values. Users can select a gene symbol to generate a scatter plot between the FPKM and Bliss synergy score. Distributions of the FPKM in groups with high (>0.13) and low (≤0.13) Bliss synergy scores are displayed in a corresponding box plot on the right.

Cell-line sensitivity analysis

This module provides the summarized synergistic/antagonistic patterns of drug combinations in a single heat map across multiple datasets within a given tissue type or a set of cell lines selected by users (Figure 3). The heat map displays cell lines from multiple datasets as columns and drug combinations as rows. The queried combination is displayed in red font (e.g. BT-combination), while compound(s) from the query are displayed in yellow or blue font. This allows users to identify which drug in a given database shows the strongest synergistic effects when one drug is administered in combination with others, or discover ones that provide greater synergism than the given combination. Finally, the distribution of synergy scores in each combination will be displayed in box plots that are aligned with each drug combination in the right panel.

Figure 3.

Comparison of synergy scores across datasets. A heat map displaying ZIP synergy scores in each cell line (X-axis) for every drug combination (Y-axis). The queried combination is displayed in red font (e.g. BT-combination), while compound(s) from the query are displayed in yellow or blue font. The box plot on the right side displays the distributions of the ZIP synergy scores sorted by the median values. Different datasets are shown using different color keys in boxes: orange and purple indicate the NCI-ALMANAC and MERCK datasets.

Tissue-specific enrichment analysis

This module displays the synergistic effect of the drug combination on specific tissue types. In order to test whether higher synergy scores were associated with a specific tissue type, we ranked the cell lines by their synergy scores and then calculated the area under the receiver operating characteristics curve (AUC) for each tissue type by comparing the tissue with the others in the combinations. AUC values were calculated using the R pROC package (35), and tissue-specific AUCs are displayed in a single table sorted by AUCs. We also used four synergy scoring methods to compare AUCs across multiple tissue types and datasets.

Consistency among datasets

This module offers users the ability to compare two synergy scoring methods for determining synergy scores on a scatter plot, with concordance quantified using the C-index, Spearman rho, and Pearson r (Figure 4). Users can select a synergy scoring metric on each axis to check whether they are consistent across methods. If users select to display the distributions in multiple datasets, a scatter plot at the top will display synergy scores in an integrated dataset, while the others display in each dataset along with statistical testing. Cell lines are indicated as dots that are color coded depending on the source of the datasets, allowing users to see trends in synergy scores across datasets. The names of the cell lines are displayed when the mouse is moved over the dots (i.e. cell lines) on the scatter plot.

Figure 4.

Comparison of the synergy scoring method for the BT-combination. Scatter plot showing the results of testing the BT-combination against 97 cell lines from NCI-ALMANAC (orange) and MERCK (purple) datasets. Three coefficients (C-index, Spearman rho, and Pearson r) were calculated for the correlation between the Bliss synergy score (X-axis) and ZIP synergy score (Y-axis). The cell line name (i.e. K-562) is displayed as a user moves the cursor over a data point on the plot.

Pharmacogenomics analysis in drug combinations

The integration of both molecular and pharmacological profiles of the preclinical model systems allows users to identify potential biomarkers in combination therapies. In the Pharmacogenomics tab, users can select one matrix from either molecular or metabolomic profiles and one synergy scoring method to perform association tests in multiple datasets, while the biomarker discovery module performs association tests in each dataset. To query the database, users need to select the following options: data type (e.g. molecular profile), feature (e.g. gene), sample(s), a pair of drugs, and one metric of synergy scores (e.g. ZIP). For example, we focused on expression-based biomarkers in the BT-combination, which dentified ATP6V1B2 as one of the potential biomarkers in the merged datasets (i.e., NCI-ALMANAC and MERCK) (N = 82, C-index = 0.634, P = 2.6E–05) (Figure 5). Because ATP6V1B2 was reported to be one of the essential basal-specific genes in breast cancer (36), this biomarker could be used to identify samples that would benefit from this drug combination in clinical settings.

Figure 5.

Pharmacogenomics analysis of the BT-combination. Association between ATP6V1B2 expression (X-axis) and ZIP synergy score (Y-axis) in 82 cell lines for which gene expression profiles were available. Selected parameters and correlation coefficients are displayed on the right side.

SUMMARY AND FUTURE DIRECTIONS

SYNERGxDB provides a unified framework to prioritize synergistic combinations by repurposing and re-evaluating existing drug combinations, allowing users to access and analyze molecular profiles to identify potential biomarkers. To our knowledge, SYNERGxDB presents the largest integrated database of both molecular and pharmacological profiles for drug combinations tested in in vitro cancer models, with rigorous curation having been applied to identifiers across datasets. Lastly, the database allows users to identify potential biomarkers, which will need further preclinical and clinical validation.

The integration of multiple profiles from different domains can provide a more-powerful and systematic approach to answering complex pharmacological and biological questions. There is a need for a computational framework that enables the integrated analysis of these complementary profiles to improve predictions of drug synergism in preclinical settings. Although machine-learning and multivariate statistical approaches are being increasingly used for predicting drug synergism, many prediction methods have used only pharmacological profiles in isolated datasets. SYNERGxDB can be further used for machine-learning analysis to identify and test predictors of drug synergies using both pharmacological and molecular profiles. We plan to integrate DrugComboRanker (37) and RACS (Ranking system of Anti-Cancer Synergy) (38) in our Web application along with our recent DNF (Drug Network Fusion) algorithm (39) to prioritize the list of drugs in combinations based on their chemical structures, transcriptomic perturbations, and the single-drug sensitivities. SYNERGxDB will ultimately provide a unique resource for investigating the mechanisms underlying variant-dependent responses to combination therapies using molecular profiles of cancer cell lines, which is a key step toward achieving precision oncology.

Click here for additional data file.