Modeling Cellular Response in Large-Scale Radiogenomic Databases to Advance Precision Radiotherapy

Radiotherapy is routinely used as curative therapy for patients with cancer. Recent technologic advances have considerably augmented the physical precision of radiotherapy, resulting in improved cure rates and less toxicity (1–3). Biologically motivated improvements (such as the addition of radiosensitizing drugs) to radiotherapy delivery have not seen such dramatic improvements despite the known differences in radiation efficacy that exist among patients with a particular tumor type. This is due in part to a lack of predictive biomarkers on which to stratify patients. Instead, the stratification of patients to different radiotherapy-containing regimens continues to be based primarily on clinical variables such as tumor stage.

The biological determinants of cellular response to radiation, referred to as radioresponse, are complex and include both genomically-based cell-intrinsic and external microenvironmental factors (4, 5). Intrinsic radiosensitivity varies among individual tumors of the same type with implications for optimal radiotherapy dosing and curability. Measurement of intrinsic radiosensitivity in molecularly-characterized cancer cell lines could provide the radiogenomic data necessary to develop radioresponse predictors. However, despite decades of research there remains no clinically utilized radiosensitivity biomarker discovered from cell culture radiogenomic studies. Reasons for this include the need for clonogenic assays when measuring intrinsic radiosensitivity in vitro, which are cumbersome and not amenable to large screens (6, 7). Furthermore, radiosensitivity varies with dose in a complex and tumor-specific manner, rendering measurements at multiple dose levels a necessity.

Most short-term cytotoxicity assays amenable for high-throughput analysis of drug response have endpoints at 72 hours. These assays are inappropriate for measuring radiosensitivity because of the delayed cellular death by mitotic catastrophe that often occurs following radiation exposure (8). To address this limitation, an extended (9-day) viability assay was developed as a surrogate for clonogenic survival that is amenable to high-throughput processing (9). This assay was recently applied to 533 cancer cell lines with multiple radiation dose levels (10), becoming the largest radioresponse dataset published by a significant margin. This increase in scale of radioresponse data could lead to robust predictive biomarkers. However, full utilization by the research community requires sophisticated analysis tools that can appropriately model cellular response to radiation and seamlessly integrate associated molecular and pharmacogenomic profiles of cell lines.

In this study, we performed a preclinical assessment of intrinsic radiosensitivity using large-scale radiogenomic datasets. We sought to (i) model dose–response data using the linear-quadratic (LQ) model (11); (ii) integrate the modeled radioresponse profiles with transcriptomic data to determine pathway- and tissue-specific determinants of radioresponse; (iii) identify mutations associated with radioresponse; (iv) estimate radioresponse under hypoxic conditions; and (v) identify classes of drugs with cytotoxic effects that correlate with radioresponse. To facilitate these and other future analyses, we developed RadioGx, a new computational toolbox enabling comparative and integrative analysis of radiogenomic datasets. Our work provides a framework for future hypothesis generation and preclinical assessments of radioresponse using appropriate biological assays and indicators.

Supplementary Table S1 presents the sensitivity and transcriptomic datasets that are used in this study. Supplementary Table S2 presents the functionality of the RadioGx package. Cell line genomic studies often lack standardized identifiers. To overcome this, we assigned a unique identifier to each cell line and radiation treatment. Within RadioGx we implemented a RadioSet (RSet), a data container storing radiation dose-response and molecular data along with experimental metadata (detailed structure provided in Supplementary Fig. S1). In addition, the RSet also enables efficient implementation of curated annotations and molecular features for cell lines, which facilitates comparisons between different datasets. We have implemented a unique set of functions that enables users to analyze radiogenomic datasets. One of the primary functions is the downloadRSet that allows users to download the RadiationSet (RSet) object. We have also incorporated a function, linearQuadraticModel, which fits the radiation cell survival data using the standard radio-biological formalism, the LQ model. This function uses a normal error distribution by default, but users also have an option to choose Cauchy distribution. For a given dataset by the end user, this function fits the dataset with the LQ model, and returns radiobiological parameters alpha and beta along with the goodness-of-fit. To extract several features from this curve, we have implemented the functions computeAUC, which enables users to compute area under the survival curve (AUC), computeSF2 function, which returns the fraction of cells that survive a radiation dose of 2 Gy, and computeD10 function, which returns the radiation dose at which only 10% of cells survive.

