Assessment of Diagnosis, Prognosis and Immune Infiltration Response to the Expression of the Ferroptosis-Related Molecule HAMP in Clear Cell Renal Cell Carcinoma

Background. Hepcidin antimicrobial peptide (HAMP) is a key factor in maintaining iron metabolism, which may induce ferroptosis when upregulated. However, its prognostic value and relation to immune infiltrating cells remains unclear. Methods. This study analyzed the expression levels of HAMP in the Oncomine, Timer and Ualcan databases, and examined its prognostic potential in KIRC with R programming. The Timer and GEPIA databases were used to estimate the correlations between HAMP and immune infiltration and the markers of immune cells. The intersection genes and the co-expression PPI network were constructed via STRING, R programming and GeneMANIA, and the hub genes were selected with Cytoscape. In addition, we analyzed the gene set enrichment and GO/KEGG pathways by GSEA. Results. Our study revealed higher HAMP expression levels in tumor tissues including KIRC, which were related to poor prognosis in terms of OS, DSS and PFI. The expression of HAMP was positively related to the immune infiltration level of macrophages, Tregs, etc., corresponding with the immune biomarkers. Based on the intersection genes, we constructed the PPI network and used the 10 top hub genes. Further, we performed a pathway enrichment analysis of the gene sets, including Huntington’s disease, the JAK-STAT signaling pathway, ammonium ion metabolic process, and so on. Conclusion. In summary, our study gave an insight into the potential prognosis of HAMP, which may act as a diagnostic biomarker and therapeutic target related to immune infiltration in KIRC.


Introduction
Clear cell renal cell carcinoma (ccRCC), also known as kidney renal clear cell carcinoma (KIRC), is among the most common malignancies worldwide, representing 80% of kidney malignancies [1]. As it is prone to metastasis, the prognosis of KIRC patients is poor, especially for patients in the late clinical stage [2,3]. One study showed that the 5-year overall survival rate (OS) of early stage KIRC could reach 96%, but it was no more than 10% for advanced stages [4]. Although advanced KIRC can be treated with molecular targeted therapy and immunotherapy, the long-term efficacy is still unsatisfactory [3]. Therefore, in order to prolong the overall survival of patients with KIRC there is a pressing need to screen potential diagnostic biomarkers, or identify the potential therapeutic targets related to tumorigenesis, metastasis and prognosis [5,6].
Notably, KIRC cells exhibit substantially higher sensitivity to ferroptosis than normal renal cells, which possess a basal level of ferroptosis sensitivity. The import, export, storage and turnover of iron impact ferroptosis sensitivity [7,8]. Moreover, previous studies have

