Development and Validation of a 6-Gene Hypoxia-Related Prognostic Signature For Cholangiocarcinoma

Cholangiocarcinoma (CHOL) is highly malignant and has a poor prognosis. This study is committed to creating a new prognostic model based on hypoxia related genes. Here, we established a novel tumor hypoxia-related prognostic model consisting of 6 hypoxia-related genes by univariate Cox regression and the least absolute shrinkage and selection operator (LASSO) algorithm to predict CHOL prognosis and then the risk score for each patient was calculated. The results showed that the patients with high-risk scores had poor prognosis compared with those with low-risk scores, which was verified as an independent predictor by multivariate analysis. The hypoxia-related prognostic model was validated in both TCGA and GEO cohorts and exhibited excellent performance in predicting overall survival in CHOL. The PPI results suggested that hypoxia-related genes involved in the model may play a central role in regulating the hypoxic state. In addition, the presence of IDH1 mutations in the high-risk group was high, and GSEA results showed that some metabolic pathways were upregulated, but immune response processes were generally downregulated. These factors may be potential reasons for the high-risk group with worse prognosis. The analysis of different immune regulation-related processes in the high- and low-risk groups revealed that the expression of genes related to immune checkpoints would show differences between these two groups. We further verified the expression of the oncogene PPFIA4 in the model, and found that compared with normal samples, CHOL patients were generally highly expressed, and the patients with high-expression of PPFIA4 had a poor prognosis. In summary, the present study may provide a valid prognostic model for bile duct cancer to inform better clinical management of patients.


INTRODUCTION
Cholangiocarcinoma (CHOL) is a malignant and aggressive disease with a poor prognosis, the median survival less than 24 months (1). Depending on the location of the tumor, CHOL is divided into intrahepatic cholangiocarcinoma (iCCA) and extra-hepatic cholangiocarcinoma (eCCA), and eCCA is divided into perihilar cholangiocarcinoma (pCCA) and distal cholangiocarcinoma (dCCA) (2). The early symptoms of CHOL are insidious and not easily detected, and most patients are diagnosed at an advanced stage (2). For patients with advanced CHOL, treatment options are limited. The median overall survival (OS) of patients treated with standard chemotherapy regimens (gemcitabine and cisplatin) is still less than 1 year (3). Therefore, it is clinically important to find the relevant factors affecting the prognosis of patients with CHOL. The American Joint Committee on Cancer (AJCC) staging manual has become a benchmark for classifying cancer patients, determining prognosis, and determining the best treatment (4). However, it has been discovered that the TNM staging approach is insufficient for determining prognosis and does not account for cancer's biological heterogeneity. There is significant variation in prognosis and treatment response even among patients with the same stage, and other factors such as age, performance status, and tumor site might affect patient survival, therefore it only provides limited information on clinical prognosis (5). As a result, developing reliable prognostic biomarkers is critical in order to improve clinical prognostic value. With the continuous improvement and development of high-throughput sequencing technology (6), it is possible to exploit potential molecular targets in the CHOL genome using bioinformatics technology. The development of promising molecular targets is essential to improve the application of targeted therapy for CHOL.
The tumor microenvironment is a complex and variable cellular microenvironment, which can not only provide the material conditions required for the growth and proliferation of tumor cells, but also change its composition through autocrine and paracrine secretion of tumor cells (7). It was reported that CHOL cells can construct their own favorable environment by secreting tumor-related regulatory mediators through the extracellular matrix and stromal cells, thus promoting the proliferation of CHOL cells and enhancing their resistance to treatment (8). Most solid tumors include hypoxia as a micro environmental trait, and developing tumors frequently live in hypoxic environments due to limited blood supply (9). The hypoxic conditions could enhance the angiogenesis, proliferation and invasive abilities of tumor cells (10,11). Similarly, in CHOL studies, hypoxia has been found to increase the aggressiveness of CHOL (12). In iCCA studies, hypoxic environment was able to promote the progression of CHOL through upregulation of Rab1a (13). It has been found that in CHOL, hypoxia-induced Sonic Hedgehog signaling pathway regulates cancer stemness, epithelial-to-mesenchymal transition and invasion (14).
In the present study, the hypoxia-related genes for CHOL were mainly screened by GEO and TCGA-CHOL dataset. The prognostic signature characterized by hypoxia-related genes were filtered according to univariate and LASSO regression analysis, and the risk score of each samples were obtained. Subsequently, CHOL samples in TCGA were classified into high-and low-risk subgroups according to the risk score, and the differences between high-and low-risk subgroups were further explored, such as survival differences, clinical characteristics differences, enrichment differences, somatic cell mutation differences, and tumor immunomodulatory gene expression differences.

