FXYD3 Expression Predicts Poor Prognosis in Renal Cell Carcinoma with Immunosuppressive Tumor Microenvironment

Simple Summary FXYD3 belongs to the protein-coding gene family associated with Na+/K+-ATPase enzymes and chloride ion channels. Recently, the biological role of FXYD3 has been reported in multiple cancers. Nevertheless, the prognostic value of FXYD3 expression has been undiscovered in clear renal cell carcinoma (KIRC). In this study, we assessed the datasets from The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) dataset (GSE29609). We found the FXYD3 high KIRC patients had distinct clinical characteristics, including hypoxia and poor overall survival. Furthermore, the algorithms discovered that FXYD3 mRNA levels were associated with tumor purity, multiple types of the tumor infiltrating lymphocytes (TILs) and several genes related to T cell exhaustion. In conclusion, FXYD3 predicts a poor prognosis associated with hypoxia, pro-tumor TILs, and T cell exhaustion in KIRC. Abstract FXYD3 is a protein-coding gene, belonging to the FXYD protein family associated with Na+/K+-ATPase enzymes and chloride ion channels. Accumulating evidence suggests the biological role of FXYD3 in multiple cancers. However, the prognostic value of FXYD3 expression in clear renal cell carcinoma (KIRC) is unclear. Therefore, we evaluated the clinical data with tumor-infiltrating lymphocytes (TILs) and immunoinhibitory gene expression data using The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) dataset (GSE29609). First, the FXYD3 high KIRC patients had distinct clinical characteristics, including age, sex, disease stage, histological grade, and hypoxia-related gene expressions. Next, FXYD3 gene expression was correlated with poor overall survival in both TCGA and GSE29609 cohorts. The ESTIMATE algorithm revealed that higher FXYD3 mRNA levels were associated with increased infiltration of immune cells and tumor purity. Moreover, the FXYD3 high KIRC tissue harbored increased TILs such as B cells, CD8+ T cells, and M1 macrophage, whereas NK cells and neutrophils were decreased. In addition, we showed FXYD3 was co-expressed with several immunoinhibitory genes related to T cell exhaustion such as LGALS9, CTLA4, BTLA, PDCD1, and LAG3. In conclusion, FXYD3 is an unfavorable prognostic biomarker associated with hypoxia, pro-tumor TILs, and T cell exhaustion.


Introduction
Renal cell carcinoma (RCC) is the most common type of kidney cancer in adults, of which clear cell RCC (KIRC) accounts for most cancer-related deaths [1]. Clinical challenges in advanced KIRC lie in managing a poor prognosis caused by resistance to radiotherapy and chemotherapy [2]. KIRC harbors an immunogenic tumor microenvironment (TME) containing various tumor-infiltrating T lymphocytes (TILs) [3]. Patients with advanced and metastatic RCC have been treated with molecularly targeted agent monotherapy such as tyrosine kinase inhibitors (TKIs). However, since the advent of immune checkpoint blockade (ICB), the treatment of RCC has reached a significant turning point. ICBs that target programmed cell death1 (PD-1), programmed cell death-ligand 1 (PD-L1), and cytotoxic T-lymphocyte-associated antigen 4 (CTLA-4) or ICB plus TKIs are now standard treatment options for RCC [4,5]. In the era of the ICB, an understanding of the immunogenic TME would be useful to find a new therapeutic strategy in KIRC management.
FXYD3 is a protein-coding gene, belonging to the FXYD protein family in association with Na+/K+-ATPase enzymes and chloride ion channels [6]. FXYD3 is expressed in various organs such as the liver, pancreas, stomach, colon, prostate, lungs, skeletal muscles, and kidneys [7]. In addition, evidence shows that FXYD3 expression can be a poor prognostic biomarker in multiple solid tumors [8][9][10][11][12][13]. However, the prognostic role of FXYD3 in KIRC has not been clarified so far. Therefore, we hypothesized FXYD3 could predict the outcome in KIRC patients. Moreover, if FXYD3 is associated with the TME of KIRC, FXYD3 could be a potential target to modulate in order to increase the response of ICB in KIRC.
Herein, we investigated the prognostic value of FXYD3 gene expressions in KIRC using The Cancer Genome Atlas (TCGA) datasets and Gene Expression Omnibus (GEO) datasets. Moreover, we investigated the tumor microenvironment (TME) by evaluating the estimated TILs [14] and immune-inhibitory gene expressions [15] in association with FXYD3 expression.

