Identification of SOX6 and SOX12 as Prognostic Biomarkers for Clear Cell Renal Cell Carcinoma: A Retrospective Study Based on TCGA Database

Background The SOX gene family has been proven to display regulatory effects on numerous diseases, particularly in the malignant progression of neoplasms. However, the molecular functions and action mechanisms of SOX genes have not been clearly elucidated in clear cell renal cell carcinoma (ccRCC). We aimed to explore the expression status, prognostic values, clinical significances, and regulatory actions of SOX genes in ccRCC. Methods RNA-sequence data and clinical information derived from The Cancer Genome Atlas (TCGA) database was used for this study. Dysregulated SOX genes between the normal group and ccRCC group were screened using the Wilcoxon signed-rank test. The Kaplan-Meier analysis and univariate Cox analysis methods were used to estimate the overall survival (OS) and disease-specific survival (DSS) differences between different groups. The independent prognostic factors were identified by the use of uni- and multivariate assays. Subsequently, the Wilcoxon signed-rank test or Kruskal-Wallis test and the chi-square test or Fisher exact probability methods were employed to explore the association between clinicopathological variables and SOX genes. Finally, CIBERSORT was applied to study the samples and examine the infiltration of immune cells between different groups. Results Herein, 12 dysregulated SOX genes in ccRCC were screened. Among them, two independent prognostic SOX genes (SOX6 and SOX12) were identified. Further investigation results showed that SOX6 and SOX12 were distinctly associated with clinicopathological features. Furthermore, functional enrichment analysis revealed that SOX6 and SOX12 were enriched in essential biological processes and signaling pathways. Finally, we found that the SOX6 and SOX12 expression levels were correlated with tumor-infiltrating immune cells (TIICs). Conclusion The pooled analyses showed that SOX6 and SOX12 could serve as promising biomarkers and therapeutic targets of patients with ccRCC.


Introduction
Clear cell renal cell carcinoma (ccRCC) is the most common aggressive kidney malignant tumor, whose incidence is increasing year by year [1]. In spite of great advances in cancer treatments, the prognosis of advanced ccRCC is still poor. Currently, the reliable prognostic biomarkers for ccRCC are still limited [2,3]; novel robust diagnostic or prognostic biomarkers are needed to improve patient prognostication.
The SOX gene family, containing more than 20 family members, are thought to be involved in diverse biological processes; deregulation of the SOX gene can disrupt gene expression, and this has been linked to many diseases, particularly in cancer progression [4]. Emerging evidence indicated that the SOX family exhibited regulatory functions in tumor development, including tumor growth, altered consistency of cancer microenvironment, and invasion [5][6][7]. Several SOX genes were proven to be able to screen several tumor specimens and predict clinical outcome of tumor patients. For instance, SOX1 has been proven to display a regulatory function in the invasion and migration of gastric cancer cells [8]. Elevated SOX2 expressions were involved in enhanced stemness of various neoplasms and acquired resistance to chemotherapeutic agents [9,10]. In breast cancer, upregulation of SOX3 suppressed proliferation and stimulate apoptosis of cancer cells [11]. The SOX4 gene was elevated in >22 types of malignant tumors, and it is widely considered as an oncogene [12]. In gastric cancer, hypoxia could lead to the activation of SOX5/Wnt pathway via inhibiting miRNA-338-3p, subsequently affecting tumor risk [13]. In general, it has been demonstrated that the SOX gene family played critical roles in the pathogenesis of various cancers. Exploring the expression patterns, prognostic values, clinical significances, and action mechanisms of the SOX gene in ccRCC offers new opportunities for targeted therapies of ccRCC. In the current study, we analyzed TCGA database to explore the role of SOX genes in ccRCC.

