Molecular Characteristics, Oncogenic Roles, and Relevant Immune and Pharmacogenomic Features of EVA1B in Colorectal Cancer

Objective: EVA1B, a protein coding gene, is a critical paralog of EVA1A gene. Herein, our study was conducted to investigate the role of EVA1B in colorectal cancer (CRC) progression and prognosis.

Methods: Pan-cancer analysis was conducted to analyze expression, genetic and epigenetic alterations, and immunological characteristics of EVA1B. Especially, immunological characteristics and mutational landscape were compared between high and low EVA1B expression groups in the combined TCGA-COAD and TCGA-READ datasets. Through random survival forest analysis, an EVA1B-derived genomic model was developed, and its prognostic value was verified in the external datasets (GSE14333, GSE39582, and GSE87211). Drug sensitivity was compared between high- and low-risk subpopulations. A nomogram was conducted through integrating independent factors.

Results: EVA1B expression presented a remarkable upregulation in most cancer types, especially CRC. EVA1B expression was significantly correlated to DNA methyltransferases, DNA mismatch repair genes, m⁶A regulators, TMB, and MSI across pan-cancer. High EVA1B expression indicated an undesirable CRC patients’ prognosis. Additionally, its upregulation was correlated to enhanced immune cell infiltration, increased stromal and immune activation, and elevated activities of cancer immunity cycle. Higher frequencies of amplification and deletion were investigated in high EVA1B expression subpopulation. Following verification, the EVA1B-derived genomic model reliably predicted patients’ prognosis and drug responses. The nomogram (age, stage, EVA1B-derived risk score) was conducted to quantify an individual’s survival probability. Furthermore, our experimental validation based on immunohistochemistry indicated that EVA1B overexpression is correlated with CRC tumorigenesis and poor outcomes in our CRC patients’ cohort.

Conclusion: Collectively, our findings provided valuable resource for guiding the mechanisms and therapeutic analysis of EVA1B in CRC.

Introduction

Colorectal cancer (CRC) involving the colon and rectum ranks the third leading cause of cancer mortality globally, with over 1.85 million cases and 850,000 deaths each year (1). Recently, CRC presents increasing incidence among younger cases (2). This disease is a heterogeneous disease with distinct pathogenesis mechanisms, involving somatic mutation, genetic fusion, genetic instability, as well as epigenetic alteration (3). Among newly diagnosed CRC, 20% of cases have occurring metastasis at diagnosis as well as 25% will metastasize following localized disease (1). Surgical resection is the major therapeutic regimen of CRC. Nevertheless, there are very few treatment regimens for metastatic patients. Although chemotherapy is usually recommended, merely few targeted therapies appropriate for cases who have specific mutation profiling are available, such as EGFR inhibitor and VEGF inhibitor (3). Hence, it is urgently required to develop novel molecule targets against CRC.

Immunotherapy is a novel alternative against CRC treatment, which utilizes cases’ own immune system to combat tumor cells (4, 5). A minority of CRC cases present microsatellite instability (MSI), a molecule predictor of defective DNA mismatch repair, and a predictive biomarker of immunotherapeutic response, but most of them are microsatellite-stable (MSS) (6). Enhanced tumor mutational burden (TMB) and neoantigen load in MSI cancers sustain the infiltrations of immune effector subpopulations, and anti-cancer immune response is strongly in relation to MSS counterparts (6). Nevertheless, a few MSS cancers present increased TMB as well as infiltrating immune cells to respond to immunotherapy. Hence, novel predictors of immunotherapeutic responses will be required. EVA1B is a protein coding gene, which is a critical paralog of EVA1A gene. It is an endoplasmic reticulum and lysosome-relevant protein involving autophagy and apoptosis (7). Previously, it triggers papillary thyroid carcinogenesis and epithelial-mesenchymal transition (EMT) through Hippo signaling (8). Flubendazole exerts an antitumor role through modulating EVA1A-mediated autophagy and apoptosis in breast carcinoma (9). MiR-125b relieves oxaliplatin resistance in liver carcinoma through negatively modulating EVA1A-mediated autophagy (10). EVA1B possesses high sequence similar to EVA1A gene, and EVA1B protein presents the similar domain to EVA1A protein (11). Nevertheless, the role of EVA1B in CRC remains indistinct. This study was conducted to investigate the correlation of EVA1B expression with prognosis, tumor immune and pharmacogenomic features in CRC.

Materials and Methods

Pan-Cancer Analyses