Radiobiologic modeling is used to allow comparisons of various clinically relevant radiotherapy treatment regimens. The most common formulation in current clinical practice is the LQ model (11), which assumes that there are 2 components to cell killing induced by radiation: one that is proportional to dose (linear, |\alpha $|⁠) and another that is proportional to the square of the dose (quadratic, |\beta $|⁠). The LQ model describes the fraction of cells that survived (S) a uniform dose D (Gy); the survival fraction of cells after irradiating with an acute dose D is given by:

The ratio \displaystyle |\frac{\alpha }{\displaystyle \beta }$| varies by the cell population or tissue that is being irradiated, and reflects the response to different fractionation schemes. Cell populations or tissues with a high value are less sensitive to the effects of fractionation than those with a low value.

Mutation data were obtained through the Cancer Cell Line Encyclopedia (CCLE) data portal (DepMap 18Q3 data release). The mutation data used in the analysis are incorporated into the RadioSet as one of the molecular profiles. Ensembl Variant Effect Predictor (VEP; Ensembl release 96 – April 2019) was used to annotate predicted functional impact of mutations. Impact ratings were categorized into 'High' and 'Not High', with the latter category encompassing “Modifier,” “Moderate,” and “Low” impact mutations (12). For genes with multiple mutations, the highest impact mutation was considered. Radiation response AUC values were compared between the 3 groups: (i) cell lines with VEP-High mutations; (ii) cell lines with mutations that were not VEP-High; and (iii) cell lines with no mutations (wildtype, WT). Wilcoxon U test compared VEP-High and WT groups. Cohen d values to assess effect size were calculated for each individual gene using the effsize package (v0.7.4). For thick forest plot creation, we used the metaviz package (v0.3.0) with default parameters.

The LQ model can also be used to model the effect of hypoxia. Hypoxia is a hallmark of many solid malignant tumors and influences tumor progression, therapy resistance, development of metastases, clinical behavior, and response to conventional treatments like radiotherapy. The survival fraction of cells due to a given radiotherapy dose is given by Eq. A under well-oxygenated, or normoxic conditions. However, the surviving fraction of cells may vary depending on the amount of oxygen concentration in the tumor, as cells in the hypoxic region are considered to be more resistant to radiation therapy. This hypoxic effect can be incorporated into the LQ model using the, “oxygen enhancement ratio (OER),” which can be normalized to yield the “oxygen modification factor (OMF)” (13). OMF is defined as follows:

where |{O_2}$| is the oxygen concentration in the system in mm Hg, |{K_m}\ = \ 3$| mm Hg, defined as the oxygen at which half of the ratio is achieved, and |{\rm{OE}}{{\rm{R}}_m}\ = \ 3$| is the maximum value at well-oxygenated condition. Therefore, the LQ model given in Eq. A can be modified to include oxygen concentration as follows:

In general, the OER can be a function of radiation dose. Some studies have suggested that the maximal oxygen enhancement varies in the range of 2.5 to 3 with differences in radiation dosage (14–16). This can be simply included into the revised LQ model by considering different OERs for the parameters α and β, that is, |{\rm{OE}}{{\rm{R}}_\alpha }$| and |{\rm{OE}}{{\rm{R}}_\beta }$|⁠. However, because we consider the normalized OER (or, OMF), the introduction of these separate terms will not produce a significant difference in the final survival fraction. Thus, we assume |{\rm{OE}}{{\rm{R}}_\alpha }= {\rm{OE}}{{\rm{R}}_\beta }$| in our mathematical framework. We assume that the system is moderately hypoxic, that is, approximately 5 mm HG for this study.

We used CTRPv2 dataset in PharmacoGx package (v1.10.3; ref. 17) that has 545 drugs to compute the association between radioresponse and drug response (defined by the AUC of the Hill function). We also performed pharmacologic enrichment analysis, an adaptation of the gene set enrichment analysis (GSEA) methodology. For this, we computed the correlation of radioresponse with each drug response, and a pharmacologic set represents a gene set. Similar to the GSEA method, a running sum is calculated, starting with the first compound-level statistic to the last. The sum is increased if a compound-level statistic belongs to the pharmacologic class of interest, otherwise, the sum is decreased. The enrichment score of the pharmacologic class of interest is defined as the maximum deviation from zero of the running sum (Supplementary Fig. S2; ref. 18).