Survival Analysis with FXYD3 mRNA Levels in the Public Database
Clinical data of patients with KIRC on TCGA datasets and GSE29609 [16] were obtained from cBioPortal (TCGA PanCancer Atlas; https://www.cbioportal.org/, accessed on 16 December 2020) [17,18], and GEO database, respectively. The patients with an upper 50% expression of FXYD3 mRNA in each dataset were classified as FXYD3 high, and the rest were classified as FXYD3 low in all the following analyses. Kaplan-Meier curves were plotted for overall survival (OS). Data of tumor purity (ESTIMATE score), the level of stromal cells present (stroma score), and the infiltration level of immune cells (immunity score) based on expression data were obtained from ESTIMATE (https: //bioinformatics.mdanderson.org/estimate/disease.html, accessed on 17 July 2022) [19]. FXYD3 expression levels with tumor grade and disease stage were achieved on TISIDB [20]. The multivariate Cox regression analysis was used on Gene Expression Profiling Interactive Analysis 2 (GEPIA2), a web tool for comprehensive analyses with TIL data on TCGA [21,22], to evaluate the prognostic significance of FXYD3 gene expression as a continuous variable with other parameters such as age, sex (male vs. female), and disease stage at initial diagnosis.

Estimated TIL Fraction
The TIL fraction data estimated using quanTIseq [14] was downloaded (http://timer. comp-genomics.org/, accessed on 17 July 2022). The composition of TILs was visualized by principal component analysis. The estimated TIL fractions between FXYD3 high vs. low patients were compared by the Mann-Whitney test. A heatmap of TIL fractions estimated using quanTIseq [14] with clinical information was generated by the R package 'ComplexHeatmap'. Each fraction of TILs was converted to a z-score. We performed hierarchical clustering to characterize the TILs, and the calculating methods of the distance for rows and columns were 'canberra' and 'euclidean', respectively. The clustering methods for rows and columns were 'ward.D2' and 'complete', respectively.

Gene Co-Expression Network Analysis
We used the R package 'igraph' to visualize the co-expression network of FXYD3 and the immune-inhibitory genes [15]. Spearman's correlation coefficients were used to generate the network graph. The nodes represented the genes and the edges corresponded to the interaction between genes. Subgraphs that show densely connected genes in the graph were created by calculating the leading non-negative eigenvector.

Statistical Analysis
For graphs and statistical analyses (principal component analysis, Mann-Whitney test), we used the R freeware (http://www.r-project.org, accessed on 17 July 2022) and GraphPad Prism (v. 9.4.0) software (GraphPad Software Inc., San Diego, CA, USA). All p-values were two-sided, and a p-value of ≤ 0.05 was considered statistically significant. Table 1 shows the clinical characteristics between FXYD3 high and low patients (N = 255, n = 255, respectively). Significant differences of age, sex, disease stage, and histologic grade were detected. In the high FXYD3 group, the median of age and the % of male patients were significantly higher. Moreover, the FXYD3 high patients harbored more hypoxia-related gene expression by Buffa hypoxia score [23]. In contrast, there was no significant difference in mutation counts.

Infiltrating Immune Cells and FXYD3 mRNA Levels
To understand the poor outcome in FXYD3 high KIRC patients, we examined the characteristics of the TME between FXYD3 high and low patients in terms of TILs, tumor purity, and abundance of stromal cells using the ESTIMATE algorithm [19]. We found FXYD3 high KIRC patients had higher levels of infiltration of immune cells and tumor purity (Figure 2A,B) whereas there was no difference of stromal cells in the TME ( Figure 2C). TILs are key factors in regulating tumor progression in KIRC [24]. Next, we compared the TILs estimated using the quanTIseq method [14] between FXYD3 high and low groups. The TIL population composition did not strongly overlap based on FXYD3 expression levels ( Figure 3A). Next, the FXYD3 high group had increased B cells, CD8+ T cells, and M1 macrophages, and decreased NK cells and neutrophils ( Figure 3B,C,E,G,H). Regulatory T cells (Treg) and M2 macrophages were not significantly different ( Figure 3D,F).
Cancers 2022, 14, x FOR PEER REVIEW 5 of 11 FXYD3 high KIRC patients had higher levels of infiltration of immune cells and tumor purity (Figure 2A,B) whereas there was no difference of stromal cells in the TME ( Figure  2C). TILs are key factors in regulating tumor progression in KIRC [24]. Next, we compared the TILs estimated using the quanTIseq method [14] between FXYD3 high and low groups. The TIL population composition did not strongly overlap based on FXYD3 expression levels ( Figure 3A). Next, the FXYD3 high group had increased B cells, CD8+ T cells, and M1 macrophages, and decreased NK cells and neutrophils ( Figure 3B,C,E,G,H). Regulatory T cells (Treg) and M2 macrophages were not significantly different ( Figure  3D,F).  The heatmap of TILs and FXYD3 mRNA levels is shown in Figure 4. Three clusters were identified. In clusters 1 and 3, most patients were FXYD3 high, while most of the FXYD3 high KIRC patients had higher levels of infiltration of immune cells and tumor purity (Figure 2A,B) whereas there was no difference of stromal cells in the TME ( Figure  2C). TILs are key factors in regulating tumor progression in KIRC [24]. Next, we compared the TILs estimated using the quanTIseq method [14] between FXYD3 high and low groups. The TIL population composition did not strongly overlap based on FXYD3 expression levels ( Figure 3A). Next, the FXYD3 high group had increased B cells, CD8+ T cells, and M1 macrophages, and decreased NK cells and neutrophils ( Figure 3B,C,E,G,H). Regulatory T cells (Treg) and M2 macrophages were not significantly different ( Figure  3D,F).  The heatmap of TILs and FXYD3 mRNA levels is shown in Figure 4. Three clusters were identified. In clusters 1 and 3, most patients were FXYD3 high, while most of the those TILs (i.e., M1/M2 macrophage, CD8+ T cells, B cells, and Treg) were relatively fewer than cluster 3. Contrary to clusters 1 and 3, the cluster 2 patients harbored increased NK cells and neutrophils and fewer M1/M2 macrophages, CD8+ T cells, B cells, and Treg. cluster 2 patients had low FXYD3 gene expression. The cluster 3 patients had high fractions of TILs, including M1/M2 macrophages, CD8+ T cells, B cells, and Treg, with fewer NK cells and neutrophils. The cluster 1 patients had similar TIL patterns to cluster 3, though all those TILs (i.e., M1/M2 macrophage, CD8+ T cells, B cells, and Treg) were relatively fewer than cluster 3. Contrary to clusters 1 and 3, the cluster 2 patients harbored increased NK cells and neutrophils and fewer M1/M2 macrophages, CD8+ T cells, B cells, and Treg.

Expression of Immunoinhibitory Genes and FXYD3 mRNA Levels
To discover which immunosuppressive gene signals contribute to a worse prognosis in FXYD3 high patients, we screened the expression of the immunoinhibitory gene set [15]. The heatmap of FXYD3 mRNA levels and immunoinhibitory gene expression data is shown in Figure 5. Three clusters were identified. In column clusters 2 and 3, the majority of patients had higher FXYD3 expression, while column cluster 1 was composed of more FXYD3 low patients. The expression levels of row cluster 1 genes (CD160, KIR2DL1, KIR2DL3, KDR, ADORA2A, IDO1, HAVCR2, NECTIN2, TGFB1, TGFBR1, VTCN1, and IL10RB) were not distinctly different across all column clusters. In contrast, the row cluster 2 gene expressions (CSF1R, IL10, LGALS9, CD244, CTLA4, BTLA, TIGIT, CD96, PDCD1, and LAG3) were higher in column clusters 2 and 3, whereas column cluster 1 rarely expressed those genes. In addition, we evaluated the correlation of FXYD3 gene expression

Discussion
The present study shows the FXYD3 mRNA level is a prognostic factor in KIRC. In the TME of KIRC, we observed the different infiltration of multiple TILs associated with FXYD3 expression levels. Moreover, KIRC patients with higher FXYD3 gene expression had distinct expression patterns of immunoinhibitory genes, suggesting T cell exhaustion was associated with FXYD3.
Our study provides a potential explanation for the poor prognosis in FXYD3 high KIRC. First, increased hypoxia-related signaling was suggested in FXYD3 high groups. Hypoxia is one of the hallmarks in various solid tumors with a critical association with tumor genetic instability and prognosis [36]. In KIRC, most tumors harbor somatic inactivation of both VHL alleles with loss of function of the VHL tumor suppressor protein (pVHL) [36]. The loss of function of pVHL leads to activation of hypoxia-inducible factor (HIF), leading to tumor progression [36]. The potential relevance of FXYD3 to hypoxiarelated genes should be explored using animal models in further studies. Second, we observed that the TIL pattern in FXYD3 high patients was characteristic of an unfavorable

Discussion
The present study shows the FXYD3 mRNA level is a prognostic factor in KIRC. In the TME of KIRC, we observed the different infiltration of multiple TILs associated with FXYD3 expression levels. Moreover, KIRC patients with higher FXYD3 gene expression had distinct expression patterns of immunoinhibitory genes, suggesting T cell exhaustion was associated with FXYD3.
Our study provides a potential explanation for the poor prognosis in FXYD3 high KIRC. First, increased hypoxia-related signaling was suggested in FXYD3 high groups. Hypoxia is one of the hallmarks in various solid tumors with a critical association with tumor genetic instability and prognosis [36]. In KIRC, most tumors harbor somatic inactivation of both VHL alleles with loss of function of the VHL tumor suppressor protein (pVHL) [36]. The loss of function of pVHL leads to activation of hypoxia-inducible factor (HIF), leading to tumor progression [36]. The potential relevance of FXYD3 to hypoxia-related genes should be explored using animal models in further studies. Second, we observed that the TIL pattern in FXYD3 high patients was characteristic of an unfavorable TME containing increased CD8+ T cells [4,37] in KIRC. In the TME of KIRC, highly infiltrated CD8+ T cells paradoxically could not contribute to a better prognosis, likely due to exhaustion of CD8+ T cells [4,38]. Our results showed that disease stage and histological grade were correlated with increased FXYD3 expression. Moreover, we found increased expression of multiple inhibitory genes involved in T cell exhaustion (e.g., PDCD1, CTLA4, TIGIT, LAG3 [39]) in FXYD3 high patients. The results of the network analysis also support the association of FXYD3 with several genes related to T cell exhaustion [39,40] in KIRC. In contrast, the FXYD3 low patients had increased NK cells and neutrophils with fewer T cell subsets such as CD8+ T cells or Treg. High NK cell infiltration in RCC is correlated with a better prognosis [41,42]. Our results suggest that the NK cells inside the TME could contribute to a better prognosis in FXYD3 low patients. A higher neutrophil-lymphocyte ratio is prognostic of poor prognosis in RCC [43]. Neutrophils in the TME have subtypes: antitumor (N1) and pro-tumor (N2) phenotypes [44]. N2 neutrophils can be induced by exposure of the neutrophils to regulatory factors like TGF-β [44]. The present study showed the low expression of TGF-β in FXYD3 low patients, suggesting that infiltrated neutrophils might not gain a pro-tumor phenotype. Taken together, FXYD3 is associated with poor prognosis in KIRC in relation to pro-tumor TME, potentially with hypoxia-related signaling, less NK cells, and exhausted CD8+ T cells.
The acknowledged limitations are that the present study was based only on the public database. In addition, this study did not examine the association between FXYD3 expression and the effects of antitumor agents, including ICBs. A preclinical model study would reveal the causative biological role of FXYD3 signaling in KIRC. Further study is necessary to demonstrate the mechanisms involving FXYD3-related signaling pathways on poor prognosis.

Conclusions
The present study demonstrates that FXYD3 is an unfavorable prognostic biomarker in KIRC with hypoxia, pro-tumor TILs, and multiple genes related to T cell exhaustion. Institutional Review Board Statement: Ethical approval was waived since we used only publicly available data.
Informed Consent Statement: Informed consent statement was waived since we used only publicly available data.

Conflicts of Interest:
The authors declare no conflict of interest.