Materials and Methods
2.1. Data Collection. RNA-seq data and clinical information were obtained from TCGA (https://portal.gdc.cancer.gov/) [14], including 530 ccRCC samples and 72 noncancerous samples' transcriptome data and detailed clinical information of the corresponding patients (Table 1).

Identification of Dysregulated SOX Genes.
The "limma" R package [15] and Wilcoxon test were applied for the preliminary identification of the dysregulated SOX genes in tumor specimens. p < 0:05 was considered as differentially expressed. The "pheatmap" R package was utilized to draw the heat map of the dysregulated SOX genes.

Survival Analysis and Uni-and Multivariate Assays.
To evaluate the prognostic value of these dysregulated SOX genes, OS and DSS were estimated using the univariate Cox analysis and Kaplan-Meier (KM) analysis methods. The SOX genes significantly associated with OS and DSS were considered as prognosis-related SOX genes. Then, combined with the RNA-seq data and clinical features, the independent prognostic SOX genes and clinical variables were determined using Cox assays. These independent prognostic SOX genes were regarded as candidate SOX genes.

Functional Enrichment Analysis.
To delve into the action mechanisms of the candidate SOX genes, we divided all samples into high-and low-expression subgroups based on the median value of the gene expression. Then, the differentially expressed genes (DEGs) between two subgroups were confirmed by the use of the "limma" R package. Genes that meet the following criteria were considered as DEGs: the absolute value log2FoldChange (FC) is greater than one, and FDR is <0.05. Metascape (http://metascape.org) was applied to study function and pathway enrichment analyses of DEGs [16].
2.5. Evaluation of Tumor-Infiltrating Immune Cells. TIIC proportions of ccRCC samples were evaluated using the CIBERSORT algorithm [17]. Then, the TIIC proportions of tumor samples were divided into two subgroups (high and low) and visualized by the use of violin plots.
2.6. Statistical Analysis. All statistical analyses were carried out using R 3.6.1 software. The OS and DSS were estimated using Cox regression assays and the Kaplan-Meier method. Univariate and multivariate assays were carried out to identify the independent prognostic SOX genes and clinical features. The chi-square test or Fisher methods were employed to explore the association between clinicopathological variables and the candidate SOX genes. p < 0:05 was regarded as statistically significant.

Results
3.1. Identification of Dysregulated SOXs. Using the "limma" R package and Wilcoxon signed-rank test, a total of 12 dysregulated SOX genes in ccRCC were identified. Figure 1(a) revealed the pheatmap of all differentially expressed SOXs. Among them, SOX4/5/6/13/15/30 were distinctly lower in ccRCC specimens than in noncancerous specimens    3.2. Identification of Prognosis-Related SOX Genes. Using the Cox regression analysis method, we found that SOX6 and SOX13 were protective factors for patients' OS (hazard ratio ðHRÞ < 1), while SOX12 and SOX15 were risk-associated factors (HR > 1) (Figure 2(a)). For patients' DSS, SOX6, SOX7, and SOX13 were protective factors, whereas SOX11, SOX12, and SOX15 were risk-associated factors (Figure 2   Disease Markers and SOX13 were protective factors for patients' OS and DSS, low SOX6 and SOX13 expression predicted worse prognosis of patients, while SOX12 was risk-associated factor patients' OS and DSS, and high SOX12 expressions were remarkably associated with adverse prognosis of patients. 3.3. Identification of Independent Prognosis Factors. As shown in Table 2, in univariate analysis, SOX6 expression (HR = 0:588, p < 0:001), SOX12 expression (p < 0:001), and SOX13 expression (p < 0:01) were associated with poorer OS of ccRCC patients. Moreover, multivariate assays suggested that SOX6 expression (HR = 0:75, p < 0:05) and SOX12 expression (HR = 1:10, p < 0:05) of patients with ccRCC were independently associated with shorter OS (Figures 4(a) and 4(b)). However, no statistical difference was found for SOX13 ( Figure 4(c)).   Table 3, the higher the histological grade, M status, T status, and clinical stage, the lower the SOX6 expression level, while the high expression level of SOX12 indicated advanced histological grade, M status, T status, and clinical stage (Table 4).

Identification of Potential Molecular Mechanisms.
To explore the action mechanisms of SOX6/12, we identified DEGs between SOX6/12 high-and low-expression subgroups based on the mean expression of SOX6/12 in all samples. Figure 6(a) showed the volcano plots of DEGs based on SOX6 expressions. Then, we performed functional enrichment assays of these DEGs to search out the potential biological processes (BP) and signaling pathways. In terms of BP, SOX6 was mainly involved in the acute inflammatory response, metal ion homeostasis, skin development, extracellular matrix organization, and negative regulation of peptidase activity ( Figure 6(b)). KEGG pathway analysis showed that SOX6 was significantly associated with several classical cancer-related signaling pathways, such as cytokine-cytokine receptor interaction, IL-17 pathway, HIF-1 signaling pathway, proteoglycans in cancer, and PI3K-Akt signaling pathway ( Figure 6(c)). Figure 6(d) showed the volcano plots of DEGs based on SOX12 expressions. SOX12 was mainly involved in endocrine system development, chemical synaptic transmission, signal release, behavior, and regulation of hormone levels ( Figure 6(e)). Moreover, SOX12 was remarkably enriched in neuroactive ligand-receptor interaction, cholinergic synapse, maturity-onset diabetes of the young, synaptic vesicle cycle, calcium signaling pathway, estrogen signaling pathway, cardiac muscle contraction, and steroid hormone biosynthesis ( Figure 6(f)).

Discussion
Herein, we found that 12 SOX genes were dysregulated in ccRCC; among them, SOX4/5/6/13/15/30 were distinctly lowly expressed in ccRCC samples. SOX7/9/11/12/18/21 were highly expressed in ccRCC. Survival assays revealed that low expressions of SOX6/13 predicted a poor OS and DSS, while higher SOX12 expression indicated a worse prognosis. Uni-and multivariate assays suggested that low SOX6 or high SOX12 expression was an independent prognostic factor for poor overall survival of patients with ccRCC. Growing evidence has demonstrated that SOX6 serve as a tumor suppressor in the onset and progression of human cancer. Chen et al. [18] demonstrated that SOX6 was downregulated in cervical cancer, and decreased SOX6 expression significantly stimulated the proliferation, migration, and invasion of tumor cells. Jiang et al. [19] showed that SOX6 was decreased in patients with pancreatic cancer, and decreased SOX6 expression could promote proliferation and metastasis of tumor cells. SOX6   Disease Markers downregulation was also found in prostate cancer, hepatocellular carcinoma, lung cancer, and breast cancer, while its downregulation was associated with the malignant phenotype of neoplasms [20][21][22][23]. SOX12 has been demonstrated to exhibit a regulatory function in carcinogenesis and cancer progression. Du et al. [24] indicated that SOX12 was highly expressed in colorectal cancer and its overexpression indicated a poorer prognosis of patients; mechanistically, SOX12 stimulated tumor cell proliferation and metastasis by modulating asparagine synthesis. Also, SOX12 was significantly increased in breast cancer, and overexpression of SOX12 predicted a poor clinical outcome of patients [25]. Zou et al. [26] found that SOX12 may serve as a cancer stem-like cell marker in hepatocellular carcinoma, which was associated with the chemore-sistance of tumor cells. In lung cancer, SOX12 was found to be upregulated in cancerous samples, and loss of SOX12 in tumor cells suppressed proliferation, migration, and invasion in vitro but stimulated apoptosis of tumor cells [27]. Although SOX6/12 play a crucial role in carcinogenicity and progression, much remains unknown about their roles in ccRCC. In our study, SOX6 was found to be downregulated in ccRCC, while SOX12 was upregulated. The low SOX6 or high SOX12 expression indicates a poor clinical outcome and prognosis. The results were consistent with what was reported in the literature, indicating that SOX6 may serve as a tumor suppressor, while SOX12 may serve as an oncogene in ccRCC. Then, to explore the molecular mechanisms involved in SOX6 and SOX12 in RCC, we performed function and      13 Disease Markers pathway enrichment analyses between high-and lowexpression groups; we found that SOX6 and SOX12 might regulate the ccRCC progression through several cancerrelated signaling pathways, such as cytokine-cytokine receptor interaction, IL-17 signaling pathway, and PI3K-Akt signaling pathway. These signaling pathways might serve as key points to figure out the potential action mechanisms of SOX6/12 in RCC. Immune cells are an essential immune microenvironment component, which are closely associated with tumor occurrence and progression [28,29]. Therefore, we estimated the association between SOX6/12 expression and TIICs in ccRCC applying the CIBERSORT algorithm for the first time. Interestingly, we found that low SOX6 expression or high SOX12 expression was distinctly associated with higher infiltrating levels of regulatory T cells (Tregs), the essential regulators of immune tolerance [30]. Our data revealed that SOX6 and SOX12 might play a specific role in the immunosuppressive tumor microenvironment by regulating regulatory T cells (Tregs).

Conclusions
Our findings validated that several SOX genes are abnormally expressed in ccRCC and distinctly associated with the outcome of ccRCC patients. Besides, SOX6 and SOX12 are expected to be the most promising therapeutic targets in ccRCC treatment.

16
Disease Markers