Pathway enrichment analysis on the gene expression data was carried out using GSEA (19) with pathways defined by QIAGEN's Ingenuity Pathway Analysis (IPA; QIAGEN Redwood City, www.qiagen.com/ingenuity). Genes were ranked based on their coefficient of correlation between the gene expressions and the radioresponse of interest (AUC or SF2). GSEA was then used to compute the enrichment score for each pathway with statistical significance calculated using a permutation test (10,000 permutations) as implemented in the piano package (20). Nominal P values obtained for each pathway are corrected for multiple testing using the false discovery approach (FDR; ref. 21).

RadioGx is implemented in R and is freely available from the Comprehensive R Archive Network (CRAN) from cran.r-project.org/web/packages/RadioGx/. The code, documentation, and detailed tutorial describing how to run our pipeline and reproduce our analysis results are open-source and publicly available through the RadioGx GitHub repository (https://github.com/bhklab/RadioGx-analysis). A virtual machine reproducing the full software environment is available on Code Ocean (DOI: 10.1101/449793). Our study complies with the guidelines outlined in refs. 22–24. All the data are available in the form of RSet objects with associated digital object identifiers (DOI).

To realize the full potential of large-scale radiogenomics datasets for robust biomarker discovery, we developed the RadioGx software package (Supplementary Fig. S1). RadioGx represents the first computational toolbox that integrates radioresponse data with radiobiologic modeling and molecular data from hundreds of cancer cell lines. Within RadioGx, datasets are standardized with comprehensive cell line annotations including the type of radioresponse assay (i.e., clonogenic assay and 9-day viability assay) and indicators used to generate dose-response data (i.e., SF2 and AUC). RadioGx enables fitting of dose–response data using established radiobiologic models, quality control in order to investigate the consistency and biological plausibility of radioresponse assays and indicators, and integration of these data with other data types and radioresponse models.

Multiple dose–response measurements from the same cell line can be incorporated into established radiobiological models to predict the effect of specific perturbations (e.g., radiotherapy fraction size or hypoxia) on radioresponse. Within RadioGx, we applied the commonly used LQ model to fit 9-day viability assay data for 533 cancer cell lines (Fig. 1A; ref. 10). The LQ model goodness-of-fit was high for the majority of cell lines (median R² = 0.958; Supplementary Fig. S3A and S3B). For 498/533 (93%) of cell lines, the model fit the experimental data reasonably well (R² ≥ 0.6); these cell lines were retained for subsequent analyses.

Using the LQ model for each cell line, we calculated AUC as a summary radioresponse indicator that is independent of a specific dose level. As expected, a range of radioresponse profiles were seen (Fig. 1B). We next compared AUC and dose-specific survival data (SF2, SF4, SF6, and SF8) from the 9-day viability assay with clonogenic survival data generated by Yard and colleagues for a subset of cell lines (Fig. 1C; Supplementary Fig. S4A–S4E; ref. 10). We observed high Pearson correlation (R ≥ 0.8) for AUC (n = 15), SF2 (n = 12), SF4 (n = 15), and SF6 (n = 15), but SF8 showed only moderate correlation (R = 0.64; n = 11), consistent with prior observations suggesting poor reproducibility of survival assays following high radiation doses (25). Taken together, the 9-day viability assay provides a robust surrogate for clonogenic survival at a range of radiation doses. Moreover, the LQ model within RadioGx allows for characterization of radioresponse and derivation of radioresponse indicators for the vast majority of cancer cell lines.

Summary indicators of radioresponse are useful for preclinical investigations. As radioresponse data within RadioGx has been fit to the LQ model, one could describe radioresponse through imputed survival across a range of dose levels (AUC) or at a specific dose level (e.g., SF2). There is currently no consensus regarding the optimal indicator for use across studies, with both AUC and SF2 frequently used (26–29). The use of SF2 as a radioresponse indicator has been bolstered by clinical observations that local tumor control following radiotherapy may be associated with SF2 measured from ex vivo tumor cells (30). Moreover, SF2 is thought to differentiate between radiosensitive and radioresistant cell types (31). However, there is currently insufficient evidence to support the routine use of SF2 or AUC when probing the molecular determinants of radioresponse.

We compared AUC and SF2 across all cell lines within RadioGx. The values were well correlated (R = 0.92; 95% CI, 0.90–0.93; P = 2.2e−16; Fig. 2A); the weakest correlations were observed among the most radioresistant cell lines, where cell death at higher doses likely contributes to the AUC value but by definition has no bearing on SF2 (Fig. 2B). We then asked whether the biological processes that govern these 2 radioresponse indicators are the same. To achieve this, we correlated the basal level gene expression data from the CCLE (32) with the radioresponse indicators (SF2 and AUC), and performed GSEA on the gene list ranked based on correlation estimates. For FDR <5%, 77 transcriptional pathways were enriched using AUC as the radioresponse indicator, out of which, 41 and 36 pathways were positively and negatively correlated with AUC, respectively (Supplementary File 1; Supplementary Fig. S5). Similarly, using SF2 as the radioresponse indicator, only 38 pathways were enriched, out of which, 19 were positively correlated with the SF2 value. All but 3 pathways enriched using SF2 were enriched using AUC (Fig. 2C and D).

The 17 pathways that were significantly correlated with radioresponse using the AUC indicator but not the SF2 indicator included biological processes known to impact radioresponse, suggesting stronger relevance for AUC. For instance, the NRF2-mediated oxidative stress response pathway was positively associated with AUC but not with SF2 (Supplementary File 1). In conditions of oxidative stress, such as following radiation, degradation of NRF2 is prevented, leading to its stabilization and translocation into the nucleus, where it activates expression of a wide variety of downstream antioxidant targets (33); this pathway has previously been described as contributing to intrinsic radioresistance (9, 34). In addition, progression through the cell cycle following radiation response is a known factor in determining cell survival vs. cell death via mitotic catastrophe. Three pathways directly related to cell-cycle progression [(i) cell cycle: G₂–M DNA damage checkpoint regulation; (ii) cell cycle: G₁–S checkpoint regulation; and (iii) cell-cycle control of chromosomal replication] were all seen exclusively when using AUC as the radioresponse indicator. Thus, as compared with SF2, AUC was able to capture more gene expression pathways putatively correlated with radioresponse. Taken together, our analyses reveal AUC and SF2 as related radioresponse indicators with AUC providing for a more comprehensive characterization of the biological processes underpinning radioresponse. As a result of these findings, we exclusively used AUC as the radioresponse indicator for subsequent analyses.

The LQ model can be used to estimate the dependence of cellular survival on radiotherapy fraction size and DNA repair. The |\alpha $| and |\beta $| values produced by the LQ model allow for comparisons among distinct cell lines or tumors, and in clinical practice the |\alpha /\beta $| ratio is used to predict cellular response to different radiotherapy fractionation schemes. Using the LQ model, we derived the |\alpha /\beta $| ratio for cancer cell lines within RadioGx. A wide range of |\alpha /\beta $| values were observed (Fig. 3A; median = 10.14; interquartile range = 4.49–28.07). As expected, the |\alpha $| component was strongly anticorrelated with AUC, whereas the |\beta $| component displayed no significant association with AUC (Fig. 3B). This result indicates that for the cell line data contained within RadioGx, dependence of cellular survival on radiotherapy fraction size is a distinct parameter that describes radioresponse and should therefore be considered alongside radiosensitivity (e.g., AUC or SF2) in preclinical investigations.

To examine the biological factors that underlie the differences between α, β, and AUC, we identified transcriptional pathways that were significantly associated with each radioresponse metric. For FDR <5%, we found 14 pathways commonly associated with all 3 metrics (Fig. 3C; Supplementary File 2). Supporting the biological relevance of these pathways, several known components of DNA damage response, signaling, and repair were represented among the 14 common pathways. For instance, pathways related to mismatch repair in eukaryotes, role of BRCA1 in DNA damage response, and cell-cycle control of chromosomal replication were each present. These results, which are consistent with fundamental tenets of radiobiology, suggest that analysis of large cell line resources within RadioGx could be performed to generate novel hypotheses and could contribute to preclinical biomarker discovery.

High fidelity DNA repair is critical for cellular survival following exposure to ionizing radiation. Among the various DNA repair processes, nonhomologous end joining (NHEJ) is the most important determinant of radiation response (35). Loss-of-function mutations within NHEJ pathway components confer radiosensitivity. We therefore used mutations in 19 genes implicated in NHEJ (Supplementary Table S3) as positive controls in the RadioGx datasets. Mutations were annotated according to predicted functional impact using VEP (12), and genes without predicted high functional impact (‘VEP-High') within our subset of cell lines were excluded. There was a strong overall association between ‘VEP High' mutations and radioresponse (Fig. 4A). Associations on an individual gene basis was limited by small sample size (Fig. 4B; Supplementary Fig. S6); however, VEP-High mutations in RAD50 were significantly associated with radioresponse. Overall, these results show that at the mutation level, the radioresponse data contained within RadioGx reflects known clinical radiosensitivity trends. Thus, the positive control analysis supports the robustness of the datasets contained within RadioGx and lends support for further studies that utilize these data.

By integrating radioresponse and molecular data, RadioGx is meant to enable new biological insights and predictions. To further demonstrate the utility of RadioGx for this purpose, we next extended the radiobiological modeling to incorporate the putative effects of oxygen availability in the tumor microenvironment on radioresponse (13).

Molecular oxygen is necessary to mediate the indirect effects of ionizing radiation to exert cell kill. Thus, cells become more resistant to radiation under oxygen-deficient conditions. We derived adjusted AUC values for the cancer cell lines within RadioGx at a range of oxygen partial pressures. As expected, reduced oxygen partial pressure resulted in a predicted increase in AUC (Fig. 5A). Cell lines from distinct cancer histologies displayed consistent increases in AUC under hypoxic conditions (P < 2.2e−16 for all, Wilcoxon test), but the magnitude of this increase differed between histologies (Fig. 5B). The largest and smallest median differences in AUC were observed for cancer cell lines from the breast and large intestine, respectively. These differences reflect a nonlinear relationship between oxygen availability and radioresponse that is dependent on |\alpha /\beta $|⁠.

Next, we evaluated the univariate association of gene expression levels measured under normoxic conditions with AUC values under normoxic and hypoxic conditions. For an FDR <5%, the numbers of genes that were significantly associated with radioresponse were 1,825 and 2,395 under normoxic and hypoxic conditions, respectively (Supplementary File 3). Moreover, 1,375 genes were negatively associated with radioresponse under normoxic condition but positively associated with radioresponse under hypoxic condition, and 471 genes were positively associated with radioresponse under normoxic condition but negatively associated with radioresponse under hypoxic condition (Supplementary Fig. S7). In keeping with these effects, we observed large changes in the ranking of strength of correlation of specific genes with radioresponse under oxic and hypoxic conditions (Fig. 5C). The gene with the greatest change, PPM1A, has been implicated in the regulation of cellular stress response and has previously been shown to have hypoxia-specific activity (36). WDR70, a gene with known roles in DNA double-strand break repair (37, 38), also displayed a large change in this analysis (Fig. 5C). One might hypothesize based on our results that WDR70 could have previously uncharacterized hypoxia-specific activities and/or expression; these findings warrant further investigation.

It is known that distinct tissues and tumor types respond differently to ionizing radiation exposure. Intrinsic radiosensitivity has been suggested as a major contributing factor to this differential response (10). We used RadioGx to interrogate radioresponse within tissue types (Fig. 6).

To examine the biological factors that may underlie suspected differences in radioresponse between tissue types, we identified 281 transcriptional pathways that were significantly associated with radioresponse within at least one tissue type (Fig. 6A; Supplementary File 4). Of these 281 pathways, 123 were statistically significant only in one tissue type (Supplementary Fig. S8). Overall, there were more statistically significant pathway associations with radiosensitivity than radioresistance (total across all tissue types: 437 and 226, respectively). Remarkably, we did not find any transcriptional pathways that were statistically significantly associated with radioresponse across all tissue types. We also observed variable |\alpha /\beta $| values among the tissue types within RadioGx (Fig. 6B), suggesting heterogeneity of DNA repair and dependence on radiotherapy fraction size.

Androgen receptor (AR) expression has emerged as a mediator of radioresistance in breast (10) and prostate cancer (39), but its effect on radioresponse in other tissue types is poorly understood. We used RadioGx to interrogate AR expression and its association with radioresponse across tissue types. As expected, among the available tissue types, breast cancer cell lines expressed AR at the highest levels and displayed the highest concordance between AR expression and radioresponse (Supplementary Fig. S9A). Notably, other tissue types including soft tissue, kidney, urinary tract, and stomach cancer showed similar concordance between AR expression and radioresponse (Supplementary Fig. S9B). Of these, soft tissue cell lines demonstrated the highest AR expression. Interestingly, AR expression has been found to be detectable in clinical samples (40) and is associated with radioresistance in rhabdomyosarcoma (41).

Datasets within RadioGx are standardized with regard to cell line annotations such that integrated analyses using other existing datasets can be easily conducted. For instance, our previously published tool, PharmacoGx (17), contains pharmacogenomic data from multiple studies and enables meta-analysis of pharmacogenomic data. We wished to identify categories of drugs with cytotoxic effects that correlate with radioresponse, so we interrogated RadioGx to compare cellular responses to ionizing radiation and chemotherapeutic agents (n = 545 distinct drugs). Drug responses were obtained from 480 cancer cell lines from the CTRPv2 pharmacogenomic dataset (Supplementary Table S1) that were in common between the datasets. We computed the correlation between drug response and radiation response across the cancer cell lines (Supplementary Fig. S10) and then classified drugs according to pharmacologic categories (i.e., by cellular targets and/or mechanisms of action). Drugs targeting the cytoskeleton, DNA replication, and mitosis displayed the strongest correlations with radioresponse (FDR < 5%; Fig. 7). Thymidylate synthetase inhibitors such as the known radiosensitizing drug, fluorourocil, also displayed cytotoxic effects that correlated with radioresponse but did not reach statistical significance. In addition to these anticipated and largely confirmatory findings, we also observed unexpected negative associations between radioresponse and cytotoxic effects of drugs targeting numerous cell signaling pathways (i.e., PI3K signaling, ERK MAPK signaling, WNT signaling, EGFR signaling, ABL signaling), although these were not statistically significant.

To date, the paradigm of precision medicine has primarily been applied to advanced incurable cancers. For early stage curable cancers for which radiotherapy is used with curative intent, there remains a need for more precise biologically-guided radiotherapy delivery. For instance, there are currently no clinically implemented molecular biomarkers that are predictive of radioresponse. This also extends to predictive insights into the response of tumors to other therapeutic agents that may be administered in combination with radiotherapy. Although molecular diagnostic tools are making their way into clinical practice in other settings, the lack of equivalent molecular indicators in the field of radiobiology has impeded translation in this domain (1, 7, 42).

Recently, large radioresponse and genomic datasets have been generated from hundreds of cancer cell lines, providing an opportunity to address this unmet need. We have developed RadioGx, an open-source software package that enables users to perform integrative analysis of radiogenomic datasets for preclinical evaluation of radioresponse determinants. RadioGx standardizes published nomenclature and annotations between datasets and integrates dose–response and molecular data. Because the RadioGx platform is developed on cell-based model systems, it suffers from certain limitations. Although cell lines have been derived from primary tumor patient samples, they evolved under different conditions temporally and may have accumulated culture-dependent genomic alterations. More importantly, these cell cultures may lack the endogenous hierarchical structure of tumors (comprising of stem cells, progenitor cells, differentiated cells, etc.) along with microenvironment components such as immune cells and other stromal components, all of which also play important roles in mediating tumor response to radiation therapy. We aimed to address many of these issues by using high quality data from large cell line factories with rigorous quality controls and by focusing on intrinsic radiosensitivity as an outcome.

We used RadioGx to evaluate the appropriateness of the 9-day viability assay for assessing radioresponse, the robustness of distinct radioresponse indicators, and the utility of applying established radiobiologic models to the data for novel hypothesis generation. We confirmed the findings from Yard and colleagues (10) that the 9-day viability assay, which is amenable to high-throughput processing and analysis, largely recapitulates the results of the more tedious gold standard clonogenic assay. We also were able to show that the radiation response data matched known biomarkers of response. We note that some prior putative intrinsic radiosensitivity gene expression signatures that were generated using cell line clonogenic survival data have failed validation using independent sets of cancer cell lines (27, 43), highlighting the need for robust reproducible methodologies for future studies. Moreover, we found that AUC derived from the LQ model might provide a more complete characterization of the biological processes underpinning radioresponse as compared with the dose-specific SF2 indicator, particularly for relatively radioresistant cell types. Based on our findings, we suggest that AUC should be the radioresponse indicator of choice for preclinical studies. Although we found that the LQ model fit the radioresponse data for the vast majority of cancer cell lines within RadioGx, a small subset was not amenable to LQ modeling. Although their use within the LQ framework should be approached with caution, it remains unclear whether these samples might reflect a distinct entity of tumors that do not abide the LQ formalism and may instead abide other useful models.

A major hurdle in the development of large-scale radioresponse datasets has been the technical and throughput challenges associated with the clonogenic assay. We demonstrated how existing data within RadioGx can be used to generate hypotheses and make predictions to inform future investigations. For instance, recognizing a dearth of large-scale radioresponse data under hypoxic conditions, we integrated radiobiological modeling with gene expression data from RadioGx, which allowed us to predict radioresponse under hypoxic conditions. Our findings suggest that the change in radioresponse under hypoxia is tissue-specific and that certain genes are either differentially associated with radioresponse under normoxic and hypoxic conditions or may have expression levels or activity that are regulated by oxygen tension; these hypotheses generated using RadioGx could be tested experimentally in future studies. In addition, by combining RadioGx with an existing pharmacogenomics analysis platform, we uncovered drugs with cytotoxic effects that are correlated or anticorrelated with radioresponse, suggestive of genomic/transcriptomic dependencies related to their mechanisms of action. We were able to confirm drug classes with therapeutic effects that overlap with ionizing radiation (e.g., mitotic inhibitors); moreover, this analysis proposed novel hypotheses regarding possible anticorrelated therapeutic effects with drugs targeting a number of cellular signaling pathways such as ABL and EGFR. Future studies may seek to examine whether members of these drug classes may make rational combination therapies with radiation as a result of reduced additive toxicity.

Predictive biomarkers of radioresponse could be applied across multiple cancer types (pan-cancer) or via a tissue-specific approach. Although conserved cellular processes are activated by radiation (44), different cell types have variable downstream responses and rates of survival (10). One factor that may have contributed to past decisions to focus on pan-cancer analysis is the limited amount of available cell-line radiosensitivity data; for instance, multiple putative radiosensitivity signatures have used radiation response data obtained from the NCI-60 cell lines (45–47). We showed that there is considerable variation in pathways associated with radiation response across tissue types within RadioGx, which supports a tissue-specific approach (48, 49). We envision that future investigations into tissue-specific radiation response biomarkers will be facilitated by the larger data sets curated in RadioGx, which will only expand with the inclusion of future data sets.

In summary, this study demonstrates the impact of combining radiogenomic datasets with established radiobiological models and other existing pharmacogenomic data. Future applications of RadioGx may include generation of biomarkers for intrinsic radiosensitivity and selection of novel combination therapies for preclinical testing. Thus, we envision that RadioGx will help to accelerate preclinical radiotherapeutic discovery pipelines and guide the selection of appropriate biological endpoints.

M.E. Abazeed reports a receiving commercial research grant from Bayer AG, Siemens Healthcare and has received speakers bureau honoraria from Bayer AG. S.V. Bratman reports receiving a commercial research grant from Nektar Therapeutics and is a co-inventor on patent licensed to Roche Molecular Diagnostics. No potential conflicts of interest were disclosed by the other authors.

Conception and design: V.S.K. Manem, M. Lambie, M.E. Abazeed, B. Haibe-Kains, S.V. Bratman

Development of methodology: V.S.K. Manem, M. Lambie, M.E. Abazeed, B. Haibe-Kains

Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): M.E. Abazeed, B. Haibe-Kains

Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): V.S.K. Manem, M. Lambie, I. Smith, M. Freeman, M. Koritzinsky, M.E. Abazeed, B. Haibe-Kains

Writing, review, and/or revision of the manuscript: V.S.K. Manem, M. Lambie, P. Smirnov, M. Koritzinsky, M.E. Abazeed, B. Haibe-Kains, S.V. Bratman

Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): I. Smith, P. Smirnov, V. Kofia, B. Haibe-Kains

Study supervision: B. Haibe-Kains, S.V. Bratman

Other (contributed to the development of the accompanying RadioGx Software Package): P. Smirnov

This work was supported by a grant from the V Foundation for Cancer Research (V2018-010) and from Canadian Institute of Health Research (PJT-162185). V.S.K. Manem was supported by the Terry Fox Research Institute. V.S.K. Manem and M. Freeman were supported the Canadian Institutes of Health Research. P. Smirnov was supported by Genome Canada and the Ontario Research Funds. S.V. Bratman and B. Haibe-Kains are supported by the Gattuso-Slaight Personalized Cancer Medicine Fund at the Princess Margaret Cancer Centre. M. Lambie was supported by a fellowship from STARS21. We also gratefully acknowledge the support from the Princess Margaret Cancer Foundation and the Princess Margaret Cancer Center Head & Neck Translational Program, with philanthropic funds from the Wharton Family, Joe's Team, and Gordon Tozer.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.