Pre-Processing and Data Collection
We utilized the R package TCGAbiolinks (15) to download the 36 CHOL RNA-seq data [log2(FPKM+1)] and relevant clinical information from the TCGA database (https://portal.gdc.cancer. gov/). The 36 samples with both available RNA-seq data and clinical information were used as the training set for hypoxiarelated markers. The GSE107943 dataset was obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/), which were used as the external validation set, containing available microarray expression data and clinical data of bile duct cancer samples. The microarray data preprocessing process was as follows: remove null probes and probes corresponding to multiple genes at the same time, and use the median expression value of these probes as the expression value of the gene when multiple probes correspond to the same gene.

Prognostic Model Construction
The 200 potential hypoxia-related genes were retrieved from MSigDB's HALLMARK HYPOXIA pathway (https://www.gseamsigdb.org/gsea/index.jsp) (v7.4 version). All prognosisassociated hypoxia genes were screened using univariate Cox hazard analysis, with a p-value of 0.1 chosen as the cutoff. Utilizing the R package glmnet, the Lasso Cox regression model was built using all prognosis-related hypoxia genes, and candidate genes were penalized to exclude redundant components and reduce the risk of overfitting (16). A 10-fold cross-validation under the constraint of low bias is used to calculate the model's penalty parameter (l).
To determine the best prognosis-associated hypoxia signature, Lasso regression analysis was used, and the regression coefficients were calculated for each gene separately, and then the sample risk scores for multivariate Cox regression were calculated using the R package survminer (17), with the formula as follows: bi * ci bi denotes the weight of each gene in Lasso Cox regression; ci denotes the expression level of each gene.
The CHOL patients were divided into high-risk and low-risk groups based on Risk score. Kaplan-Meier curves of survival differences between samples in the high-risk and low-risk groups were plotted, and significant differences in OS between groups were assessed based on log-rank tests.

Correlation Analysis of Prognostic Models and Clinical Characteristics
The univariate and multivariate Cox hazard regression analysis were performed separately for risk scores and patient clinical characteristics (including age, gender, and stage) to determine whether the predictive power of the prognostic model was independent of clinical characteristics. Risk ratios (HR) and 95% confidence intervals were calculated for each candidate prognostic factor, with a p-value < 0.05 as the threshold of significance. In addition, differences in risk scores across age, sex, and stage subgroups were assessed based on the Wilcoxon ranksum test.

Gene Interactions Analysis
The expression correlation between genes in the prognostic model was assessed based on Pearson correlation analysis. In addition, the interactions of proteins encoding hypoxia-related genes were obtained based on the STRING database (http:// www.string-db.org/), and protein-protein interactions (PPI) network of hypoxia-related genes were constructed, and the network graphs were visualized using Cytoscape (v3.9.1) to assess the interactions of the genes in the prognostic model.

Functional Enrichment Analysis
Gene set enrichment analysis (GSEA) was used to identify biological processes that differed significantly between highand low-risk group samples (18). We examined biological processes from the GO database and KEGG signaling pathways for significant variations between high-and low-risk group samples using the R tool ClusterProfiler (19). To screen the collection of substantially linked genes, a corrected p-value of less than 0.05 was employed as a criterion, and the p-value correction technique Benjamini Hochberg was applied.

Somatic Cell Mutations
TCGA provided the tumor mutational data. The mutation genes were calculated using tumor mutational burden (TMB). The mutational data in both the high-risk and low-risk groups were analyzed using the R maftools package (20). P-values less than 0.05 were deemed statistically significant.

Expression of Immune-Negative Regulatory Genes
The cancer immune cycle has emerged as the primary research framework for cancer immunotherapy. It describes a cyclical process in which the immune system eradicates cancer: Antigen Processing and Presentation (Category 1), Atimicrobials (Category 2), Chemokines (Category 3), Cytokine Recepors (Category 4), Cytokines (Category 5), Matrix remodeling (Category 6), and Nature Killing Cell Cytoxicity (Category 7) (21). In the low and high hypoxia risk groups, we looked at the expression of genes that negatively regulate these processes. Tracking Tumor Immunophenotype (http://biocc.hrbmu.edu. cn/TIP/index.jsp) was used to identify immunonegative regulation-related genes (22).

Tissue Specimens
From September 2021 to March 2022, 10 paired fresh cholangiocarcinoma tissues and neighboring non-tumor tissues were taken from patients who underwent radical surgery at Xi'an Jiaotong University's First Affiliated Hospital. In addition, between December 2013 and April 2016, 62 patients with cholangiocarcinoma who had a curative procedure were enrolled in a cohort. The protocol was authorized by the First Affiliated Hospital of Xi'an Jiaotong University's Institutional Medical Ethics Committee (No. XJTU1AF2021LSK-261).

RNA Isolation and RT-qPCR
The RNeasy Mini Kit (Qiagen) was used to extract total RNA, and the reverse transcription was done with the High Capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific). A LightCycler480 PCR equipment was used to perform quantitative PCR using SYBR Green dye (Roche). GAPDH was used to standardize the relative fold change in expression.

Statistical Analysis
All statistical analyses were performed using R software (version 4.1.2; https://www.R-project.org). Wilcoxon rank-sum test was used to compare significant differences in continuous variables between the two groups. Pearson correlation analysis was used to analyze the expression correlation between genes. Differentially expressed genes between high and low risk groups were identified based on the R package limma (23) with a differentially expressed gene threshold of |log2FC| > 0.58 and p-values < 0.05. The univariate Cox regression analysis was performed using the R package survival, and multi-factor Cox regression analysis was performed on prognostic model scores to explore the the independent prognostic value of the score. Kaplan-Meier curves were plotted using the R package survminer, and the significance of differences between groups was assessed based on the log-rank test. And p-values < 0.05 were used as the threshold of significance for statistical analysis.

Establishment of a Hypoxia-Related Prognostic Model for CHOL
We acquired 200 hypoxia-related genes from the MSigDB database to study the predictive efficiency of hypoxia-related genes in CHOL. EGFR, PFKL, CHST2, ADORA2B, EDN2, HMOX1, CP, PPFIA4, PLAC8, IDS, ALDOC, STBD1, BCL2 were among 13 hypoxia-related genes found to be substantially linked with overall survival (OS) in the TCGA CHOL cohort using univariate Cox hazard analysis (p-value < 0.1, Figure 1A). Based on the median value of gene expression (EGFR, PFKL, CHST2, ADORA2B, EDN2, HMOX1), we separated the samples into two groups, and Kaplan-Meier curve analysis revealed that gene expression affected patient survival (log-rank test p-value < 0.05, Figure 1B).
Although we initially observed genes with prognostic efficacy in CHOL patients by univariate Cox regression analysis and logrank test, we further performed Lasso regression analysis based on these 13 genes in order to integrate these informative genes to obtain a more optimal prognostic model, while removing redundant factors and controlling the risk of overfitting. Tenfold cross-validation was performed under optimal conditions to determine the penalty parameter (l) of the model. We screened out 6 most predictive factors for OS ( Figures 1C, D). PPFIA4 and CP were retained as a valid prognostic risk factor, while BCL2, PLAC8, PFKL and EGFR were retained as a prognostic protective factor, Together, they constituted a prognostic risk model associated with hypoxia in CHOL. Based on the linear combination of the expression levels of these genes and the corresponding weights, we could assess the hypoxia-related prognostic risk score for each patient ( Figure 1E). The corresponding association between the expression of 6 hypoxia-related signature in the prognostic model and the vital status and survival of patients were visualized in Figure 1F.

Independent Validation of the Hypoxia-Related Prognostic Model
The risk score for each patient with CHOL in the TCGA was determined using the hypoxia-related prognostic model for CHOL that was constructed. After obtaining each patient's risk value, the samples were sorted from low to high in terms of risk score, and scatter plots were used to display the results. The median of all sample risk scores was used as the cutoff, and the samples were divided into two groups: high risk and low risk, the low-risk group's OS was likewise considerably lower than the high-risk group's (log-rank test p-value < 0.001, Figure 2A). We also calculated heat maps for display to investigate the expression of different genes in the risk model across all genes in different samples, and we discovered that the expression of six hypoxiarelated genes in the model differed considerably in the high-and low-risk groups ( Figure 2B).
To assess the robustness and generalizability of the hypoxiarelated prognostic model, we collected CHOL samples from GSE107943 as an independent validation cohort. Patients were divided into high-risk and low-risk groups according to the risk score of the prognostic model. We obtained similar results to the training cohort through an external validation cohort ( Figures 2C, D).

Hypoxia-Related Risk Score Is an Independent Prognostic Factor
We then explored the independence of the prognostic model in patients with CHOL from TCGA. The univariate Cox regression analysis revealed that the hypoxia-related risk score was the only prognostic factor of OS (HR = 5.94, p-value < 0.001, Figure 3A). Further, the multivariate Cox regression analysis confirmed the independent prognostic value of the hypoxia-related risk score (HR = 74.71, p-value = 0.004, Figure 3B).
In the GEO cohort, univariate Cox analysis showed that hypoxia-related risk score (HR = 2.65, p-value = 0.046), stage (HR = 3.48, p-value = 0.015) were prognostic factors in patients with CHOL ( Figure 3C). Hypoxia-related risk score (HR = 3.06, p-value = 0.032) and stage (HR = 3.88, p-value = 0.017) were found to be independent predictive markers in patients with CHOL in a multifactorial Cox analysis ( Figure 3D), which is consistent with the TCGA cohort's findings. The prognostic risk scores established in our study might be used as independent risk factors affecting patient prognosis based on the above investigations. To further corroborate the independence of prognostic risk scores from clinical characteristics, we compared the differences in the distribution of risk scores by grouping patients based on their age, gender and stage, respectively. We found no significant differences in risk scores by age and tumor stage grouping (pvalue > 0.05, Figures 3E, G), but slightly higher in male patients than in female patients (p-value < 0.05, Figure 3F).

The Interaction Mechanism Between Prognostic Genes
Although the 6 signature genes we obtained from the prognostic model were all involved in the tumor hypoxic state hallmark, but their interactions with other hypoxic genes were not clear, so we examined the co-expression of these genes. Overall, there was no strong co-expression among the hypoxia-related genes involved in the prognostic model (Pearson |r| < 0.34, Figure 4A), indicating that the redundant genes with strong associations had been filtered out through LASSO Cox analysis. In addition, we looked at how the hypoxia-related genes in the model interacted with other hypoxia-related genes. We built a PPI network connecting the genes in the prognostic model and other hypoxia-related genes using the interactions of gene-encoded proteins retrieved from the String database online tool ( Figure 4B).
We found extensive interactions between EGFR, CP, PFKL genes and other hypoxia-associated genes in the model. The model includes hypoxia-related genes, which might play a key role in tumor hypoxia regulation.

Hypoxia-Related Risk Xcores Reveal Differences in Molecular Tumor Characteristics
To explore the underlying molecular mechanisms of hypoxiarelated prognostic models, we first examined somatic mutations in different risk groups in the TCGA-CHOL cohort. Due to the small sample size of the TCGA-CHOL dataset, the comparison of differences in mutation frequencies failed to present significant results ( Figure 5A, p-value < 0.05). However, we found four patients in the high-risk group enriched for the IDH1 mutation, while only one patient in the low-risk group was enriched for this mutation ( Figure 5A).

Hypoxia-Related Risk Scores Reveal Differences in Tumor Immune Microenvironment Characteristics
Tumor-immune interactions are crucial in controlling the immune milieu because they incorporate the effects of several regulators. We looked examined how immune regulatory genes were expressed differently in high-risk and low-risk people. We discovered that genes involved in immune modulation were expressed differently in high and low risk groups ( Figure 6A). We looked examined how divergent genes were expressed in distinct immune regulatory groups. The genes HLA-DPB1, HLA-DQA1, HLA-DRA, HLA-DRB1 and HLA-DRB5 in the "Antigen Processing and Presentation" category were considerably downregulated in the high-risk group. PLA2G2A, KNG1, HRG, and other genes associated to "antimicrobials" were considerably downregulated in the high-risk group. CCL17, a gene associated to "chemokines," was considerably downregulated in the high-risk group. The gene EGFR, which belongs to the "Cytokine Receptors" category, was considerably downregulated in the high-risk group. The genes CTF1, GDNF, HAMP, IGF2, SAA1, SAA2, TGFA, which belong to the "cytokines" category, were considerably downregulated in the high-risk group. CA9 and PLOD2, which are connected to "matrix remodeling," were considerably downregulated in the high-risk group. The gene PLCG2, which is connected to "Natural Killer Cell Cytotoxicity," was considerably downregulated in the high-risk group. The above results are shown in figures (Figures 6A, B). These findings show that immune checkpoint expression differs under the effect of CHOL's hypoxic condition, affecting the immunological microenvironment's control.

The Expression of PPFIA4 Is Dramatically Increased in CHOL Tissues
To further corroborate these findings, PPFIA4 mRNA levels were measured in 10 matched fresh CHOL tissues with surrounding non-tumor tissues. When compared to nont u m o r t i s s u e s , t h e a m o u n t o f P P F I A 4 m R N A i n cholangiocarcinoma tissues became much higher ( Figure 7A). In addition, we investigated PPFIA4 expression in a cohort of patient with cholangiocarcinoma samples from our hospital by IHC assay, which demonstrated predominant expression of PPFIA4 within tumor cells ( Figure 7C, D). The expression of PPFIA4 was found in 26 of the 62 carcinoma samples (41.9%), whereas the remaining 36 samples (58.1%) were negative. Table 1 outlined the link between PPFIA4 expression and clinicopathological characteristics. The presence of PPFIA4 was significantly linked to large tumor size, advanced TNM stage, and lymph node metastasis, while there was no correlation of PPFIA4 expression with gender, age, differentiation and tumor stage. Additionally, individuals with high PPFIA4 expression had a poorer 5-year OS prognosis ( Figure 7B, p-value = 0.03).

DISCUSSION
The molecular mechanisms underlying the development of CHOL are still unclear. It was reported that hypoxia enhances the aggressiveness of CHOL cells (12) and the malignant behavior of tumors can be regulated through hypoxia-related gene (13,24). Gene signatures have been utilized to predict the prognosis of many malignancies in recent years, and they have proven to be more accurate than TNM staging and histological diagnosis (25,26). Hypoxia has been identified as a factor impacting patient prognosis that is independent of clinically established prognostic variables such as tumor stage, lymph node status, and tumor grade (27). Although approaches such as nitroimidazole, PET imaging, and IHC biomarker expression have been used to determine the degree of hypoxia in patients' tumors, the precise proven method remains unknown (28,29). For a number of malignancies, such as prostate cancer, lung cancer, and breast cancer, prognostic models based on hypoxiarelated gene expression have been described (30)(31)(32). In this study, we used a Cox and Lasso regression model to analyze 200 hypoxia-related genes, and then screened 6 prognosis-related hypoxia-related genes and constructed a prognostic prediction model for CHOL, and the OS of the high-risk group was significantly lower than that of the low-risk group, and we validated our model using the GEO CHOL date set. The prognostic model we constructed is an ideal model with better sensitivity and specificity.
Phosphofructokinase L (PFKL) includes Binding locations for hypoxia-inducible factor 1 (HIF-1) (33). HIF-1 is a nuclear transcription factor that is active in the oxygen state where the body is adapted to hypoxia or hypoxia. HIF-1a has a certain transcriptional activity and a broad target gene spectrum (34). In hypoxia, the transcriptional and translational levels of HIF-1a can increase exponentially with prolonged hypoxia (35). Dfructose 6-phosphate is converted to D-fructose 1,6bisphosphonate, a crucial step in glucose metabolism, by an enzyme encoded by PFKL (glycolysis) (36). In the tumor hypoxic region, tumor cells translocate extracellular glucose into the cell via glucose transporter proteins on the cell membrane to produce pyruvates and ATP via the glycolytic pathway, and LDH-A acts to convert it to lactate (37). The production of lactate constitutes a local acidic environment for tumor cell proliferation and activates the cascade reaction of protein hydrolases to convert pre-matrix metalloproteinases into matrix metalloproteinases, which promote the degradation of the extracellular matrix of cancer cells and facilitate the invasion of surrounding tissues by tumor cells (38). A number of cancers have been identified to have PFKL as a component of the glycolysis process, including lung cancer (39), esophageal cancer (40), neuroblastoma (41) and hepatocellular carcinoma (42). The Bcl-2 protein family is one of the core regulatory mechanisms of apoptosis, which can receive and transmit intrinsic intracellular signals or external environmental stress signals, such as nutrient or hypoxic stress, DNA damage, oncogene over-activation, endoplasmic reticulum stress, etc (43). It was found that Bcl-2-associated transcription factor BCLAF1 was highly expressed under hypoxia with a timedependent hypoxic effect and positively correlated with HIF1a expression (44). PTPRF Interacting Protein Alpha 4 (PPFIA4) was reported to be glycolysis-associated signature, involved in glycolysis of tumor microenvironment (45). Proliferation and migration of tumor cells in colon cancer are facilitated by PPFIA4's ability to increase tumor glycolysis (46); PPFIA4 also increases mitochondrial metabolism via MTHFD2 in desmoplastic-resistant prostate cancer (47). In our study, we found that PPFIA4 was oncogenic in CHOL and its high expression was associated with poor prognosis of patients. We confirmed the high expression of PPFIA4 in the tumor tissues of CHOL patients through fresh tumor tissues and their corresponding paraneoplastic tissue and previous pathological slides collected clinically. Likewise, we confirmed that high PPFIA4 expression was significantly associated with poor patient prognosis. We hypothesized that the expression status of PPFIA4 in CHOL tumor tissues could reflect the hypoxic status in tumor tissues. It has been discovered that PLAC8 enhances human interstitial extravillous trophoblast cell invasion and migration, and hypoxia dramatically upregulates PLAC8 expression in this cell population (48). In hypoxic conditions, Ceruloplasmin (CP) induces generates oxygen radicals (49), while inducing activation of HIF-1a (50), Overexpression of CP in clear cell renal carcinoma is linked to oncogenic pathways and worse survival rates (51). Epidermal Growth Factor Receptor (EGFR) is a primary regulator of HIF-1 under hypoxic conditions, and mutant EGFR is a prominent regulator of HIF-1a in EGFR mutant NSCLC cells (52). In addition, EGFR was reported that it can blocked anti-cancer miRNA maturation under hypoxic status (53). Somatic mutations between high-and low-risk groups were compared in this study. Among the high-risk individuals, we observed that Isocitrate Dehydrogenase (NADP(+)) 1 (IDH1) mutations were more common than expected. IDH1 has been linked to tumor hypoxic reprogramming, according to a study (54). Mutations in IDH1 reduce IDH1's catalytic activity and activate HIF-1 a in glioblastoma cells (55). Wild-type IDH1 suppresses tumor development in renal cell carcinoma via degrading HIF-a (56). By using GSEA to examine the variations in cellular activities and pathways between high and low risk groups, we found that metabolic processes and pathways associated with hypoxic state were upregulated, but immune response processes were generally downregulated. Prognostic models developed by our team may be able to respond well to tumor microenvironment hypoxia.
Genes associated to the negative regulation of the tumor immune, such as checkpoints and suppressor-related genes, were examined between high and low risk groups to study the hypoxia/immunity association. Immune checkpoints are targets of immunotherapy, and researchers have found that patients with tumors co-expressing HIF-1a, a subtype of HIF, and PD-L1, an immune checkpoint, have a high risk of tumor recurrence and metastasis as well as lethality (57). In CHOL, hypoxic condition altered the regulatory impact on the immune milieu by modifying the expression of immunological checkpoints. In the high-risk scenario, the hypoxic condition could boost the ability of tumors to escape and invade.
It is possible that our research contained some flaws. In the first place, the research is a retrospective one with a relatively limited sample size. Therefore, the sparse mutation landscape was not suitable for further investigation of genetic factors associated with the hypoxia-related risk score. For example, only three samples showed KRAS and TP53 mutations. In addition, the insufficient clinical information of patients and lack of diversity also hindered the further in-depth analysis. In the present study, most (32/36) of the obtained TCGA-CHOL samples were intrahepatic cholangiocarcinoma, and cholangiocarcinoma at different histological sites may present different molecular features, therefore, study with large sample size would be needed in the future to further analyze cholangiocarcinoma samples from different sites. Nevertheless, the 6-gene hypoxia-related signature identified in our study demonstrated the value to better predict the prognostic status of patients and respond to the hypoxic status of tumor cells.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/ supplementary material.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by First Affiliated Hospital of Xi'an Jiaotong University's Institutional Medical Ethics Committee. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
QS, HW, and GW undertook the data analysis and wrote the manuscript. BX and DX helped design the study and edit the manuscript. HW and BX aided in the methodology and data duration. BX and DX assisted in the conception and data representation. All the authors reviewed the data and analysis, read and approved the manuscript.