RNA sequencing data (FPKM), somatic mutational data [Mutation Annotation Format (MAF)], copy number variation (CNV) data, and clinical information of 33 cancer types [adrenocortical carcinoma (ACC); bladder urothelial carcinoma (BLCA); breast invasive carcinoma (BRCA); cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC); cholangiocarcinoma (CHOL); colon adenocarcinoma (COAD); lymphoid neoplasm diffuse large B-cell lymphoma (DLBC); esophageal carcinoma (ESCA); glioblastoma multiforme (GBM); head and neck squamous cell carcinoma (HNSC); kidney chromophobe (KICH); kidney renal clear cell carcinoma (KIRC); kidney renal papillary cell carcinoma (KIRP); acute myeloid leukemia (LAML); lower-grade glioma (LGG); liver hepatocellular carcinoma (LIHC); lung adenocarcinoma (LUAD); lung squamous cell carcinoma (LUSC); mesothelioma (MESO); ovarian serous cystadenocarcinoma (OV); pancreatic adenocarcinoma (PAAD); pheochromocytoma and paraganglioma (PCPG); prostate adenocarcinoma (PRAD); rectum adenocarcinoma (READ); sarcoma (SARC); skin cutaneous melanoma (SKCM); stomach adenocarcinoma (STAD); testicular germ cell tumors (TGCT); thyroid carcinoma (THCA); thymoma (THYM); uterine corpus endometrial carcinoma (UCEC); uterine carcinosarcoma (UCS); uveal melanoma (UVM)] were curated from The Cancer Genome Atlas (TCGA) through the GDC data portal (https://portal.gdc.cancer.gov/). Tumor Immune Estimation Resource (TIMER; version 2.0) web tool (http://timer.cistrome.org/) was adopted to analyze the expression of EVA1B across 33 cancer types (12). For TCGA dataset, FPKM value was transformed to TPM value. The expression profiling of four major DNA methyltransferases (DNMT1, DNMT2, DNMT3A, and DNMT3B), DNA mismatch repair genes (MLH1, MSH2, MSH6, PMS2, and EPCAM), and m⁶A regulators (CBLL1, VIRMA, METTL3/14, RBM15/15B, WTAP, ZC3H13, ALKBH5, FTO, ELAVL1, FMR1, HNRNPA2B1, HNRNPC, IGF2BP1/2/3, LRPPRC, YTHDC1/2, and YTHDF1/2/3) was extracted from pan-cancer specimens. TMB and MSI of pan-cancer were also harvested from TCGA project. Spearman correlation analysis was conducted to investigate the correlation of EVA1B with above factors across pan-cancer.

Estimation of Immunological Features

The abundance of tumor-infiltrating immune cells was inferred with single-sample gene set enrichment analysis (ssGSEA) derived from gene set variation analysis (GSVA) package (13). In total, 122 immunomodulatory factors containing MHCs, receptors, chemokines, and immune stimulators were collected from Charoentong et al. (14). Additionally, immune checkpoint molecules were retrieved from the study of Auslander et al. (15). Overall infiltration of stromal and immune cells was estimated in CRC tissues with estimation of stromal and immune cells in malignant tumors using expression data (ESTIMATE) algorithm in accordance with mRNA expression data (16). Cancer immunity cycle was curated from previous research (17), and the activities of all steps were estimated with ssGSEA (18).

Collection of CRC Datasets

The GSE14333 (19), GSE39582 (20), and GSE87211 (21) datasets were collected from the Gene Expression Omnibus (GEO; https://www.ncbi.nlm.nih.gov/geo/) project. Batch effects were adjusted with ComBat algorithm derived from sva package (22). Somatic mutational data were visualized with maftools package (23). GISTIC2.0 was adopted to analyze amplification and deletion using CNV data (24). Raw “CEL” files of microarray profiling from Affymetrix platform were downloaded, followed by background correction and quantile normalization via robust multiarray averaging algorithm utilizing affy and simpleaffy packages (25, 26). Meanwhile, normalized matrix files of microarray profiling from other platforms were directly curated.

Curation of Gene Sets of Known Biological Processes

Gene sets of known biological processes containing CD8+ T effector, DNA damage repair, pan-fibroblast TGF-β response signature (pan-F-TBRS), antigen processing machinery, immune checkpoints, EMT1-3, FGFR3-relevant gene signatures, angiogenesis, KEGG discovered histones, Fanconi anemia, cell cycle, DNA replication, nucleotide excision repair, homologous recombination, mismatch repair, WNT target, as well as cell cycle regulators were curated from previous research (27–29). The activities of biological processes were quantified with ssGSEA algorithm.

Functional Enrichment Analysis

GSEA was conducted to investigate the differences in signaling pathways activated in two subpopulations, with the gene set “c2.cp.kegg.v6.2.symbols.gmt” as the reference. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of EVA1B-relevant genes was conducted with clusterProfiler package (30). GO categories contained biological process (BP), cellular component (CC), and molecular function (MF). Hallmark gene sets were harvested from the Molecular Signatures Database (31), and their activities were quantified with ssGSEA algorithm.

Screening EVA1B-Relevant Genes

To identify EVA1B-relevant genes, CRC patients were classified into high and low EVA1B expression subpopulations with the mean value of EVA1B expression. Through empirical Bayesian method from limma package (32), differentially expressed genes between subpopulations were determined. The significance criteria of EVA1B-relevant genes were set as |fold-change| >1.5 and adjusted p-value <0.01.

Exploitation of an EVA1B-Derived Genomic Model

Through univariate Cox regression models, interaction of EVA1B-relevant genes with CRC prognosis was estimated in TCGA cohort. EVA1B-relevant genes with p-value <0.05 were determined as prognostic factors, which were input into random survival forest analysis (33). The number of Monte Carlo iterations was set as 100, and the number of steps forward was 5. Thereafter, their relative importance was ranked. In accordance with relative importance >0.28, characteristic EVA1B-relevant genes were determined, which were input into a multivariate Cox regression model. The formula of EVA1B-relevant genomic model was conducted as follows: risk score = ${\Sigma }_{k-1}^n{\text{ Exp}}_i\ast e^{\text{HR}}i$, in which n indicated the number of characteristic EVA1B-relevant genes, Exp_i indicated the expression of characteristic genes, and e^HRi indicated the regression coefficients of genes derived from the multivariate Cox regression analysis. In accordance with this formula, risk score of each CRC patient was calculated. Thereafter, CRC patients were clarified into high- and low-risk subpopulations following the mean value of risk score. The risk score distribution, survival state, and heatmap of the expression of characteristic EVA1B-relevant genes were drawn following patients’ risk score. Kaplan–Meier curves of overall survival (OS), disease-specific survival (DSS), and progression-free survival (PFS) were depicted between high- and low-risk subpopulations using survival and survminer packages, followed by log-rank test. Thereafter, receiver operating characteristic (ROC) curves of 1-, 3-, and 5-year survival were established utilizing survivalROC package, and area under the curve (AUC) was calculated to assess the discrimination (34). Additionally, prognostic value of the EVA1B-derived genomic model was externally verified in the GSE14333, GSE39582, and GSE87211 cohorts.

Assessment of Drug Sensitivity

The Genomics of Drug Sensitivity in Cancer (GDSC) project (https://www.cancerrxgene.org/) offers drug sensitivity information of 138 anticancer agents across approximately 75,000 experiments in 700 cancer cell lines (35). The half-maximal inhibitory concentration (IC50) that represented drug response was estimated with pRRophetic package (36).

Nomogram Establishment

Uni- and multivariate Cox regression models were established to evaluate the associations of clinical indicators (age, gender, stage, T, N, and M) and EVA1B-relevant risk score with CRC patients’ OS. ROC curves were conducted to evaluate the prognostic value of this nomogram. Additionally, calibration curves were utilized to assess the consistency between the actual and nomogram-predicted survival probabilities.

Clinical Specimens and Data Source

Formalin-fixed, paraffin-embedded specimens, including primary carcinoma specimens (n = 131), corresponding noncancerous normal tissues (n = 19) and liver metastasis specimens (n = 19) for immunohistochemistry (IHC) were obtained from 131 CRC patients who underwent surgery between 2014 and 2015. CRC primary carcinomas were assessed according to the 8th edition of the American Joint Committee on Cancer (AJCC) staging system. CRC patient data and tissue samples were obtained from the Liaoning Cancer Hospital. Patients included in the current study did not receive preoperative radiotherapy or chemotherapy prior to the study. The current study (20210804GP) was approved by the Ethics Committee of the Liaoning Cancer Hospital, Shenyang, China.

Immunohistochemistry

EVA1B expression levels in samples were determined by IHC. Tissue sections were incubated with EVA1B antibody (Manufacturer: Atlas Antibodies; Catalog number: HPA043537). The level of EVA1B expression was determined by counting the percentages of positively stained immunoreactive cells and evaluating cell staining intensity. EVA1B IHC staining was scored as “−” ~ “+” and “++” ~ “+++,” which represented negative and positive, respectively. All samples were reviewed by two independent, experienced pathologists in who were blinded to the identity of the samples.

Statistical Analysis

Statistical analysis was implemented with R software (v3.4.1; https://www.r-project.org/) and its appropriate packages. Comparison between groups was conducted utilizing Student’s t- or Wilcoxon rank-sum test. Spearman correlation test was adopted to determine the interactions between variables. Chi-square test was performed to analyze the correlation between EVA1B and clinicopathologic characteristics. Kaplan–Meier curves and the log-rank test were used to compare the OS between EVA1B expression level. The prognostic ability of the predictors for OS was evaluated by ROC curves and the AUC values. Univariate and multivariate Cox regression analyses were utilized to evaluate the independent prognostic value of EVA1B regarding OS. p-value <0.05 was set as the threshold.

Results

Analysis of Expression, Genetic and Epigenetic Alterations, and Immunological Characteristics of EVA1B Across Pan-Cancer

Figure 1 depicts the workflow of this study. EVA1B, as known as FAM176B, displayed remarkably increased expression in most cancer types in comparison with corresponding normal tissues, containing BLCA, BRCA, CHOL, COAD, ESCA, GBM, HNSC, KIRC, LIHC, PRAD, READ, STAD, and THCA (Figure 2A). Differently, downregulated EVA1B was investigated in CESC, KICH, KIRP, and UCEC. DNA methyltransferases exert critical roles in changing chromatin structure and gene expression. We noted that EVA1B presented negative interactions with four major DNA methyltransferases containing DNMT1, DNMT2, DNMT3A, and DNMT3B in most cancer types (Figure 2B). Additionally, EVA1B was negatively correlated to DNA mismatch repair genes (MLH1, MSH2, MSH6, PMS2, and EPCAM) across cancer types (Figure 2C). M⁶A methylation represents the most common RNA modification, which affects the complexity of cancer progression (37). Especially, we focused on the interactions of EVA1B with m⁶A regulators in CRC. With the mean value of EVA1B expression, CRC patients were clustered into high and low EVA1B expression subpopulations. As depicted in Figure 2D, downregulated EVA1B was correlated to most m⁶A regulators like CBLL1, VIRMA, METTL14, and RBM15. TMB becomes a reliable biomarker of sensitivity to immune checkpoint blockage (38). In Figure 2E, EVA1B presented remarkably positive interaction with TBM across THYM, and LGG. Oppositely, it was negatively correlated with TMB across BRCA, CESC, ESCA, HNSC, STAD, and UCEC. MSI is a hypermutation phenotype caused by frequent polymorphisms in short repetitive DNA sequence as well as single nucleotide substitution due to DNA MMR defects (39). We noted the positive interaction of EVA1B with MSI in THCA, BLCA, DLBC, HNSC, KIRC, and PRAD but negative interaction with MSI in COAD, READ, STAD, and UCEC (Figure 2F). Thereafter, we analyzed the immunological characteristics of EVA1B across pan-cancer. Our data suggested that EVA1B was remarkably correlated to immune cell infiltration (Figure 2G), MCH molecules (Figure 2H), chemokines (Figure 2I), immunostimulatory factors (Figure 2J), receptors (Figure 2K), and immune checkpoint molecules (Figure 2L) across pan-cancer.

FIGURE 1

Figure 1 The workflow of this study.

FIGURE 2

Figure 2 Analysis of expression, genetic and epigenetic alterations, and immunological characteristics of EVA1B across pan-cancer. (A) TIMER analysis identifies the difference in expression of EVA1B between diverse cancer types and matched normal specimens. (B) Circle diagram visualizes the interaction of EVA1B with four major DNA methyltransferases: DNMT1, DNMT2, DNMT3A, and DNMT3B across pan-cancer. The first outer ring, cancer types; the second ring, four DNA methyltransferases; the third ring, correlation coefficients; the fourth ring, p-values; and numbers in the inner ring, correlation coefficients and p-values. (C) Heatmap visualizes the relationship of EVA1B with five DNA mismatch repair genes: MLH1, MSH2, MSH6, PMS2, and EPCAM across diverse cancer types. For each relationship, the top left triangle is colored to indicate correlation coefficients, while the bottom right triangle is colored to indicate p-values. (D) Box plots depict the difference in expression of m⁶A regulators in high and low EVA1B expression groups across CRC specimens. (E) Radar chart shows the association of EVA1B with TMB in each cancer type. (F) Radar chart shows the association of EVA1B with MSI across pan-cancer. (G–L) Heatmaps visualize the association of EVA1B with (G) tumor-infiltrating immune cells, (H) MHC molecules, (I) chemokines, (J) immunostimulatory factors, (K) receptors, and (L) immune checkpoint molecules across cancer types. For each relationship, the top left triangle is colored to indicate p-values, while the bottom right triangle is colored to indicate correlation coefficients. ACC, adrenocortical carcinoma; BLCA, bladder urothelial carcinoma; BRCA, breast invasive carcinoma; CESC, cervical squamous cell carcinoma and endocervical adenocarcinoma; CHOL, cholangiocarcinoma; COAD, colon adenocarcinoma; DLBC, lymphoid neoplasm diffuse large B-cell lymphoma; ESCA, esophageal carcinoma; GBM, glioblastoma multiforme; HNSC, head and neck squamous cell carcinoma; KICH, kidney chromophobe; KIRC, kidney renal clear cell carcinoma; KIRP, kidney renal papillary cell carcinoma; LAML, acute myeloid leukemia; LGG, lower-grade glioma; LIHC, liver hepatocellular carcinoma; LUAD, lung adenocarcinoma; LUSC, lung squamous cell carcinoma; MESO, mesothelioma; OV, ovarian serous cystadenocarcinoma; PAAD, pancreatic adenocarcinoma; PCPG, pheochromocytoma and paraganglioma; PRAD, prostate adenocarcinoma; READ, rectum adenocarcinoma; SARC, sarcoma; SKCM, skin cutaneous melanoma; STAD, stomach adenocarcinoma; TGCT, testicular germ cell tumors; THCA, thyroid carcinoma; THYM, thymoma; UCEC, uterine corpus endometrial carcinoma; UCS, uterine carcinosarcoma; UVM, uveal melanoma. (^*p < 0.05; ^**p < 0.01; ^***p < 0.001).

Prognostic Significance of EVA1B and Its Association With Drug Sensitivity in CRC

We also focused on the biological role and significance of EVA1B in CRC. We combined the COAD and READ datasets from TCGA project and removed batch effects (Figures 3A, B). Survival analysis was indicative that patients with high EVA1B expression possessed more undesirable OS outcomes than those with low EVA1B expression (Figure 3C). Targeted therapies play a key role in CRC management, and genetic alterations that attribute to CRC heterogeneity is correlated to the responses to targeted therapies (40). Reliable predictive biomarkers for targeted therapies remain greatly lacking. Thus, we evaluated whether EVA1B expression was correlated to drug sensitivity. As a result, EVA1B presented positive association with sensitivity to bleomycin, lenvatinib, zoledronate, floxuridine, mitoxantrone, simvastatin, topotecan, and dasatinib (Figure 3D). Differently, EVA1B expression was negatively correlated to sensitivity to nilotinib, tamoxifen, panobinostat, docetaxel, vinorelbine, crizotinib, ethinyl estradiol, and paclitaxel in CRC. Therefore, EVA1B might be a predictive marker of above agents.

FIGURE 3

Figure 3 Prognostic significance of EVA1B and its association with drug sensitivity in CRC. (A, B) PCA diagrams depict the combined COAD and READ datasets from TCGA project before and after removal of batch effects. (C) Kaplan–Meier curves of OS for high and low EVA1B expression subpopulations. (D) Scatter plots show the correlation of EVA1B expression with drug sensitivity in CRC.

Immunological and Biological Significance of EVA1B in CRC

Our ssGSEA results demonstrated the remarkably increased infiltration levels of most immune cells in high EVA1B expression subpopulation (Figure 4A). Additionally, high EVA1B expression was correlated to increased stromal and immune score as well as reduced tumor purity in CRC (Figure 4B). There were prominently enhanced activities of most steps within cancer immunity cycle in high EVA1B expression subpopulation (Figure 4C). We also investigated the increased activities of immune activation processes (like CD8+ T effector, and immune checkpoint) and stromal activation processes (like EMT, pan-F-TBRS, FGFR3-related genes, and angiogenesis) in high EVA1B expression subpopulation (Figure 4D). Additionally, the interactions of EVA1B expression with cancer immunity cycle and known biological processes were evaluated across CRC specimens. In Figure 4E, EVA1B expression presented positive correlations to most steps within cancer immunity cycle as well as immune and stromal activation processes.

FIGURE 4

Figure 4 Immunological and biological significance of EVA1B in CRC. (A) Comparison of the infiltration of immune cells between high and low EVA1B expression subpopulations. (B) Evaluation of the difference in stromal score, immune score, and tumor purity in high and low EVA1B expression subpopulations. (C) Quantification of the activities of cancer immunity cycle in two subpopulations. (D) Comparing the activation of known biological signatures in two subpopulations. (E) Associations of EVA1B expression with the activities of cancer immunity cycle and known biological signatures across CRC specimens. (^*p < 0.05; ^**p < 0.01; ^***p < 0.001).

Association of EVA1B With Mutational Landscape in CRC

In Figure 5A, we evaluated the prevalence of somatic mutation in high and low EVA1B expression subpopulations. APC, TP53, and TTN ranked the first three mutational genes. Nevertheless, no prominent difference in somatic mutation was investigated in high and low EVA1B expression subpopulations. The GISTIC2.0 results demonstrated that amplification and deletion displayed higher frequencies in high EVA1B expression subpopulation (Figures 5B, C) compared with low EVA1B expression subpopulation (Figures 5D, E). Additionally, we calculated the G-score in accordance with the amplitude of the aberrations and the frequency of their incidence in CRC specimens. Following comparison, mutations occurred in more regions in high EVA1B expression subpopulation (Figures 5F, G).

FIGURE 5

Figure 5 Association of EVA1B with mutational landscape and signaling pathways in CRC individuals. (A) Landscape of somatic mutation in high and low EVA1B expression CRC subpopulations. Genes are ranked following the mutational frequency. Upper portion of the chart displays TMB score of each patient. (B, C) Landscape of (B) amplification and (C) deletion in high EVA1B expression subpopulation. (D, E) Landscape of (D) amplification and (E) deletion in low EVA1B expression subpopulation. The genome is oriented vertically from top to bottom; GISTIC2.0 q-value at each locus is depicted from left to right. The green line indicates the cutoff value of q-value = 0.25. (F, G) Detection and comparison of amplification and deletion of copy number in high and low EVA1B expression subpopulations. (H) GSEA for the signaling pathways activated in high EVA1B expression subpopulation. (I) GSEA for the signaling pathways activated in low EVA1B expression subpopulation.

Signaling Pathways Associated With EVA1B

We further analyzed the signaling pathways involving EVA1B via GSEA. In Figure 5H, high EVA1B expression was positively correlated to cell adhesion molecules (CAMs), chemokine signaling pathways, cytokine–cytokine receptor interaction, Leishmania interaction, leukocyte trans-endothelial migration, and VEGF signaling pathways. Meanwhile, low EVA1B expression was negatively correlated to cell cycle, homologous recombination, mismatch repair, nucleotide excision repair, oocyte meiosis, and spliceosome (Figure 5I). The above data showed that EVA1B was involved in several oncogenic pathways, indicating the oncogenic role of EVA1B.

Identification of EVA1B-Derived Genes and Their Biological Significance

Through the cutoffs of |fold-change| >1.5 and adjusted p-value <0.01, we determined 602 EVA1B-derived genes in CRC individuals (Supplementary Table S1). Their biological significance was further investigated. In Figure 6A, EVA1B-derived genes were mainly correlated to immune response in accordance with GO annotation results. Additionally, there were remarkable interactions of EVA1B-derived genes with tumorigenic pathways (such as PI3K-Akt signaling pathway, proteoglycans in cancer, and ECM-receptor interaction) and immune activation pathways (such as IL-17 signaling pathway, Th17-cell differentiation, cytokine–cytokine receptor interaction, Th1- and Th2-cell differentiation, antigen processing and presentation, intestinal immune network for IgA production, complement and coagulation cascades; Figure 6B). Overall, EVA1B-derived genes might exert remarkable roles in CRC progression.

FIGURE 6

Figure 6 Exploitation of an EVA1B-derived genomic model for CRC prognosis. (A, B) GO and KEGG enrichment analyses of EVA1B-derived genes. (C) Random survival forest analysis for calculating the relative importance of EVA1B-derived genes. (D) Determination of characteristic EVA1B-derived genes with relative importance >0.28. (E) Distribution of high- and low-risk CRC subpopulations with the mean value (vertical dashed line) of EVA1B-derived risk score. (F) Kaplan–Meier curves of OS for high and low CRC subpopulations. (G) Distribution of survival status in two CRC subpopulations. (H) Heatmap for the expression of IGFBP3, EGFL7, SULT1B1, MMP1, CYP2W1, as well as ISYNA1 in two CRC subpopulations. (I) ROC curves at 1-, 3-, and 5-year OS outcomes in accordance with EVA1B-derived risk score.

The EVA1B-Derived Genomic Model for CRC Prognosis

In total, 66 EVA1B-derived genes were prominently correlated to CRC patients’ prognosis utilizing univariate-Cox regression models (Table 1). With random survival forest analysis, we ranked the relative importance of above prognostic EVA1B-derived genes (Figure 6C). In accordance with relative importance >0.28, six characteristic EVA1B-derived genes were determined, containing IGFBP3, EGFL7, SULT1B1, MMP1, CYP2W1, as well as ISYNA1 (Figure 6D). On the basis of them, a multivariate-Cox regression model was conducted and risk scoring system was developed following the formula: risk score = (−0.117800092) * ISYNA1 expression + 0.067856101 * CYP2W1 expression + (−0.145078003) * MMP1 expression + (−0.107319254) * SULT1B1 expression + 0.427836868 * EGFL7 expression + 0.188926758 * IGFBP3 expression. With the mean value, we classified CRC patients into high- and low-risk subpopulations (Figure 6E). Survival analysis uncovered that low-risk CRC individuals displayed a remarkable survival advantage (Figure 6F). Nevertheless, no distinct difference in survival status was noted between high- and low-risk subpopulations (Figure 6G). Additionally, the heterogeneity in expression of IGFBP3, EGFL7, SULT1B1, MMP1, CYP2W1, as well as ISYNA1 was displayed in two subpopulations (Figure 6H). ROC curves confirmed that the EVA1B-derived genomic model possessed the well potency in estimating 1-, 3-, and 5-year OS probabilities (Figure 6I).

TABLE 1

Table 1 Univariate-Cox regression models determine prognostic EVA1B-derived genes in CRC.

Verification of Prognostic Significance of EVA1B-Derived Genomic Model in CRC

Except for OS, our analysis demonstrated that high-risk subpopulation was indicative of more undesirable DSS and PFS outcomes in comparison with low-risk subpopulation in TCGA cohort (Figures 7A, B). The prognostic significance of the EVA1B-derived genomic model was verified in external cohorts. Our data confirmed that high-risk score was predictive of unfavorable OS outcome in the GSE14333, GSE39582, and GSE87211 datasets (Figures 7C–E).

FIGURE 7

Figure 7 Verification of prognostic significance of EVA1B-derived genomic model in CRC. (A, B) Kaplan–Meier curves of DSS and PFS for high and low CRC subpopulations in TCGA cohort. (C–E) Kaplan-Meier curves of OS for high and low CRC subpopulations in the (C) GSE14333, (D) GSE39582, and (E) GSE87211 datasets.

The EVA1B-Derived Genomic Model Predicts Drug Responses of CRC Patients

Further analysis was conducted to uncover whether this EVA1B-derived genomic model was predictive of drug responses of CRC individuals. In Figure 8A, we noted that high-risk group had significantly lower estimated IC50 of sunitinib and docetaxel relative to low-risk group, indicating that high-risk subpopulations presented higher sensitivity to sunitinib and docetaxel. Meanwhile, the low-risk group showed significantly lower estimated IC50 of sorafenib and gemcitabine, indicating that low-risk subpopulations were more likely to respond to sorafenib and gemcitabine. We also evaluated the interactions of EVA1B-derived genes with drug responses. As a result, MMP1 was negatively correlated to IC50 of mithramycin, actinomycin D, depsipeptide, and homoharringtonine; EGFL7 presented positive correlation to IC50 of raltitrexed, gemcitabine, cytarabine, and TFdU; IGFBP3 was positively correlated to IC50 of lenvatinib and JNJ-42756493 but negatively correlated to IC50 of LEE-011, cobimetinib (isomer 1), and selumetinib; CYP2W1 displayed positive relationship to IC50 of tegafur and fluorouracil (Figure 8B). These data indicated that EVA1B-derived genes might be correlated to sensitivity to above agents.

FIGURE 8

Figure 8 EVA1B-derived genomic model predicts drug responses of CRC patients. (A) Comparison of the estimated IC50 values of sunitinib, sorafenib, gemcitabine, and docetaxel between high- and low-risk CRC subpopulations. (B) Associations of EVA1B-relevant genes with IC50 values of small molecule compounds in CRC.

Establishment of a Reliable Nomogram for Prediction of CRC Prognosis

Univariate-Cox regression analysis uncovered the remarkable association of age, stage, T, N, M, and EVA1B-derived risk score with CRC prognosis (Figure 9A). Among them, age, T, and EVA1B-derived risk score served as independent prognostic indicators of CRC (Figure 9B). By integrating above three independent prognostic indicators, this study exploited a nomogram to estimate CRC patients’ survival outcomes (Figure 9C). This EVA1B-derived risk score occupied the most contribution to the prediction of 1-, 3-, and 5-year OS duration. ROC curves confirmed the desirable efficacy in predicting patients’ survival outcomes (Figure 9D). Additionally, we evaluated the prediction performance of the nomogram with calibration curves. Our data demonstrated that 1-, 3-, and 5-year predicted by this nomogram was close to the actual survival duration (Figures 9E–G). Above data indicated the superior predictive capacity of this nomogram.

FIGURE 9

Figure 9 Establishment of a reliable nomogram for prediction of CRC prognosis. (A, B) Uni- and multivariate Cox regression models are conducted to uncover the association of clinical features and EVA1B-derived risk score with CRC survival outcome. (C) A prognostic nomogram is exploited through integrating independent prognostic indicators (age, T stage, and EVA1B-derived risk score) to estimate 1-, 3-, and 5-year survival probability. (D) Predictive efficacy of this nomogram is verified through ROC curves at 1-, 3-, and 5-year survival. (E–G) Calibration plots show the association of predicted 1-, 3-, and 5-year OS with actual survival duration.

Molecular Mechanisms Involving the EVA1B-Derived Genomic Signature

We further investigated the molecular mechanisms involving the EVA1B-derived genomic signature. As depicted in Figure 10A, unfolded protein response, E2F targets, G2M checkpoint, mTORC1 signaling, MYC targets, DNA repair, and oxidative phosphorylation presented higher activities in low-risk subpopulation. Meanwhile, tumorigenic pathways were remarkably activated in high-risk subpopulation, such as P53 pathway, TGF-beta signaling, apoptosis, Hedgehog signaling, EMT, Notch signaling, and angiogenesis, indicative of high-risk subpopulation’s unfavorable survival outcomes. Moreover, immune activation pathways such as allograft rejection, inflammatory response, and complement were significantly activated in high-risk subpopulation.

FIGURE 10

Figure 10 Molecular mechanisms involving EVA1B-derived genomic signature as well as prognostic value of EVA1B-derived genes and their interaction with m⁶A regulators and immune cell infiltration in CRC. (A) Heatmap visualizes the activities of hallmark gene signatures in high- and low-risk subpopulations. (B) Forest diagram depicts the interaction of EVA1B-derived genes with CRC clinical outcomes in accordance with univariate-Cox regression models. (C) Heatmap depicts the association of EVA1B-derived genes with m⁶A regulators across CRC specimens. (D) Heatmap visualizes the interaction of EVA1B-derived genes with immune cell infiltration in CRC. (E) Association among drug sensitivity, predictive immune responses, and risk status of CRC individuals. (^*p < 0.05; ^**p < 0.01).

Prognostic Significance of EVA1B-Derived Genes and Their Interaction With m⁶A Methylation and Immune Cell Infiltration

Figure 10B demonstrates that TRADD, MXD4, NR1H2, EVA1B, and C11orf68 acted as risk factors of CRC outcomes. Additionally, we investigated that above EVA1B-derived genes presented positive correlations to m⁶A regulator ALKBH5 and RBM15B but negative correlations to other m⁶A regulators in CRC (Figure 10C). In Figure 10D, these EVA1B-derived genes were negatively correlated to activated CD4+ T cell, effector memory CD4+ T cell, memory B cell, and type 2 T helper cell but positively correlated to other tumor-infiltrating immune cells in CRC. Figure 10E shows the association among drug sensitivity, predictive immune responses, and risk status of CRC individuals.

EVA1B Overexpression Is Correlated With CRC Tumorigenesis and Poor Outcomes in CRC Patients

EVA1B IHC was performed on 19 pairs of primary CRC, corresponding noncancerous and matching liver metastasis specimens from 19 liver metastasis CRC patients. The findings showed that EVA1B was overexpressed in primary CRC tissues compared with the corresponding noncancerous normal tissues (p < 0.01). In addition, expression level of EVA1B was similar in the liver metastases samples to the expression level in matching primary CRC tissues without statistical significance (Figures 11A, B).

FIGURE 11

Figure 11 EVA1B overexpression is correlated with CRC tumorigenesis and poor outcomes in CRC patients. (A, B) Comparison of EVA1B expression in noncancerous mucosa, paired primary CRC, and liver metastases by IHC. (C) Kaplan–Meier plots of CRC specimens with negative and positive EVA1B expression. (D) Time-dependent ROC analysis of EVA1B showing the overall survival of patients with CRC. (E, F) Univariate and multivariate Cox regression of EVA1B expression and clinical parameters in CRC. ***p < 0.001; ns, no significance.

In addition, the relationship between EVA1B expression and clinicopathologic characteristics of 131 patients with CRC was explored (Table 2). Expression of EVA1B was highly correlated with pN stage (p < 0.001), pM stage (p = 0.005), AJCC stage (p = 0.001), CEA (p = 0.003), and CA199 (p = 0.002). However, the finding showed no significant correlation between EVA1B expression level and age, gender, tumor location, differentiation, and pT stage of CRC patients (Table 2).

TABLE 2

Table 2 EVA1B staining expression associations with clinicopathologic and survival characteristics in CRC patients.

Kaplan–Meier curve was generated to explore the prognosis value of EVA1B in these 131 CRC patients. The findings showed that CRC patients with high expression levels of EVA1B presented with poor OS compared with those with low expression levels of EVA1B (p < 0.001, Figure 11C). AUC values for EVA1B expression (0.630, 0.679, and 0.720 for the 2-, 4-, and 6-year OS, respectively) showed that EVA1B accurately discriminated between the high- and low-expression groups (Figure 11D). Univariate and multivariate Cox analyses were performed, indicating that EVA1B was an independent as a prognostic factor by adjusting for stage, gender, location, differentiation, stage, CEA, and CA199 (Figures 11E, F).

Discussion

This study presented an integrative analysis of molecular characteristics, oncogenic roles, and relevant immune and pharmacogenomic features of EVA1B across pan-cancer, especially CRC. Our data demonstrated the remarkable upregulation of EVA1B expression in most cancer types. Additionally, its upregulation was significantly correlated to DNA methyltransferases, DNA mismatch repair genes, m⁶A regulators, TMB, and MSI across pan-cancer. Tumorigenesis is a multistep process where normal cells acquire genetic and epigenetic alterations contributing to tumor initiation and progression (41, 42). The pathological mechanisms of CRC initiation and progression include chromosomal instability, high MSI and methylation, contributing to oncogene, tumor suppressor gene, and mismatch repair-relevant gene mutations (43). Our data were indicative of the critical role of EVA1B in rectal and colon carcinogenesis.

High EVA1B expression was indicative of undesirable CRC patients’ clinical outcomes. Nevertheless, we did not note the distinct difference in survival status between high- and low-risk subpopulations due to this retrospective cohort. The role of EVA1B in survival status of CRC should be investigated in a prospective cohort in our future studies. Previous research proposed the upregulation of EVA1B in glioma (11). At present, the therapeutic regimen is surgical resection plus chemotherapy in more advanced or inoperable patients (44). Nevertheless, immunotherapy to supplement curative and palliative therapy has been undergoing clinical trials (45). Our further analysis demonstrated that EVA1B upregulation was correlated to enhanced immune cell infiltration, increased stromal and immune activation, and elevated activities of cancer immunity cycle. Additionally, higher frequencies of amplification and deletion were noted in high EVA1B expression subpopulation. Our data were indicative that EVA1B might be utilized as a promising immunotherapeutic predictor.

We further determined 602 EVA1B-derived genes in CRC that were mainly correlated to immune response, tumorigenic pathways (such as PI3K-Akt signaling pathway, proteoglycans in cancer, and ECM-receptor interaction), and immune activation pathways (such as IL-17 signaling pathway, Th17 cell differentiation, cytokine-cytokine receptor interaction, Th1 and Th2 cell differentiation, antigen processing and presentation, intestinal immune network for IgA production, complement and coagulation cascades), indicative of their remarkable roles in CRC progression. Utilizing random survival forest analysis, we conducted an EVA1B-derived genomic model, composed of IGFBP3, EGFL7, SULT1B1, MMP1, CYP2W1, and ISYNA1. Following verification in diverse external datasets, the genomic model reliably and independently predicted patients’ prognosis and relapse. Additionally, it possessed the potential in estimating drug responses, like sunitinib, sorafenib, gemcitabine, and docetaxel. There was remarkable correlation of characteristic EVA1B-derived genes with responses to small molecular compounds. TRADD, MXD4, NR1H2, EVA1B, and C11orf68 served as remarkable risk factors of CRC outcomes. Previously, modulation of TRADD restores cellular homeostasis as well as alleviates apoptosis (46). It triggers tumor inhibition through regulation of ULF-dependent p19Arf ubiquitylation (47). CDK4/6 inhibitor upregulates MXD4 expression that negatively modulates MYC in CD8+ T cells (48). NR1H2 modulates cholesterol homeostasis within human cells, controlling fitness and function of activated T cells (49). Nevertheless, more experimental evidence should be conducted to uncover the functions of above EVA1B-relevant genes in CRC.

A nomogram acts as a powerful tool to quantify an individual’s risk in a clinical setting through integrating diverse risk factors (50). Herein, this study utilized the nomogram in prediction of 1-, 3-, and 5-year OS probabilities via incorporating age, T stage, as well as EVA1B-derived genomic model for CRC individuals. Each indicator was assigned a score on the basis of its contribution to survival risk. Thereafter, ROC curves demonstrated that the nomogram presented the favorable efficiency in predicting an individual patient’ OS outcome. Additionally, calibrated curves confirmed the actual survival duration agreed with the nomogram-estimated survival duration. A few limitations of our research will be pointed out. The functional role of EVA1B in tumor immunity requires in-depth experimental verification. Moreover, prognostic value of EVA1B needs to be verified in larger CRC cohorts.

Conclusion

Overall, our study conducted an integrative analysis to uncover molecular characteristics, oncogenic roles, and relevant immune and pharmacogenomic features of EVA1B in CRC. Our findings indicated that EVA1B acted as a convincing prognostic marker as well as a predictor of therapeutic responses.

Data Availability Statement

Publicly available datasets were analyzed in this study. These data can be found here: https://portal.gdc.cancer.gov/ and GEO, https://www.ncbi.nlm.nih.gov/geo/, under accession numbers GSE14333, GSE39582, and GSE87211.

Ethics Statement

The studies involving human participants were reviewed and approved by the human ethics committees of the Liaoning Cancer Hospital. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

BM and YML conceived and designed this study. YML and KW performed the bioinformatics analyses and visualization. YL collected data and performed the statistical analysis. BM, YML and QM wrote the original draft. BM and YML revised the manuscripts. All authors revised and approved the final manuscript.

Funding

This study was supported by the National Natural Science Foundation of China (81902383); the Revitalizing Liaoning Talents Program (XLYC1907004); the Young and Middle-aged Scientific and Technological Innovation Talent Support Plan of Shenyang City (RC200223); Cultivation Program of National Science Foundation of Liaoning Cancer Hospital (2021-ZLLH-03); Cultivation Program of National Science Foundation of Liaoning Cancer Hospital (2020-ZLLH-44); Natural Guiding Plan Foundation of Liaoning Province (2019-ZD-0584); and Beijing Xisike Clinical Oncology Research Foundation (Y-QL202101-0039).

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/fimmu.2022.809837/full#supplementary-material

Supplementary Table S1 | The information of 602 EVA1B-derived genes in CRC patients.

Abbreviations

CRC, colorectal cancer; IHC, immunohistochemistry; MSI, microsatellite instability; MSS, microsatellite-stable; TMB, tumor mutational burden; CNV, copy number variation; TCGA, The Cancer Genome Atlas; TIMER, Tumor Immune Estimation Resource; ssGSEA, single-sample gene set enrichment analysis; GSVA, gene set variation analysis; ESTIMATE, estimation of stromal and immune cells in malignant tumors using expression data; EMT, epithelial-mesenchymal transition; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; OS, overall survival; DSS, disease-specific survival; PFS, progression-free survival; ROC, receiver operating characteristic; AUC, area under the curve; GDSC, Genomics of Drug Sensitivity in Cancer; IC50, half-maximal inhibitory concentration.

References

1. Biller LH, Schrag D. Diagnosis and Treatment of Metastatic Colorectal Cancer: A Review. Jama (2021) 325(7):669–85. doi: 10.1001/jama.2021.0106

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Nguyen M, Tipping Smith S, Lam M, Liow E, Davies A, Prenen H, et al. An Update on the Use of Immunotherapy in Patients With Colorectal Cancer. Expert Rev Gastroenterol Hepatol (2021) 15(3):291–304. doi: 10.1080/17474124.2021.1845141

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Xu H, Liu L, Li W, Zou D, Yu J, Wang L, et al. Transcription Factors in Colorectal Cancer: Molecular Mechanism and Therapeutic Implications. Oncogene (2021) 40(9):1555–69. doi: 10.1038/s41388-020-01587-3

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Zhao J, Zhong S, Niu X, Jiang J, Zhang R, Li Q. The MHC Class I-LILRB1 Signalling Axis as a Promising Target in Cancer Therapy. Scand J Immunol (2019) 90(5):e12804. doi: 10.1111/sji.12804

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Kong X, Fu M, Niu X, Jiang H. Comprehensive Analysis of the Expression, Relationship to Immune Infiltration and Prognosis of TIM-1 in Cancer. Front Oncol (2020) 10:1086. doi: 10.3389/fonc.2020.01086

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Picard E, Verschoor CP, Ma GW, Pawelec G. Relationships Between Immune Landscapes, Genetic Subtypes and Responses to Immunotherapy in Colorectal Cancer. Front Immunol (2020) 11:369. doi: 10.3389/fimmu.2020.00369

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Zhao S, Wang H. EVA1A Plays an Important Role by Regulating Autophagy in Physiological and Pathological Processes. Int J Mol Sci (2021) 22(12):6181. doi: 10.3390/ijms22126181

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Lin BY, Wen JL, Zheng C, Lin LZ, Chen CZ, Qu JM. Eva-1 Homolog A Promotes Papillary Thyroid Cancer Progression and Epithelial-Mesenchymal Transition via the Hippo Signalling Pathway. J Cell Mol Med (2020) 24(22):13070–80. doi: 10.1111/jcmm.15909

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Zhen Y, Zhao R, Wang M, Jiang X, Gao F, Fu L, et al. Flubendazole Elicits Anti-Cancer Effects via Targeting EVA1A-Modulated Autophagy and Apoptosis in Triple-Negative Breast Cancer. Theranostics (2020) 10(18):8080–97. doi: 10.7150/thno.43473

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Ren WW, Li DD, Chen X, Li XL, He YP, Guo LH, et al. MicroRNA-125b Reverses Oxaliplatin Resistance in Hepatocellular Carcinoma by Negatively Regulating EVA1A Mediated Autophagy. Cell Death Dis (2018) 9(5):547. doi: 10.1038/s41419-018-0592-z

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Qu S, Liu J, Wang H. EVA1B to Evaluate the Tumor Immune Microenvironment and Clinical Prognosis in Glioma. Front Immunol (2021) 12:648416. doi: 10.3389/fimmu.2021.648416

PubMed Abstract | CrossRef Full Text | Google Scholar

12. 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(W1):W509–w514. doi: 10.1093/nar/gkaa407

PubMed Abstract | CrossRef Full Text | Google Scholar

13. 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

CrossRef Full Text | Google Scholar

14. Sokolov A, Carlin DE, Paull EO, Baertsch R, Stuart JM. Pathway-Based Genomics Prediction Using Generalized Elastic Net. PloS Comput Biol (2016) 12(3):e1004790. doi: 10.1371/journal.pcbi.1004790

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Auslander N, Zhang G, Lee JS, Frederick DT, Miao B, Moll T, et al. Robust Prediction of Response to Immune Checkpoint Blockade Therapy in Metastatic Melanoma. Nat Med (2018) 24(10):1545–9. doi: 10.1038/s41591-018-0157-9

PubMed Abstract | CrossRef Full Text | Google Scholar

16. 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

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Chen DS, Mellman I. Oncology Meets Immunology: The Cancer-Immunity Cycle. Immunity (2013) 39(1):1–10. doi: 10.1016/j.immuni.2013.07.012

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Xu L, Deng C, Pang B, Zhang X, Liu W, Liao G, et al. TIP: A Web Server for Resolving Tumor Immunophenotype Profiling. Cancer Res (2018) 78(23):6575–80. doi: 10.1158/0008-5472.Can-18-0689

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Jorissen RN, Gibbs P, Christie M, Prakash S, Lipton L, Desai J, et al. Metastasis-Associated Gene Expression Changes Predict Poor Outcomes in Patients With Dukes Stage B and C Colorectal Cancer. Clin Cancer Res (2009) 15(24):7642–51. doi: 10.1158/1078-0432.Ccr-09-1431

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Marisa L, de Reyniès A, Duval A, Selves J, Gaub MP, Vescovo L, et al. Gene Expression Classification of Colon Cancer Into Molecular Subtypes: Characterization, Validation, and Prognostic Value. PloS Med (2013) 10(5):e1001453. doi: 10.1371/journal.pmed.1001453

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Hu Y, Gaedcke J, Emons G, Beissbarth T, Grade M, Jo P, et al. Colorectal Cancer Susceptibility Loci as Predictive Markers of Rectal Cancer Prognosis After Surgery. Genes Chromosomes Cancer (2018) 57(3):140–9. doi: 10.1002/gcc.22512

PubMed Abstract | CrossRef Full Text | Google Scholar

22. 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(6):882–3. doi: 10.1093/bioinformatics/bts034

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Mayakonda A, Lin DC, Assenov Y, Plass C, Koeffler HP. Maftools: Efficient and Comprehensive Analysis of Somatic Variants in Cancer. Genome Res (2018) 28(11):1747–56. doi: 10.1101/gr.239244.118

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Mermel CH, Schumacher SE, Hill B, Meyerson ML, Beroukhim R, Getz G. GISTIC2.0 Facilitates Sensitive and Confident Localization of the Targets of Focal Somatic Copy-Number Alteration in Human Cancers. Genome Biol (2011) 12(4):R41. doi: 10.1186/gb-2011-12-4-r41

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Gautier L, Cope L, Bolstad BM, Irizarry RA. Affy–Analysis of Affymetrix GeneChip Data at the Probe Level. Bioinformatics (2004) 20(3):307–15. doi: 10.1093/bioinformatics/btg405

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Wilson CL, Miller CJ. Simpleaffy: A BioConductor Package for Affymetrix Quality Control and Data Analysis. Bioinformatics (2005) 21(18):3683–5. doi: 10.1093/bioinformatics/bti605

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Rosenberg JE, Hoffman-Censits J, Powles T, van der Heijden MS, Balar AV, Necchi A, et al. Atezolizumab in Patients With Locally Advanced and Metastatic Urothelial Carcinoma Who Have Progressed Following Treatment With Platinum-Based Chemotherapy: A Single-Arm, Multicentre, Phase 2 Trial. Lancet (2016) 387(10031):1909–20. doi: 10.1016/s0140-6736(16)00561-4

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Şenbabaoğlu Y, Gejman RS, Winer AG, Liu M, Van Allen EM, de Velasco G, et al. Tumor Immune Microenvironment Characterization in Clear Cell Renal Cell Carcinoma Identifies Prognostic and Immunotherapeutically Relevant Messenger RNA Signatures. Genome Biol (2016) 17(1):231. doi: 10.1186/s13059-016-1092-z

PubMed Abstract | CrossRef Full Text | Google Scholar

29. 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(7693):544–8. doi: 10.1038/nature25501

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Yu G, Wang LG, Han Y, He QY. Clusterprofiler: An R Package for Comparing Biological Themes Among Gene Clusters. Omics (2012) 16(5):284–7. doi: 10.1089/omi.2011.0118

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) Hallmark Gene Set Collection. Cell Syst (2015) 1(6):417–25. doi: 10.1016/j.cels.2015.12.004

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. Limma Powers Differential Expression Analyses for RNA-Sequencing and Microarray Studies. Nucleic Acids Res (2015) 43(7):e47. doi: 10.1093/nar/gkv007

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Wang H, Zhou L. Random Survival Forest With Space Extensions for Censored Data. Artif Intell Med (2017) 79:52–61. doi: 10.1016/j.artmed.2017.06.005

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Heagerty PJ, Lumley T, Pepe MS. Time-Dependent ROC Curves for Censored Survival Data and a Diagnostic Marker. Biometrics (2000) 56(2):337–44. doi: 10.1111/j.0006-341x.2000.00337.x

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Yang W, Soares J, Greninger P, Edelman EJ, Lightfoot H, Forbes S, et al. Genomics of Drug Sensitivity in Cancer (GDSC): A Resource for Therapeutic Biomarker Discovery in Cancer Cells. Nucleic Acids Res (2013) 41(Database issue):D955–961. doi: 10.1093/nar/gks1111

PubMed Abstract | CrossRef Full Text | Google Scholar

36. 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(9):e107468. doi: 10.1371/journal.pone.0107468

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Ma S, Chen C, Ji X, Liu J, Zhou Q, Wang G, et al. The Interplay Between M6a RNA Methylation and Noncoding RNA in Cancer. J Hematol Oncol (2019) 12(1):121. doi: 10.1186/s13045-019-0805-7

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Chalmers ZR, Connelly CF, Fabrizio D, Gay L, Ali SM, Ennis R, et al. Analysis of 100,000 Human Cancer Genomes Reveals the Landscape of Tumor Mutational Burden. Genome Med (2017) 9(1):34. doi: 10.1186/s13073-017-0424-2

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Baretti M, Le DT. DNA Mismatch Repair in Cancer. Pharmacol Ther (2018) 189:45–62. doi: 10.1016/j.pharmthera.2018.04.004

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Messersmith WA. NCCN Guidelines Updates: Management of Metastatic Colorectal Cancer. J Natl Compr Canc Netw (2019) 17(5.5):599–601. doi: 10.6004/jnccn.2019.5014

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Chen L, Wang G, Qiao X, Wang X, Liu J, Niu X, et al. Downregulated miR-524-5p Participates in the Tumor Microenvironment of Ameloblastoma by Targeting the Interleukin-33 (IL-33)/Suppression of Tumorigenicity 2 (ST2) Axis. Med Sci Monit (2020) 26:e921863. doi: 10.12659/msm.921863

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Niu X, Xu J, Liu J, Chen L, Qiao X, Zhong M. Landscape of N(6)-Methyladenosine Modification Patterns in Human Ameloblastoma. Front Oncol (2020) 10:556497. doi: 10.3389/fonc.2020.556497

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Li J, Liang L, Yang Y, Li X, Ma Y. N(6)-Methyladenosine as a Biological and Clinical Determinant in Colorectal Cancer: Progression and Future Direction. Theranostics (2021) 11(6):2581–93. doi: 10.7150/thno.52366

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Alexander PG, McMillan DC, Park JH. A Meta-Analysis of CD274 (PD-L1) Assessment and Prognosis in Colorectal Cancer and its Role in Predicting Response to Anti-PD-1 Therapy. Crit Rev Oncol Hematol (2021) 157:103147. doi: 10.1016/j.critrevonc.2020.103147

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Qiu XT, Song YC, Liu J, Wang ZM, Niu X, He J. Identification of an Immune-Related Gene-Based Signature to Predict Prognosis of Patients With Gastric Cancer. World J Gastrointest Oncol (2020) 12(8):857–76. doi: 10.4251/wjgo.v12.i8.857

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Xu D, Zhao H, Jin M, Zhu H, Shan B, Geng J, et al. Modulating TRADD to Restore Cellular Homeostasis and Inhibit Apoptosis. Nature (2020) 587(7832):133–8. doi: 10.1038/s41586-020-2757-z

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Chio II, Sasaki M, Ghazarian D, Moreno J, Done S, Ueda T, et al. TRADD Contributes to Tumour Suppression by Regulating ULF-Dependent p19Arf Ubiquitylation. Nat Cell Biol (2012) 14(6):625–33. doi: 10.1038/ncb2496

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Heckler M, Ali LR, Clancy-Thompson E, Qiang L, Ventre KS, Lenehan P, et al. Inhibition of CDK4/6 Promotes CD8 T-Cell Memory Formation. Cancer Discovery (2021) 11(10):2564–81. doi: 10.1158/2159-8290.Cd-20-1540

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Michaels AJ, Campbell C, Bou-Puerto R, Rudensky AY. Nuclear Receptor Lxrβ Controls Fitness and Functionality of Activated T Cells. J Exp Med (2021) 218(4):e20201311. doi: 10.1084/jem.20201311

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Miao Y, Li Q, Wang J, Quan W, Li C, Yang Y, et al. Prognostic Implications of Metabolism-Associated Gene Signatures in Colorectal Cancer. PeerJ (2020) 8:e9847. doi: 10.7717/peerj.9847

PubMed Abstract | CrossRef Full Text | Google Scholar