Data Acquisition
We downloaded the gene expression profiles and clinical information from patients with renal clear cell carcinoma from the Cancer Genome Atlas (TCGA, https://tcga-data. nci.nih.gov/tcga/, accessed on 21 July 2022), a public repository of high-throughput experimental data. The type of data was RNAseq-FPKM, including 539 tumor tissues and 72 adjacent tumorous tissue samples.

Oncomine Database Analysis
The Oncomine database is one of the popular data-mining platforms, which we used to analyze the differential expression of HAMP in different cancer types [24]. In this study, the thresholds of the p-value, fold change and gene rank were, respectively, set as 0.001, 2 and the top 10%.

Timer Database Analysis
The Timer (Tumor Immune Estimation Resource) database is a comprehensive tool mainly used to analyze immune infiltrates by RNA-seq [25]. The website contains over 10,000 samples in 32 different cancer types from TCGA [26]. We studied the differential expression of HAMP between tumor and normal tissues through use of the "Diff Exp" module. In addition, we used the "Correlation" module to analyze the expression scatterplots between HAMP and the immune marker genes in KIRC. Moreover, the "SCNA" module was used to evaluate the differences in tumor infiltration with different copy number alterations for HAMP.

Ualcan
Ualcan is an online resource with an interface that allows the user to analyze the relative expression levels of a query protein or set of proteins across specific tumor subgroups [27]. In the CPTAC confirmatory/discovery datasets, we analyzed the total protein expression levels of HAMP in KIRC.

GeneMANIA Database Analysis
GeneMANIA is a website for generating hypotheses about genes' functions using the available genomic and proteomic data [28]. In our study, HAMP was submitted to the website to select the 50 most closely related genes.

STRING Databases Analysis
STRING is designed to analyze gene-gene and protein-protein interactions, which are shown in a PPI network [29]. We used the STRING database's "multiple proteins" function to construct a PPI network of HAMP and its correlated genes.

GEPIA Database Analysis
GEPIA is a web server that provides expression analysis functions for TCGA and GTEx data to predict genes' correlations in cancer [30]. To pick the closest neighboring genes, we inputted HAMP in KIRC to find and download the 100 most similar genes, which we later used to make an intersection with another group selected by R.

Cytoscape Software
Cytoscape is a software package which is often used to visualize and construct functional networks and make further explorations [31]. We inputted the file downloaded from STRING and formed a PPI, then used the "cytoHubba" app to select the top 10 hub genes.

Gene Set Enrichment Analyis
Gene set enrichment analysis (GSEA) is a computational method that determines whether a predefined set of genes has statistically significant, concordant differences between two biological states [32]. To evaluate whether the selected gene, HAMP, was expressed significantly differentially in KIRC, we divided the samples into the high and low groups according to the differences in the expression and selected the enriched pathways.

Statistical Analysis
The data in this study were mostly analyzed by R 3.6.3, a package downloaded from Bioconductor. We compared the expression of HAMP between tumors and normal tissues in paired and non-paired samples by Wilcoxon's rank-sum test, and the difference showed statistical significance. The differences in the expression levels between Stages T1-T2 and Stages T3-T4 were evaluated through Wilcoxon's rank-sum test. Various factors including sex, stage, age, and so on, were included to analyze the probability of survival by univariate and multivariate Cox regression. Results with p < 0.05 were considered to be significant.
To analyze the prognosis of KIRC, we used Kaplan-Meier curves to analyze the disease-specific survival, overall survival and progression-free survival under low and high expression levels of HAMP. Furthermore, we used the ROC and time-based ROC to evaluate the sensitivity of HAMP in diagnosis. The "Survival" and "survminer" packages in R were used in the survival analysis, the pROC and timeROC packages were used for the ROC analysis, and the ggplot2 package was used to visualize the results.
Next, we used the data from CIBERSORT to create a violin plot with the R package vioplot. We also compared the 22 types of immune cells with each other to make a heatmap to select the most correlated cells with the R package pheatmap. To analyze the Pearson correlations between genes and immune cells, the R programming ssGSEA was used to calculate the correlations between the expression levels of HAMP and tumor-infiltrating immune cells. 4 of 17 To find the correlated genes and construct a PPI network, we first applied a singlegene correlation screening, and created an intersection with the similar genes selected from GEPIA to find the superposition, which established a foundation for further study. Meanwhile, with the coincident gene group, we drew a co-expression heatmap with the R.ggplot 2 package.

Differential Expression of HAMP in KIRC
Given that HAMP may be a potential biomarker for KIRC, we first used the Oncomine database to analyze the mRNA expression levels of HAMP in multiple types of cancer and normal tissues. Significantly, higher expression levels of HAMP were found in brain, breast, color, kidney and lung cancers compared with the adjacent normal tissues ( Figure 1A).

The Clinical Correlation and Prognostic Value of HAMP in KIRC
To assess the correlation between the expression of HAMP and clinical pathologic outcomes in tumors, we explored the pathological stage of KIRC patients. The results revealed that the gene HAMP was significantly upregulated in Stages T3-T4 compared with Stages T1-T2 ( Figure 2A).
In addition, we constructed a nomogram to visualize the prognostic model of the Cox regression analysis and found that the concordance was 0.762 (0.738-0.786) ( Figure 2B). Further, we performed univariate and multivariate Cox regressions to find the correlations between HAMP and age, stage, etc. (Supplementary Table S1).
Furthermore, as Figure 1C shows, the median expression level of HAMP in tumors is much higher than in normal tissue. To eliminate other mixing factors we analyzed paired samples, and the trend was consistent and maintained a significant difference ( Figure 1D).
Since epigenetic alterations have been proven to play a role in cancer biology, we also explored the promoter methylation level of HAMP in KIRC and normal samples with UALCAN. Interestingly, the expression level in primary tumors was significantly lower compared with normal tissues, and the normal-vs.-primary statistic was 1.62 × 10 −12 ( Figure 1E).

The Clinical Correlation and Prognostic Value of HAMP in KIRC
To assess the correlation between the expression of HAMP and clinical pathologic outcomes in tumors, we explored the pathological stage of KIRC patients. The results revealed that the gene HAMP was significantly upregulated in Stages T3-T4 compared with Stages T1-T2 ( Figure 2A).
In addition, we constructed a nomogram to visualize the prognostic model of the Cox regression analysis and found that the concordance was 0.762 (0.738-0.786) ( Figure 2B). Further, we performed univariate and multivariate Cox regressions to find the correlations between HAMP and age, stage, etc. (Supplementary Table S1).
To evaluate the prognostic significance of HAMP in KIRC, we compared the patients with high and low expression levels in terms of overall survival (OS), disease-specific survival (DSS) and progression-free interval (PFI), and the results were HR = 1.53 (1.55-3.51; p < 0.001) for DSS, HR = 1.82 (1.33-2.48; p < 0.001) for OS and HR = 1.68 (1.22-2.32; p = 0.001) for PFI. The outcomes indicate that HAMP is a dangerous factor for KIRC, which plays an important role in prognosis ( Figure 2C-E).
Furthermore, we used the ROC to construct a model to predict the probability of survival and found that the AUC was 0.911 (0.879-0.944) in the ROC curve, which showed that HAMP is a significant predictive index in diagnosis ( Figure 2F). Similarly, the AUC in the ROC for OS was 0.649 (0.601-0.698) and the AUC of the time-dependent ROC predicted that the 1-year, 3-year and 5-year outcomes were, respectively, 0.686, 0.631 and 0.627 ( Figure 2G,H).

Correlation between Immune Cells and HAMP Expression Levels in KIRC
As we know, tumor-infiltrating immune cells are one of the representative cellular components of the host's anti-tumor immune responses and tumor immune escape [33]. Therefore, we explored the correlation between HAMP expression levels and immune infiltration levels to confirm the effect of HAMP in prognosis in KIRC. According to the HAMP expression levels, we divided the samples into the tumor and normal expression groups and analyzed the differences in 22 immune cells. The results revealed that CD8T cells, CD4 memory resting T cells and M2 macrophage cells were actively expressed, and regulatory T cells (Treg) and neutrophils were closely related to the HAMP expression levels. Interestingly, regulatory T cells (p = 0.012) were much higher in the normal tissue group than in the tumor group, while neutrophils (p = 0.015) showed a completely opposite pattern ( Figure 3A). Furthermore, we used the ROC to construct a model to predict the probability of survival and found that the AUC was 0.911 (0.879-0.944) in the ROC curve, which showed that HAMP is a significant predictive index in diagnosis ( Figure 2F). Similarly, the AUC in the ROC for OS was 0.649 (0.601-0.698) and the AUC of the time-dependent ROC predicted that the 1-year, 3-year and 5-year outcomes were, respectively, 0.686, 0.631 and 0.627 ( Figure 2G,H).  groups and analyzed the differences in 22 immune cells. The results revealed that CD8T cells, CD4 memory resting T cells and M2 macrophage cells were actively expressed, and regulatory T cells (Treg) and neutrophils were closely related to the HAMP expression levels. Interestingly, regulatory T cells (p = 0.012) were much higher in the normal tissue group than in the tumor group, while neutrophils (p = 0.015) showed a completely opposite pattern ( Figure 3A).  In addition, a correlation heatmap was used to visualize the relational degree within the subpopulation of immune cells. CD8 T cells were the most negatively related to CD4 memory resting T cells and closely positively related to regulatory T cells ( Figure 3B).
We also used a lollipop chart to clearly compare the degree of correlation among the 24 immune cells and selected the most representative immune cells for further exploration. Moreover, we used the TIMER database to analyze the correlation between the expression levels of HAMP and the representative immune markers ( Table 1). As shown in Figure 4A, the top 8 immune cells selected in Figure 3C with a Pearson correlation of >0.4 all had a linear relationship with the expression levels of HAMP ( Figure 3C). To further find the mechanism of how HAMP expression acts on the immune cells, we drew scatterplots between HAMP expression levels and the immune marker genes of M1/M2 macrophages, T cells and B cells in KIRC based on the TIMER database. The results showed statistical significance ( Figure 4B-E). Moreover, the correlations between HAMP and the related marker genes of B cells, T cells and macrophages in KIRC are listed in Table 2.
After we had analyzed the somatic copy number alterations (SCNAs) in KIRC, there were significant differences in B cells, CD8+ T cells, CD4+ cells, neutrophils and dendritic cells, which confirmed the prediction that HAMP may regulate the tumor process through immune infiltration ( Figure 4G). We also assessed the correlations between the most common immune checkpoints and the expression levels of HAMP to further strengthen the potential mechanism between the expression levels of HAMP and immune cells ( Figure 4F).

The Network of HAMP Expression in KIRC
In cancer research, constructing a PPI network is a useful method for revealing coexpressed and related genes [34]. We constructed a PPI network of HAMP and another 50 interacting proteins in GeneMANIA. (Figure 5A) Moreover, we used the GEPIA database and R programming to find the top 100 similar genes of HAMP and constructed a Venn diagram to make an intersection ( Figure 5C). There were 39 co-expressed genes, as shown in Supplementary Table S2. The interactions of the 39 genes and HAMP were visualized with STRING.  Table 2.   In addition, we plotted a heatmap of HAMP's co-expressed genes to analyze the correlations, and the results showed a highly consistent trend, which means that these 39 co-expressed genes have significant positive correlations with the key gene, HAMP ( Figure 5D).
We used the "cytoHubba" module on the PPI network to calculate the nodes' scores and selected the top 10 hub nodes as ranked by MCC. The key genes extracted and the shortest paths between the hub genes are displayed in Figure 5E. lation and T cell receptor signaling pathways, and six negatively related pathways: aut immune thyroid disease, the B cell receptor signaling pathway, Huntington's disease, t JAK-STAT signaling pathway, non-small cell lung cancer, pancreatic cancer and small c lung cancer ( Figure 6C). Two GO items, namely the ammonium ion metabolic process an beta catenin binding, showed an upregulation of HAMP according to the FDP q-valu and NOM p-values ( Figure 6D).

Enrichment Analysis of HAMP in KIRC
As GSEA has been extensively utilized in the context of differential expression analysis, we also analyzed the intersection genes described above to uncover the likely functional interpretations [35]. Gene Ontology (GO) is commonly used for describing our knowledge of the biological domain with respect to three aspects: molecular function, biological processes, and location in the cellular component [36]. In total, 269 GO terms (251 terms for BP, 10 terms for CC and 8 terms for MF) were enriched, with an adjusted p-value of <0.05 and a q-value of <0.2. For biological process, regulation of lymphocyte activation, synapse pruning, and phagocytosis were enriched. For the cellular component, tertiary granules, collagen trimers and secretory granule membranes were the main locations. In the MF category, amyloid-beta binding, immunoglobulin binding and IgG binding were enriched. In addition, KEGG enrichment analysis was applied to find some pathways where the genes were enriched, such as Fc gamma R-mediated phagocytosis, Staphylococcus aureus infection, osteoclast differentiation, and so on ( Figure 6A). Based on the GO and KEGG pathway enrichment analyses, the networks of the genes, terms and their interactions were visualized ( Figure 6B).
In addition, based on the groups with low and high expression levels of HAMP in KIRC, we explored the HAMP-related signaling pathways in the GO and KEGG enrichment analyses, which showed significant differences (FDR < 0.05 and p < 0.05). The selected pathways in the groups with high and low expression are listed in Supplementary Tables S3 and S4. The KEGG pathway analysis revealed two pathways that were positively correlated with high expression levels of HAMP, namely the oxidative phosphorylation and T cell receptor signaling pathways, and six negatively related pathways: autoimmune thyroid disease, the B cell receptor signaling pathway, Huntington's disease, the JAK-STAT signaling pathway, non-small cell lung cancer, pancreatic cancer and small cell lung cancer ( Figure 6C). Two GO items, namely the ammonium ion metabolic process and beta catenin binding, showed an upregulation of HAMP according to the FDP q-values and NOM p-values ( Figure 6D).

Discussion
HAMP is closely related to ferroptosis, which is a key factor of the iron overload disease hemochromatosis [37]. It is reported ubiquitously expressed in cancers and can, to some extent, be a cancer-driving gene [38]. We found that HAMP was highly expressed in KIRC, KIRP, LUAD, SKCM, and so on, which corresponds to Liu's former study [13]. Significant differences between normal and tumor tissues in KIRC were explored. In a previous study, KIRC patients had a poor prognosis, especially those in advanced T stages [3]. Our study has proven the overexpression of HAMP in later stages of the disease. We also observed that the expression level of HAMP is closely related to the probability of survival and the higher the expression, the worse the DDS, OS and PFS. In addition, our results also showed the mechanism of how HAMP affects prognosis, namely the interactions among the expression level of HAMP and the degree of immune infiltration and different immune markers. Based on all these results, HAMP can be considered as a potential diagnostic and prognostic biomarker and could be an immune-related treatment target.

Discussion
HAMP is closely related to ferroptosis, which is a key factor of the iron overload disease hemochromatosis [37]. It is reported ubiquitously expressed in cancers and can, to some extent, be a cancer-driving gene [38]. We found that HAMP was highly expressed in KIRC, KIRP, LUAD, SKCM, and so on, which corresponds to Liu's former study [13]. Significant differences between normal and tumor tissues in KIRC were explored. In a previous study, KIRC patients had a poor prognosis, especially those in advanced T stages [3]. Our study has proven the overexpression of HAMP in later stages of the disease. We also observed that the expression level of HAMP is closely related to the probability of survival and the higher the expression, the worse the DDS, OS and PFS. In addition, our results also showed the mechanism of how HAMP affects prognosis, namely the interactions among the expression level of HAMP and the degree of immune infiltration and different immune markers. Based on all these results, HAMP can be considered as a potential diagnostic and prognostic biomarker and could be an immune-related treatment target.
We first analyzed the mRNA expression levels of HAMP in different cancers with Oncomine and R programming based on the TCGA database. There were discrepancies in the differential expression levels of HAMP in different cancer types among the different databases, which may reflect the differences in the data gathering methods and the fundamental mechanisms in terms of biological characteristics. With the same results in different methods, we can safely draw the conclusion that HAMP plays a significant role in cancers. In addition, it can be said that methylation in the promoter genes related to tumors is significantly connected to the clinical behavior of cancers [39,40]. Our UALCAN results showed the lower promoter methylation level of HAMP in KIRC, further confirming HAMP as a dangerous factor for epigenetic alteration.
Moreover, since the mRNA and protein expression levels of HAMP varied between normal and tumor tissues in KIRC, it showed high prognostic potential. In terms of prognosis, the T stage is one of the well-known prognostic factors, and a higher T stage indicates poor prognosis [41]. Upregulated expression levels of HAMP are related to the clinical stage and are highly likely to be a key factor in a lower probability of survival. Our study also found that a high expression level of HAMP is related to the diagnosis of KIRC. However, further steps are needed to explore the specific mechanisms and find strong support.
As we know, the immune system plays a major role in the control and progression of cancer [37]. Some studies have also reported that immune infiltration has an influence on the therapeutic responsiveness and prognosis of KIRC patients [42]. Thus, we have provided insight into the relationship between the expression levels of HAMP and immune cell infiltration in KIRC. In this study, the results showed that the expression level of HAMP has a converse correlation with the infiltration level of regulatory T cells and neutrophil. Former studies have shown that the neutrophil recruitment/activation results in tumor cell death, while Treg cells are recruited into the tumor microenvironment to mediate immune suppression, which is the reason for their positive and negative relationship, respectively, with the expression level of HAMP [43,44]. The results also showed that the expression level of HAMP has a strong positive correlation with the extent of infiltration of macrophages, Treg, Th1 cells, T cells and B cells. Moreover, the relationship between the expression level of HAMP and various immune cells' gene markers indicates that HAMP could regulate the tumor immune microenvironment. The gene markers of M1 macrophages (NOS2, IRF5 and PTGS2) had a weak relationship with HAMP expression, whereas those of M2 macrophages were strongly correlated with HAMP, which predicted the regulatory function of HAMP in TAM polarization. In addition, the infiltration of a large number of Treg cells into tumor tissues is often associated with poor prognosis [45]. The increased expression of HAMP was found to have a significant correlation with Treg cells' gene markers in KIRC. These findings suggested that HAMP plays a role in the recruitment and functioning on the immune infiltrating cells in KIRC.
After we constructed a PPI network of HAMP's co-expressed genes in STRING and found their intersections in R, the collected 48 genes showed a strong relationship with HAMP. In addition, the selected 10 hub genes revealed a high degree of interaction, thus indicating potential biomarkers of KIRC. Based on the GO and KEGG pathway analyses, we found that HAMP was significantly associated with several immune response pathways and cancer pathways. Our study showed that HAMP is related to the biological processes of regulating lymphocyte activation, synapse pruning and phagocytosis, which correspond to the tumor immune environment. Moreover, the GO and KEGG pathways revealed the differentially enriched gene sets between the groups with high and low expression levels. The JAK-STAT signaling pathway and the related pathways were associated with increased mortality rates in renal cancer [46].
However, there are some defects and imperfections in our experiment. Although we used multiple databases for a comprehensive analysis, there were differences in the algorithms among the databases which may have led to deviations in the results. Similarly, as the data came from TCGA, we could not exclude whether the patients' survival and other indicators were interfered with by external factors such as drug use. In addition, the specific mechanisms underlying the role of genes in cancer diagnosis remain to be explored, as the relationship between immune infiltration markers and gene expression is not sufficient to support this hypothesis. Nevertheless, our experiment is the first to investigate the effect of a high expression level of HAMP in KIRC on the prognosis of tumor patients, and it provides a possible target for clinical treatment.

Conclusions
In summary, our study has provided an insight into the potential prognostic use of HAMP, and thus it may act as a diagnostic biomarker and therapeutic target related to immune infiltration in ccRCC.
We found through bioinformatic analyses that the ferroptosis-related gene HAMP is a key gene in the KIRC. The expression of HAMP significantly increased in KIRC, and the evaluated expression levels of HAMP played a potentially important role in the diagnosis and prognosis of KIRC. Additionally, HAMP expression is positively related to immune infiltration in KIRC. Moreover, we found 39 genes that are similar to HAMP and selected 10 hub genes, which function in several KEGG and GO pathways.