- 1Department of Urology, Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, China
- 2Department of Gastroenterology, Union Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, China
- 3Department of Geriatrics, Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, China
- 4Department of Oncology, Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, China
Clear cell renal cell carcinoma (ccRCC) is a prevalent urinary malignancy. Despite the recent development of better diagnostic tools and therapy, the five-year survival rate for individuals with advanced and metastatic ccRCC remains dismal. Unfortunately, ccRCC is less susceptible to radiation and chemotherapy. Consequently, targeted therapy and immunotherapy play a crucial role in the treatment of ccRCC. Enhancer RNAs (eRNAs) are noncoding RNAs transcribed by enhancers. Extensive research has shown that eRNAs are implicated in a variety of cancer signaling pathways. However, the biological functions of eRNAs have not been systematically investigated in ccRCC. In this study, we conducted a comprehensive investigation of the role of eRNAs in the onset and management of ccRCC. Patient prognosis-influencing eRNAs and target genes were chosen to construct a predictive signature. On the basis of the median riskscore, ccRCC patients were split into high- and low-risk subgroups. The prediction efficiency was assessed in several cohorts, and multi-omics analysis was carried out to investigate the differences and underlying mechanisms between the high- and low-risk groups. In addition, we investigated its potential to facilitate clinical treatment choices. The riskscore might be used to forecast a patient’s response to immunotherapy and targeted therapy, giving a revolutionary method for selecting treatment regimens with pinpoint accuracy.
Introduction
Renal cancer is one of the most common malignant tumors, and the estimated numbers of new deaths and new cases of renal cancer in the United States in 2021 were 13,780 and 76,080, respectively (1). Renal cell carcinoma (RCC) represents 90% of renal malignancies, with clear cell renal cell carcinoma (ccRCC) being the most prevalent subtype (2). Currently, the most effective therapy for ccRCC is surgery. However, following surgery, relapses and metastases are prevalent (3). Over 30% of ccRCC patients are reportedly affected by metastatic disease (4). Despite advancements in cancer detection and treatment, the 5-year survival rate for individuals with metastatic ccRCC remains around 20% (5). Since ccRCC does not respond well to radiation and chemotherapy, targeted therapy and immunotherapy become particularly important in the treatment of ccRCC (6, 7). Evidently, a customized strategy is crucial. On the basis of both the features of carcinomas and the situations of patients, the treatment regimen should be developed comprehensively.
Immunotherapies, particularly immune checkpoint blockade (ICB), have changed cancer therapy in recent years (8). Certain individuals with ccRCC exhibit astonishing clinical improvements when immune suppression is eliminated by ICB, and in the treatment of advanced ccRCC, checkpoint blockade in combination with other anticancer drugs is presently the first-line therapy (9). However, a substantial proportion of patients do not qualify for checkpoint blockade. This necessitates the use of reliable biomarkers to determine the optimal therapy strategy.
Enhancers are short clusters of regulatory DNA elements with transcription factor recognition sequences. They are either close to or far from the promoters of the target genes and control their expression (10, 11). Enhancers have been shown to play crucial roles in cancer formation and tumor response (12–14). Beyond their typical role of recruiting transcription factors to gene promoters, enhancer elements are also transcribed into noncoding RNAs known as enhancer RNAs (eRNAs) (15). eRNAs appear to be engaged in several cancer signaling pathways by modulating the expression of their target genes, including clinically actionable genes and immune checkpoints (16). By interacting with transcription factors, eRNAs may actually boost therapeutic anticancer therapies (16–18). eRNAs may serve as diagnostic indicators and tissue-specific treatment targets due to their tissue-specific expression across various cancer types. The risk model has been constructed in glioma (19), prostate cancer (20), and hepatocellular carcinoma (21). However, to the best of our knowledge, how eRNAs affect immune function and influence survival outcomes in ccRCC patients remains to be investigated.
In this study, we conducted an integrated analysis to investigate the functions of eRNAs in the development and management of ccRCC. Patient prognosis-influencing eRNAs and target genes were chosen to generate a predictive signature. Multi-omics analysis was performed to evaluate the differences and underlying processes between the high- and low-risk groups. Prediction accuracy was assessed in several cohorts. In addition, we explored the correlations of patient responsiveness to immunotherapy and targeted therapies with riskscore to evaluate its potentiality to facilitate treatment choices in clinical practice.
Materials and methods
Data collection and preprocessing
From the TCGA database, transcriptional data, somatic mutation data, and clinical features were obtained (Table S1). Then, we obtained the following datasets from the GEO database: GSE53757 (22) (N=144), GSE46699 (23) (N=130), GSE36895 (24) (N=76), and GSE22541 (25) (N=68). By using the “ComBat” algorithm of the “sva” package (26), we eliminated the batch effect and pooled the four datasets as a validation cohort. As further confirmation, we also retrieved RNA-seq data from the ICGC-RECA-EU database. The E-MTAB-1980 cohort (27) from the ArrayExpress database was also used to verify the predictive value.
Identification of prognostic eRNAs in ccRCC
PreSTIGE (Predicting Specific Tissue Interactions of Genes and Enhancers) was used to investigate lncRNAs produced from active tissue-specific enhancers and their putative target genes. We obtained lncRNA−target gene association data from a previous study for further analysis (28). To recognize the prognostic eRNAs, univariate Cox regression was used to determine candidate eRNAs and their respective target genes. After analysing the connection between eRNAs and target genes (Cor > 0.30, P-value< 0.01), we applied LASSO regression analysis for further screening. Finally, a multivariate Cox regression analysis (P-value< 0.05) was utilized to further identify eRNAs playing important roles in the prognosis of patients.
Establishment and evaluation of the prognostic signature
The prognostic model was constructed according to the results of multivariate Cox regression analysis. The riskscore for each patient was calculated by applying the following formula:
*β represents the regression coefficient, exp represents the expression value of each eRNA
On the basis of the median riskscore, ccRCC patients in the TCGA and E-MTAB-1980 cohorts were split into high- and low-risk subgroups. Using Kaplan-Meier curves, survival disparities between the high- and low-risk groups were displayed. Moreover, the “timeROC” package (29) was applied to assess the prognostic efficacy of the signature. To uncover the differences between the high- and low-risk groups, differential expression analysis was performed with the parameters of |logFC| > 1 and P< 0.05. Furthermore, GO and KEGG enrichment analyses were carried out by the R package “clusterProfiler” (30) to unearth the potential mechanism (FDR< 0.05 and q value< 0.05).
Correlation between riskscore and clinical parameters
To elucidate the influence of the riskscore on cancer development, the difference in riskscore across patients stratified by clinical parameters was calculated. We examined the percentage of patients in various stages between the high- and low-risk groups and investigated the predictive validity of the riskscore in early-stage and late ccRCC patients, given that pathologic stage is crucial for therapy selection. In addition, univariate and multivariate Cox regression analyses were conducted to assess the independence of the prognostic model from other clinicopathological characteristics.
Tumor mutation burden analysis
Tumor mutation burden (TMB) has been shown to be a predictive biomarker for determining whether cancer patients would react favorably to immune checkpoint inhibitors (31). Therefore, we investigated the relationship between TMB and riskscore in this study. Furthermore, patients were separated into four categories based on the median riskscore and TMB. The Kaplan-Meier method was used to evaluate differences in the survival distributions between groups.
Correlation between immune infiltration and riskscore
Using the “GSVA” R package, we compared the infiltration levels of 28 immune cells and immunological functioning across high- and low-risk groups (32). The Immune score and Stromal score of each sample were calculated by the “ESTIMATE” package (33). The correlations between immune infiltration and riskscore were then calculated after downloading the infiltration estimate data predicted by CIBERSORT, TIMER, xCell, quanTIseq, MCP-counter, and EPIC from TIMER 2.0 (34). In addition, the immunological subtypes of patients in the high- and low-risk groups were compared. Finally, a comparison was made between the riskscore and the human leukocyte antigen (HLA) gene family.
Benefits of the riskscore to aid treatment decision
Immunotherapy for cancer has altered the standard of treatment for certain patients with advanced disease (35). ICB has received attention in ccRCC patients as a promising anticancer treatment. Here, we first evaluated the differential expression of several immune checkpoint genes across high- and low-risk groups. Furthermore, we evaluated the correlations between CTLA4 and PD-1 expression and riskscore in different stages. Additionally, we retrieved the immunophenoscore (IPS) from The Cancer Immunome Database (TCIA) for ccRCC patients. The patient’s IPS was determined objectively by taking into account the four categories of immunogenicity-determining genes: effector cells, immunosuppressor cells, MHC molecules, and immune modulators (36). Immunogenicity is positively linked with higher IPS scores (37). In addition, we downloaded TCGA patient analysis findings from the Tumor Immune Dysfunction and Exclusion (TIDE). Three cell types believed to inhibit T-cell infiltration in tumors were studied across high- and low-risk groups: cancer-associated fibroblasts (CAFs), myeloid-derived suppressor cells (MDSCs), and the M2 subtype of tumor-associated macrophages (TAMs). Due to the absence of publicly available data on ccRCC cohorts adopting immunotherapy, we employed an immunotherapeutic cohort (IMvigor210 cohort) as a validation cohort (38). Comparisons of riskscore were made between patients with various clinical states after therapy.
Given that VEGFR-targeted therapy remains the first-line therapeutic choice for advanced ccRCC, we compared the susceptibility to sorafenib, sunitinib, pazopanib, and axitinib across high- and low-risk groups. The half maximal inhibitory concentrations (IC50) of these drugs were compared by applying the pRRophetic (39) package. In addition to the TCGA cohort, the findings were verified in the GEO and ICGC cohorts as well.
Statistical analysis
All of the analyses were conducted using RStudio 4.0.4. The Spearman method was used to calculate correlations. Two groups were compared using Student’s t test (mean with SD) and the Wilcoxon test. Multiple groups were analyzed using the Kruskal–Wallis test and one-way ANOVA. Statistical significance was defined as P-value< 0.05.
Results
Identification of prognosis-associated eRNAs and target genes in ccRCC
Figure 1A illustrates the research design for this study. We obtained 10066 regulatory relationships of eRNA and its target gene from previous research (28). After univariate Cox regression analysis and correlation analysis, 85 paired eRNAs and their target genes were identified (Table S2), including 69 eRNAs and 80 target genes. The corresponding relationship between the top 20 eRNAs associated with survival and their target genes is shown in Figure 1B. Moreover, we displayed the top 20 eRNAs or target genes with the highest mutation frequency (Figure 1C).
Figure 1 Flow chart and identification of prognostic eRNAs: (A) Diagram showing the methodology of this study. (B) Top 20 correlated eRNAs and their target genes relevant to survival. (C) The top 20 mutated eRNAs or target genes associated with patient prognosis.
Construction of the predictive signature and the prognostic value
The 69 eRNAs were filtered by using LASSO regression analysis (Figures 2A, B). The number of potential eRNAs was reduced to eight. After multivariate Cox regression analysis, five eRNAs (EMX2OS, GNG12-AS1, ZFHX4-AS1, AFG3L1P, LINC01271) remained to construct the predictive signature. The correlations between the five eRNAs are shown in Figure S1. The survival analysis and differential analysis between tumor and normal tissues performed well in terms of the prognostic genes (Figures 2C, D). Thus, we calculated the riskscore of each patient as follows: riskscore = (-0.1284 × EMX2OS expression) + (0.3209 × AFG3L1P expression) + (0.3546 × ZFHX4-AS1 expression) + (-0.8311 × GNG12-AS1 expression) + (0.9646 × LINC01271 expression) (Table 1).
Figure 2 Identification of prognosis-associated eRNAs: (A, B) LASSO regression analysis of the eRNAs correlated with OS in univariate Cox regression analysis. Survival analysis (C) and differential expression analysis (D) of the five eRNAs constructing the predictive signature. (E, F) GO and KEGG analyses between the high- and low-risk groups. (**p < 0.0, ****p < 0.0001).
After generating the riskscore, we conducted assessment and validation analyses. First, ccRCC patients from the TCGA-KIRC and E-MTAB-1980 cohorts were separated into high- and low-risk groups (Table S3), and the survival status of each patient in the cohort assessed by riskscore is shown in Figure 3A, B. The high-risk group had a greater mortality rate than the low-risk group, and Kaplan-Meier analysis indicated substantial differences in overall survival between the high-risk and low-risk groups (Figure 3C, D). Additionally, the predicted signature’s accuracy was also assessed using a time-dependent receiver operating characteristic curve analysis. The AUCs for 1-year, 2-year, 3-year, 4-year, and 5-year survival in the TCGA-KIRC dataset were 0.766, 0.714, 0.719, 0.710, and 0.743, respectively (Figure 3E). The E-MTAB-1980 cohort had AUCs of 0.732, 0.771, 0.730, 0.708, and 0.714, respectively, for survival rates at one, two, three, four, and five years (Figure 3F).
Figure 3 Validation of the predictive signature and performance analysis: (A, B) The riskscore distribution, survival status, and the expression of five eRNAs of patients in the TCGA cohort and E-MTAB-1980 cohort. (C, D) Kaplan-Meier survival curves of OS and the proportion of the survival rate between the high- and low-risk groups of the two datasets. (E, F) Time-dependent receiver operating characteristic analysis of the two datasets.
To investigate the biological differences between the high- and low-risk groups. GO and KEGG enrichment analyses were carried out to identify the biological functions of differentially expressed genes (Figures 2E, F). We discovered that the majority of the biological processes were immune system functions, such as humoral immune response, complement activation, and phagocytosis. Analysis of KEGG pathways indicated that differentially expressed genes were mostly engaged in neuroactive ligand-receptor interaction and cytokine-cytokine receptor interaction.
Correlation between the predictive signature and clinical parameters
To investigate the link between the predictive signature and clinical parameters, we evaluated the riskscore verified by clinical parameters, and the results indicated that the riskscore was markedly increased in patients with advanced pathologic stage, histologic grade and TNM stage (Figure 4A). Due to the importance of the pathologic stage in determining the optimal therapy, we examined the percentage of patients at various stages between the high- and low-risk groups. We discovered that the percentage of advanced patients in the high-risk category was dramatically increased (Figure 4B). (Figure 4B). We also evaluated the predictive validity of the riskscore in patients with early and advanced ccRCC (Figure 4C). We observed that the predictive signature could differentiate patient prognosis in both early and advanced stages well.
Figure 4 Association between the predictive signature and clinical parameters: (A) Correlations between riskscore and progression of ccRCC. (B) The proportion of patients in different stages between the high- and low-risk groups. (C) Kaplan-Meier survival curves of OS in early and advanced patients. (D) Univariate and (E) multivariate Cox regression analyses of correlations between the riskscore and clinical parameters.
To determine whether the predictive signature was an independent prognostic indicator in ccRCC, univariate and multivariate analyses were performed. The HR of the riskscore was 3.09 (95% CI: 2.51-3.81) and 2.24 (95% CI: 1.64-3.06) (Figures 4D, E), respectively, indicating that the riskscore was an independent prognostic factor in ccRCC. Interestingly, in the multivariate analyses, pathologic stage, histologic grade and TNM stage showed no significant differences when combined with riskscore, which also demonstrated a tight correlation between the predictive signature and clinical parameters.
Relationship between riskscore and TMB
TMB is emerging as a promising biomarker for predicting patients’ immune checkpoint inhibitor responses (40). Here, we assessed the relationship between riskscore and TMB. Patients in the high-risk group had a higher TMB than those in the low-risk group (Figure 5C), and patients with a higher TMB had a shorter PFS (Figure 5D). In addition, the correlation analysis revealed that the riskscore was positively related to TMB (Figure 5E). Then, patients were separated into four groups based on the riskscore and TMB. Patients in riskscorehighTMBhigh had the shortest median PFS, whereas those in riskscorelowTMBlow had the best prognosis (Figure 5F). Finally, a graphic representation of the mutation status of genes that had high mutation rates in the high- and low-risk categories was visualized (Figures 5A, B).
Figure 5 The correlation between TMB and riskscore: Genomic mutation status in the low-risk (A) and high-risk (B) groups is shown via a waterfall graphic. (C) The difference in TMB between groups with high- and low-risk. (D) The PFS survival curve between the high- and low-TMB groups. (E) Correlation analysis of riskscore and TMB. (F) Graphs showing Kaplan-Meier survival curves for four subgroups stratified by risk scores and TMB.
Riskscore was associated with ccRCC immune infiltration
We then explored the possible association between the riskscore and the tumor immunity of patients. In the high-risk and low-risk groups, there was a significant difference between the levels of most immune cells infiltrated (Figure S2) and the functions of these cells (Figure S3). Furthermore, we compared the levels of immune cell infiltration with the riskscore, and the findings indicated that the riskscore was closely associated with the majority of the immune cells (Figure 6B). Immune subtypes are characteristics of the tumor microenvironment (TME) that cut beyond conventional cancer classifications to generate groups and imply that some therapeutic methods may be independent of histologic type. We retrieved immune subtypes of TCGA from a previous work (41) and compared the prevalence of various immunological subtypes in high- and low-risk groups. The findings demonstrated striking changes in immunological properties between the high- and low-risk groups (Figure 6A). As shown in Figure 6C, the Immune Score of the high-risk patients was significantly higher than that of the low-risk patients. In addition, the expression levels of the majority of HLA-related genes were elevated in the high-risk groups (Figure 6D).
Figure 6 Relationship between immune infiltration and riskscore: (A) The proportion of patients with different immune subtypes between the high- and low-risk groups. (B) The correlation of the riskscore and immune cell infiltration levels. (C) Comparison of Immune score and Stromal score between the high- and low-risk subgroups. (D) The expression of the HLA gene family between the high- and low-risk groups. (*p < 0.05, **p < 0.01, ***p < 0.001). ns, no significance.
Riskscore predicts therapeutic benefits
Anticancer immunotherapies using inhibitors of immune checkpoints have been developed as novel treatment regimens (42). The tumor microenvironment (TME) has been shown to be closely connected to tumor growth and metastasis (43) and may dampen the treatment response, hence influencing the clinical outcome (44). To further investigate the relationships between immunotherapy response and riskscore, we analyzed the association between riskscore and immune checkpoint inhibitor genes. The expression of ICI genes, including CTLA4 and PD-1, was elevated in the high-risk groups relative to the low-risk groups (Figure 7C). The correlation analysis also suggested that the riskscore was positively correlated with immune checkpoint gene expression (Figure 7D). Considering that immunotherapies are typically prescribed to patients with advanced ccRCC, we also estimated the correlation between riskscore and immune checkpoint gene expression across different stages. It was observed that the correlation coefficient was increased in advanced status, especially in stage iv (Figure 7E). IPS analysis was used to assess the immunogenicity of the two prognostic groupings. Patients in the high-risk group had higher IPS, IPS-CTLA4, IPS-PD1, and IPS-PD1-CTLA4 scores, suggesting that immunotherapy may have a greater response in this group (Figure 7B). CAFs, MDSCs and M2-TAMs have been shown to inhibit T-cell infiltration into cancers, and we calculated their prevalence in prostate cancer. High-risk individuals had considerably fewer CAFs and M2-TAMs (Figure 7A). However, the high-risk group had a larger percentage of MDSCs. Finally, we examined the riskscore performance in the IMvigor210 cohort. We compared the riskscore of patients with different immunotherapy responses. Patients with the greatest riskscore were those who had a complete response (Figure 7F). According to these findings, patients in the high-risk category may benefit from immunotherapy.
Figure 7 Predictive value of riskscore for immunotherapy: (A) The infiltration levels of CAFs, MDSCs, and M2-TAMs between the high- and low-risk groups. (B) The correlation between IPS and riskscore. (C) The expression of immune checkpoint genes in subgroups with high- and low-risk. (D) The association between riskscore and immune checkpoint gene expression. (E) The correlation between CTLA4, PD-1 expression and the riskscore in different pathologic stages. (F) The difference in riskscore among individuals in the IMvigor210 cohort performing different responses to immunotherapy. (*p < 0.05, **p < 0.01, ***p < 0.001). ns, no significance.
Given that VEGFR-targeted therapy remains the first-line therapeutic choice for advanced ccRCC, we compared the sensitivity of sorafenib, sunitinib, pazopanib, and axitinib across high- and low-risk groups. In the TCGA cohort, we discovered that the high-risk group was more likely to respond to sunitinib and axitinib, while the low-risk group responded better to sorafenib and pazopanib (Figure 8A). The validation cohort was subsequently used to verify these findings. The GEO cohort confirmed the sensitivities to sunitinib and axitinib for patients with a high riskscore (Figure 8B), while the ICGC cohort validated the sensitivities to sunitinib and resistance to sorafenib (Figure 8C) (Table S4).
Figure 8 Riskscore predicts therapeutic benefits of targeted therapies: The estimated IC50 values of sunitinib, axitinib, sorafenib and pazopanib in the TGCA cohort (A), GEO cohort (B), and ICGC cohort (C). (*p < 0.05, **p < 0.01, ***p < 0.001). ns, no significance.
Discussion
The most prevalent type of kidney cancer is ccRCC. Due to the absence of accurate diagnostic biomarkers and the constraints of efficient screening detection, about 30% of RCC patients first present with metastatic disease (45). Recurrence occurs in around 30% of individuals after full excision of the initial tumor (46). Traditional chemotherapy and radiation therapy are generally unsuccessful against all subtypes of RCC (47). Lack of susceptibility to chemotherapy and radiation therapy has motivated study into other therapeutic methods. In recent years, the introduction of targeted therapies, such as multi-targeted tyrosine kinase inhibitors and mTOR inhibitors, has been a significant advancement in ccRCC treatment (7).
Recent studies indicates that the immune microenvironment generated by tumor immune cells may control the development of cancer (48–50). Immune checkpoint inhibitors have emerged as viable therapy choices for advanced ccRCC in recent years (51–53). Immunotherapies such as nivolumab have already been included in metastatic ccRCC therapy regimens (2). However, in the setting of high tumor heterogeneity in terms of cell biological characteristics and genetics, metastatic or advanced RCC patients still have a poor prognosis due to the absence of effective therapeutic approaches that may produce long-lasting responses (54, 55). Therefore, identifying and validating biomarkers are essential for enhancing first-line selection and treatment sequences (56).
Recent breakthroughs in genome-wide investigations have shown that the dysregulation of distal gene regulatory elements, such as enhancers, occurs in a number of pathophysiological diseases (57). They generate eRNAs, which directly contribute to tumorigenesis (58, 59). Through their ability to resolve intratumor heterogeneity with specificity of cell-type enhancers, eRNAs provide additional explanations for cancer phenotypes beyond those provided by mRNA expression (60). In addition to regulating gene expression in cancer, eRNAs also maintain genome stability, which is associated with chromosome rearrangements and genome instability (61). Moreover, it has been shown that eRNAs are expressed across human cancer tissues, indicating that they could be useful as biomarkers and therapeutic targets (62). Cancer resistance is a well-defined phenomenon that arises when cancer cells develop tolerance to a certain therapeutic dosage. Several molecular pathways, including genetic or epigenetic changes and enhanced drug efflux, contribute to therapeutic resistance (57). Deregulation of enhancer transcription has been linked to various malignancies and their treatments during the last decade (63, 64). However, the prognostic significance of eRNAs in ccRCC is yet unknown.
In this study, active tissue-specific enhancers and their predicted potential target genes were obtained from previous research. The prognostic eRNA and its target gene were filtered using univariate Cox regression analysis and correlation analysis. Next, LASSO regression analysis and multivariate Cox regression analysis were used to build the predictive signature. Finally, five eRNAs (EMX2OS, GNG12-AS1, ZFHX4-AS1, AFG3L1P, and LINC01271) remained to construct the predictive signature. The prognostic value of riskscore was demonstrated in both the TCGA and E-MTAB-1980 cohorts. In addition, the clinical correlation study revealed a positive association between the riskscore and T stage, N stage, M stage, histologic grade, and pathologic stage. In the high-risk group, the percentage of advanced patients was much higher. GO and KEGG analyses were carried out to identify the likely underlying mechanism. We observed that the biological processes enriched by differentially expressed genes were highly connected with the immune response, suggesting that there may be substantial variations between the high- and low-risk groups in terms of immunological features.
Multiple recent studies have shown that TMB correlates with immunotherapy response because it reflects the total neoantigen load (65, 66). Low TMB is a good prognostic marker but predicts the adverse predictive effectiveness of ICI therapy (65). Herein, we discovered a substantial positive association between the riskscore and TMB. Patients in the high-risk group had a higher TMB, indicating a more favorable treatment response. In addition, the expression levels of the majority of HLA-related and immune checkpoint genes were considerably elevated in the high-risk groups. Interestingly, we found that the correlation coefficient was higher in individuals with advanced disease, particularly those diagnosed with stage IV. Therefore, high-risk individuals with advanced disease may have a higher expression of CTLA4 and PD-1. Higher scores on the IPS are related to increased immunogenicity (36). In the high-risk group, the IPS, IPS-CTLA4, IPS-PD1, and IPS-PD1-CTLA4 values were considerably higher. Thus, patients in the high-risk group may exhibit a stronger immunotherapy response. Due to the rarity of open-access data on ccRCC cohorts receiving immunotherapy, we utilized IMvigor210 cohort patients for preliminary validation. Consistent with the aforementioned hypotheses, we discovered that patients who had a complete response had the highest risk score.
Given the importance of targeted therapy for patients with advanced or metastatic ccRCC, we compared the IC50 values of the first-line agents between the high- and low-risk groups. In the TCGA cohort, we found that the high-risk group was more likely to be responsive to sunitinib and axitinib, while the low-risk group was more sensitive to sorafenib and pazopanib. These findings were confirmed in the GEO and ICGC cohorts. Obviously, the combination of targeted therapies and personalized immunotherapy is more suited for patients with advanced ccRCC, since it is unreasonable to expect that a single treatment can significantly improve their prognosis. In fact, according to the recommendations of the European Association of Urology, the combination of nivolumab and ipilimumab or pembrolizumab and axitinib are the new first-line therapy for patients at intermediate and poor risk (67). It has been observed that pembrolizumab plus axitinib therapy resulted in substantially higher overall survival and progression-free survival than sunitinib treatment among patients with advanced renal-cell carcinoma who had not been treated before (68). Another trial also recognized that progression-free survival was considerably longer with avelumab plus axitinib compared to sunitinib in patients with advanced renal cell carcinoma (69). Actually, targeted metastatic RCC treatments have been demonstrated to have immunomodulatory effects, such as raising tumor cell antigenicity and encouraging T-cell infiltration (70). These results sparked an interest in exploring the possibility of combining targeted antiangiogenic medicines with immunotherapies to maximize any potential synergies (71). In this study, it was shown that high-risk individuals respond better to immunotherapy and are more sensitive to sunitinib and axitinib. On the basis of these results, we hypothesize that the riskscore might be used to aid in the selection of treatment regimens in practical practice. The combination of sunitinib, axitinib, and immunotherapy may be beneficial for high-risk individuals, while sorafenib and pazopanib are more effective for low-risk patients.
Although we constructed a predictive signature to aid treatment decision in ccRCC, several limitations to this study need to be acknowledged. First, eRNA is a subtype of noncoding RNA produced by enhancers, and it is not or is rarely detected in many RNA-seq datasets. Therefore, validation in the GEO cohort and E-MTAB-1980 cohort could not include all the eRNAs in the signature, which may reduce the efficiency of the signature. In addition, due to the paucity of open-access immunotherapy cohort data for ccRCC, a preliminary validation was conducted in the IMvigor210 cohort of bladder cancer. The potential of this signature to predict the immunotherapy response of patients in ccRCC groups receiving treatment requires additional confirmation.
Conclusions
In summary, we conducted an exhaustive investigation of the involvement of eRNAs in the development and treatment of ccRCC. The prognosis-influencing eRNAs and target genes were chosen to generate a predictive signature. Multi-omics analysis was performed to examine the differences and underlying mechanisms between the high- and low-risk groups, and the predictive accuracy was assessed in multiple cohorts. In addition, relationships between patient response to immunotherapy, targeted therapies, and riskscore were assessed. Our findings indicated that the combination of sunitinib, axitinib, and immunotherapy may be beneficial for high-risk individuals, while sorafenib and pazopanib are more effective for low-risk patients.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.
Author contributions
QL and XX designed the study and wrote the manuscript. CL, YL and BL provided administrative support. BC, GS, KZ, BNL and JM analyzed the data and modified the manuscript. All authors contributed to the article and approved the submitted version.
Funding
This work was supported by grants from the National Natural Science Foundation of China (grant number 81902619 and 82173068) and Wuhan Shuguang Project (grant number 2022020801020447).
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2022.964838/full#supplementary-material
References
1. Siegel RL, Miller KD, Fuchs HE, Jemal A. Cancer statistics, 2021. CA Cancer J Clin (2021) 71:7–33. doi: 10.3322/caac.21654
2. Hsieh JJ, Purdue MP, Signoretti S, Swanton C, Albiges L, Schmidinger M, et al. Renal cell carcinoma. Nat Rev Dis Primers (2017) 3:17009. doi: 10.1038/nrdp.2017.9
3. Cancer Genome Atlas Research Network. Comprehensive molecular characterization of clear cell renal cell carcinoma. Nature (2013) 499:43–9. doi: 10.1038/nature12222
4. Terrematte P, Andrade DS, Justino J, Stransky B, de Araújo DSA, Dória Neto AD. A novel machine learning 13-gene signature: Improving risk analysis and survival prediction for clear cell renal cell carcinoma patients. Cancers (Basel) (2022) 14:2111. doi: 10.3390/cancers14092111
5. Coppin C, Kollmannsberger C, Le L, Porzsolt F, Wilt TJ. Targeted therapy for advanced renal cell cancer (RCC): a cochrane systematic review of published randomised trials. BJU Int (2011) 108:1556–63. doi: 10.1111/j.1464-410X.2011.10629.x
6. Tacconi EMC, Tuthill M, Protheroe A. Review of adjuvant therapies in renal cell carcinoma: Evidence to date. Onco Targets Ther (2020) 13:12301–16. doi: 10.2147/OTT.S174149
7. Makhov P, Joshi S, Ghatalia P, Kutikov A, Uzzo RG, Kolenko VM. Resistance to systemic therapies in clear cell renal cell carcinoma: Mechanisms and management strategies. Mol Cancer Ther (2018) 17:1355–64. doi: 10.1158/1535-7163.MCT-17-1299
8. Zou W, Wolchok JD, Chen L. PD-L1 (B7-H1) and PD-1 pathway blockade for cancer therapy: Mechanisms, response biomarkers, and combinations. Sci Transl Med (2016) 8:328rv4. doi: 10.1126/scitranslmed.aad7118
9. Díaz-Montero CM, Rini BI, Finke JH. The immunology of renal cell carcinoma. Nat Rev Nephrol (2020) 16:721–35. doi: 10.1038/s41581-020-0316-3
10. Li W, Notani D, Rosenfeld MG. Enhancers as non-coding RNA transcription units: recent insights and future perspectives. Nat Rev Genet (2016) 17:207–23. doi: 10.1038/nrg.2016.4
11. Natoli G, Andrau J-C. Noncoding transcription at enhancers: general principles and functional models. Annu Rev Genet (2012) 46:1–19. doi: 10.1146/annurev-genet-110711-155459
12. Bahr C, von Paleske L, Uslu VV, Remeseiro S, Takayama N, Ng SW, et al. A myc enhancer cluster regulates normal and leukaemic haematopoietic stem cell hierarchies. Nature (2018) 553:515–20. doi: 10.1038/nature25193
13. Takeda DY, Spisák S, Seo J-H, Bell C, O’Connor E, Korthauer K, et al. A somatically acquired enhancer of the androgen receptor is a noncoding driver in advanced prostate cancer. Cell (2018) 174:422–32.e13. doi: 10.1016/j.cell.2018.05.037
14. Mack SC, Pajtler KW, Chavez L, Okonechnikov K, Bertrand KC, Wang X, et al. Therapeutic targeting of ependymoma as informed by oncogenic enhancer profiling. Nature (2018) 553:101–5. doi: 10.1038/nature25169
15. Lee J-H, Xiong F, Li W. Enhancer RNAs in cancer: regulation, mechanisms and therapeutic potential. RNA Biol (2020) 17:1550–9. doi: 10.1080/15476286.2020.1712895
16. Zhang Z, Lee J-H, Ruan H, Ye Y, Krakowiak J, Hu Q, et al. Transcriptional landscape and clinical utility of enhancer RNAs for eRNA-targeted therapy in cancer. Nat Commun (2019) 10:4562. doi: 10.1038/s41467-019-12543-5
17. Liang J, Zhou H, Gerdt C, Tan M, Colson T, Kaye KM, et al. Epstein-Barr Virus super-enhancer eRNAs are essential for MYC oncogene expression and lymphoblast proliferation. Proc Natl Acad Sci U.S.A. (2016) 113:14121–6. doi: 10.1073/pnas.1616697113
18. Ye Y, Xiang Y, Ozguc FM, Kim Y, Liu C-J, Park PK, et al. The genomic landscape and pharmacogenomic interactions of clock genes in cancer chronotherapy. Cell Syst (2018) 6:314–28.e2. doi: 10.1016/j.cels.2018.01.013
19. Tian W, Chen K, Yan G, Han X, Liu Y, Zhang Q, et al. A novel prognostic tool for glioma based on enhancer RNA-regulated immune genes. Front Cell Dev Biol (2022) 9:798445. doi: 10.3389/fcell.2021.798445
20. Fan S, Wang Z, Zhao L, Zhao C, Yuan D, Wang J. A robust prognostic gene signature based on eRNAs-driven genes in prostate cancer. Front Genet (2021) 12:676845. doi: 10.3389/fgene.2021.676845
21. Cai S, Hu X, Chen R, Zhang Y. Identification and validation of an immune-related eRNA prognostic signature for hepatocellular carcinoma. Front Genet (2021) 12:657051. doi: 10.3389/fgene.2021.657051
22. von Roemeling CA, Radisky DC, Marlow LA, Cooper SJ, Grebe SK, Anastasiadis PZ, et al. Neuronal pentraxin 2 supports clear cell renal cell carcinoma by activating the AMPA-selective glutamate receptor-4. Cancer Res (2014) 74:4796–810. doi: 10.1158/0008-5472.CAN-14-0210
23. Eckel-Passow JE, Serie DJ, Bot BM, Joseph RW, Cheville JC, Parker AS. ANKS1B is a smoking-related molecular alteration in clear cell renal cell carcinoma. BMC Urol (2014) 14:14. doi: 10.1186/1471-2490-14-14
24. Peña-Llopis S, Vega-Rubín-de-Celis S, Liao A, Leng N, Pavía-Jiménez A, Wang S, et al. BAP1 loss defines a new class of renal cell carcinoma. Nat Genet (2012) 44:751–9. doi: 10.1038/ng.2323
25. Wuttig D, Zastrow S, Füssel S, Toma MI, Meinhardt M, Kalman K, et al. CD31, EDNRB and TSPAN7 are promising prognostic markers in clear-cell renal cell carcinoma revealed by genome-wide expression analyses of primary tumors and metastases. Int J Cancer (2012) 131:E693–704. doi: 10.1002/ijc.27419
26. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics (2012) 28:882–3. doi: 10.1093/bioinformatics/bts034
27. Sato Y, Yoshizato T, Shiraishi Y, Maekawa S, Okuno Y, Kamura T, et al. Integrated molecular analysis of clear-cell renal cell carcinoma. Nat Genet (2013) 45:860–7. doi: 10.1038/ng.2699
28. Vučićević D, Corradin O, Ntini E, Scacheri PC, Ørom UA. Long ncRNA expression associates with tissue-specific enhancers. Cell Cycle (2015) 14:253–60. doi: 10.4161/15384101.2014.977641
29. Blanche P, Dartigues J-F, Jacqmin-Gadda H. Estimating and comparing time-dependent areas under receiver operating characteristic curves for censored event times with competing risks. Stat Med (2013) 32:5381–97. doi: 10.1002/sim.5958
30. Yu G, Wang L-G, Han Y, He Q-Y. clusterProfiler: an r package for comparing biological themes among gene clusters. OMICS (2012) 16:284–7. doi: 10.1089/omi.2011.0118
31. Merino DM, McShane LM, Fabrizio D, Funari V, Chen S-J, White JR, et al. Establishing guidelines to harmonize tumor mutational burden (TMB): in silico assessment of variation in TMB quantification across diagnostic platforms: phase I of the friends of cancer research TMB harmonization project. J Immunother Cancer (2020) 8:e000147. doi: 10.1136/jitc-2019-000147
32. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinf (2013) 14:7. doi: 10.1186/1471-2105-14-7
33. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun (2013) 4:2612. doi: 10.1038/ncomms3612
34. Li T, Fu J, Zeng Z, Cohen D, Li J, Chen Q, et al. TIMER2.0 for analysis of tumor-infiltrating immune cells. Nucleic Acids Res (2020) 48:W509–14. doi: 10.1093/nar/gkaa407
35. Heinhuis KM, Ros W, Kok M, Steeghs N, Beijnen JH, Schellens JHM. Enhancing antitumor response by combining immune checkpoint inhibitors with chemotherapy in solid tumors. Ann Oncol (2019) 30:219–35. doi: 10.1093/annonc/mdy551
36. Charoentong P, Finotello F, Angelova M, Mayer C, Efremova M, Rieder D, et al. Pan-cancer immunogenomic analyses reveal genotype-immunophenotype relationships and predictors of response to checkpoint blockade. Cell Rep (2017) 18:248–62. doi: 10.1016/j.celrep.2016.12.019
37. Gui C-P, Wei J-H, Chen Y-H, Fu L-M, Tang Y-M, Cao J-Z, et al. A new thinking: extended application of genomic selection to screen multiomics data for development of novel hypoxia-immune biomarkers and target therapy of clear cell renal cell carcinoma. Brief Bioinform (2021) 22:bbab173. doi: 10.1093/bib/bbab173
38. Mariathasan S, Turley SJ, Nickles D, Castiglioni A, Yuen K, Wang Y, et al. TGFβ attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature (2018) 554:544–8. doi: 10.1038/nature25501
39. Geeleher P, Cox N, Huang RS. pRRophetic: an r package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS One (2014) 9:e107468. doi: 10.1371/journal.pone.0107468
40. Chan TA, Yarchoan M, Jaffee E, Swanton C, Quezada SA, Stenzinger A, et al. Development of tumor mutation burden as an immunotherapy biomarker: utility for the oncology clinic. Ann Oncol (2019) 30:44–56. doi: 10.1093/annonc/mdy495
41. Thorsson V, Gibbs DL, Brown SD, Wolf D, Bortone DS, Ou Yang T-H, et al. The immune landscape of cancer. Immunity (2018) 48:812–30.e14. doi: 10.1016/j.immuni.2018.03.023
42. O’Donnell JS, Teng MWL, Smyth MJ. Cancer immunoediting and resistance to T cell-based immunotherapy. Nat Rev Clin Oncol (2019) 16:151–67. doi: 10.1038/s41571-018-0142-8
43. Brassart-Pasco S, Brézillon S, Brassart B, Ramont L, Oudart J-B, Monboisse JC. Tumor microenvironment: Extracellular matrix alterations influence tumor progression. Front Oncol (2020) 10:397. doi: 10.3389/fonc.2020.00397
44. Wu T, Dai Y. Tumor microenvironment and therapeutic response. Cancer Lett (2017) 387:61–8. doi: 10.1016/j.canlet.2016.01.043
45. Motzer RJ, Bukowski RM, Figlin RA, Hutson TE, Michaelson MD, Kim ST, et al. Prognostic nomogram for sunitinib in patients with metastatic renal cell carcinoma. Cancer (2008) 113:1552–8. doi: 10.1002/cncr.23776
46. Nerich V, Hugues M, Paillard MJ, Borowski L, Nai T, Stein U, et al. Clinical impact of targeted therapies in patients with metastatic clear-cell renal cell carcinoma. Onco Targets Ther (2014) 7:365–74. doi: 10.2147/OTT.S56370
47. Goyal R, Gersbach E, Yang XJ, Rohan SM. Differential diagnosis of renal tumors with clear cytoplasm: clinical relevance of renal tumor subclassification in the era of targeted therapies and personalized medicine. Arch Pathol Lab Med (2013) 137:467–80. doi: 10.5858/arpa.2012-0085-RA
48. Gajewski TF, Schreiber H, Fu Y-X. Innate and adaptive immune cells in the tumor microenvironment. Nat Immunol (2013) 14:1014–22. doi: 10.1038/ni.2703
49. Hinshaw DC, Shevde LA. The tumor microenvironment innately modulates cancer progression. Cancer Res (2019) 79:4557–66. doi: 10.1158/0008-5472.CAN-18-3962
50. Pitt JM, Marabelle A, Eggermont A, Soria J-C, Kroemer G, Zitvogel L. Targeting the tumor microenvironment: removing obstruction to anticancer immune responses and immunotherapy. Ann Oncol (2016) 27:1482–92. doi: 10.1093/annonc/mdw168
51. Braun DA, Hou Y, Bakouny Z, Ficial M, Sant’ Angelo M, Forman J, et al. Interplay of somatic alterations and immune infiltration modulates response to PD-1 blockade in advanced clear cell renal cell carcinoma. Nat Med (2020) 26:909–18. doi: 10.1038/s41591-020-0839-y
52. Miao D, Margolis CA, Gao W, Voss MH, Li W, Martini DJ, et al. Genomic correlates of response to immune checkpoint therapies in clear cell renal cell carcinoma. Science (2018) 359:801–6. doi: 10.1126/science.aan5951
53. Au L, Hatipoglu E, Robert de Massy M, Litchfield K, Beattie G, Rowan A, et al. Determinants of anti-PD-1 response and resistance in clear cell renal cell carcinoma. Cancer Cell (2021) 39:1497–518.e11. doi: 10.1016/j.ccell.2021.10.001
54. Mitchell AP, Hirsch BR, Harrison MR, Abernethy AP, George DJ. Deferred systemic therapy in patients with metastatic renal cell carcinoma. Clin Genitourin Cancer (2015) 13:e159–166. doi: 10.1016/j.clgc.2014.12.017
55. Park I, Lee J-L, Ahn J-H, Lee D-H, Lee K-H, Jeong IG, et al. Active surveillance for metastatic or recurrent renal cell carcinoma. J Cancer Res Clin Oncol (2014) 140:1421–8. doi: 10.1007/s00432-014-1680-9
56. Atkins MB, Tannir NM. Current and emerging therapies for first-line treatment of metastatic clear cell renal cell carcinoma. Cancer Treat Rev (2018) 70:127–37. doi: 10.1016/j.ctrv.2018.07.009
57. Adhikary S, Roy S, Chacon J, Gadad SS, Das C. Implications of enhancer transcription and eRNAs in cancer. Cancer Res (2021) 81:4174–82. doi: 10.1158/0008-5472.CAN-20-4010
58. Murakawa Y, Yoshihara M, Kawaji H, Nishikawa M, Zayed H, Suzuki H, et al. Enhanced identification of transcriptional enhancers provides mechanistic insights into diseases. Trends Genet (2016) 32:76–88. doi: 10.1016/j.tig.2015.11.004
59. Chen H, Li C, Peng X, Zhou Z, Weinstein JN, Cancer Genome Atlas Research Network, et al. A pan-cancer analysis of enhancer expression in nearly 9000 patient samples. Cell (2018) 173:386–99.e12. doi: 10.1016/j.cell.2018.03.027
60. Chen H, Liang H. A high-resolution map of human enhancer RNA loci characterizes super-enhancer activities in cancer. Cancer Cell (2020) 38:701–15.e5. doi: 10.1016/j.ccell.2020.08.020
61. Aguilera A, García-Muse T. R loops: from transcription byproducts to threats to genome stability. Mol Cell (2012) 46:115–24. doi: 10.1016/j.molcel.2012.04.009
62. Sartorelli V, Lauberth SM. Enhancer RNAs are an important regulatory layer of the epigenome. Nat Struct Mol Biol (2020) 27:521–8. doi: 10.1038/s41594-020-0446-0
63. Schmitt AM, Chang HY. Long noncoding RNAs in cancer pathways. Cancer Cell (2016) 29:452–63. doi: 10.1016/j.ccell.2016.03.010
64. Hamdan FH, Johnsen SA. Perturbing enhancer activity in cancer therapy. Cancers (Basel) (2019) 11:E634. doi: 10.3390/cancers11050634
65. Liu L, Bai X, Wang J, Tang X-R, Wu D-H, Du S-S, et al. Combination of TMB and CNA stratifies prognostic and predictive responses to immunotherapy across metastatic cancer. Clin Cancer Res (2019) 25:7413–23. doi: 10.1158/1078-0432.CCR-19-0558
66. Rizvi NA, Hellmann MD, Snyder A, Kvistborg P, Makarov V, Havel JJ, et al. Cancer immunology. mutational landscape determines sensitivity to PD-1 blockade in non-small cell lung cancer. Science (2015) 348:124–8. doi: 10.1126/science.aaa1348
67. Motzer RJ, Tannir NM, McDermott DF, Arén Frontera O, Melichar B, Choueiri TK, et al. Nivolumab plus ipilimumab versus sunitinib in advanced renal-cell carcinoma. N Engl J Med (2018) 378:1277–90. doi: 10.1056/NEJMoa1712126
68. Rini BI, Plimack ER, Stus V, Gafanov R, Hawkins R, Nosov D, et al. Pembrolizumab plus axitinib versus sunitinib for advanced renal-cell carcinoma. N Engl J Med (2019) 380:1116–27. doi: 10.1056/NEJMoa1816714
69. Motzer RJ, Penkov K, Haanen J, Rini B, Albiges L, Campbell MT, et al. Avelumab plus axitinib versus sunitinib for advanced renal-cell carcinoma. N Engl J Med (2019) 380:1103–15. doi: 10.1056/NEJMoa1816047
70. Das R, Verma R, Sznol M, Boddupalli CS, Gettinger SN, Kluger H, et al. Combination therapy with anti-CTLA-4 and anti-PD-1 leads to distinct immunologic changes. vivo. J Immunol (2015) 194:950–9. doi: 10.4049/jimmunol.1401686
Keywords: clear cell renal cell carcinoma, drug sensitivity, enhancer RNA, immunotherapy, treatment decision
Citation: Li Q, Xiao X, Chen B, Song G, Zeng K, Li B, Miao J, Liu C, Luan Y and Liu B (2022) A predictive signature based on enhancer RNA associates with immune infiltration and aids treatment decision in clear cell renal cell carcinoma. Front. Oncol. 12:964838. doi: 10.3389/fonc.2022.964838
Received: 22 June 2022; Accepted: 26 September 2022;
Published: 12 October 2022.
Edited by:
Yin Sun, University of Rochester Medical Center, United StatesReviewed by:
Jesús Beltrán-García, University of California, San Diego, United StatesZeyan Li, Shandong University, China
Copyright © 2022 Li, Xiao, Chen, Song, Zeng, Li, Miao, Liu, Luan and Liu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Bo Liu, Ym9saXU4ODhAaG90bWFpbC5jb20=; Yang Luan, eWFuZ2x1YW5AdGpoLnRqbXUuZWR1LmNu; Chaofan Liu, bGNmMjAxMDE5OTVAc2luYS5jb20=