Prognostic Value of SLC4A4 and its Correlation with Immune Infiltration in Colon Adenocarcinoma

Background

The prevalence of colorectal cancer is high, with estimates indicating that over 100,000 patients in the United States were diagnosed with colon cancer in 2019 [1]. Many risk factors affect the development and progression of colorectal cancer, including smoking, obesity, a low-fiber diet, and alcohol [2]. In addition, the incidence rate of colon cancer is higher in older than in younger people. Despite significant advances in individualized treatment, colon cancer is highly invasive and can metastasize to other organs, resulting in high mortality rates [3,4]. Moreover, colon cancer can bypass lymph nodes, making it more difficult to assess distant metastasis [5].

Most colorectal cancers are adenocarcinomas, with this disease showing significant heterogeneity. Cancer immunotherapy has been shown effective in colon cancer patients with DNA mismatch repair defects [6]. In addition, the immune microenvironment has been found to affect the occurrence and development of colon cancer [7,8]. Traditional categorical variables, such as tumor stage, make predicting the prognosis of patients difficult, as these categories lack individualization [9]. Several risk gene signatures have been shown effective in predicting the prognosis of patients with colon adenocarcinoma [10–12]. These risk genes were only used to build prognostic models, and their relationship with immune components in colon adenocarcinoma is unclear. Therefore, identifying new prognostic biomarkers and their association with immune characteristics is critical for better understanding of colon adenocarcinoma.

The solute carrier family 4 member 4 (SLC4A4) gene encodes a Na+/HCO3− cotransporter, which is mainly involved in the secretion and absorption of sodium bicarbonate. This process is highly important in maintaining the dynamic pH equilibrium within cells. SLC4A4 and other solute carrier family members have been found to contribute to tumorigenesis and development. For example, SLC4A4 or SLC4A9 have been reported to reduce intracellular pH and suppress tumor growth. The levels of expression of SLC4 and SLC26 family transporters have been shown to differ in colon cancer and normal colon tissue, suggesting that the dysregulation of bicarbonate transporters may be diagnostic of tumor development and may be a therapeutic target in cancer treatment [13]. MiR-223-3p decreases the expression of SLC4A4 in clear cell renal cell, enhancing tumor cell proliferation and metastasis [14]. SLC4A4 expression was also shown to be downregulated in thyroid cancer, suggesting diagnostic efficacy [15,16]. However, SLC4A4 is not downregulated in all tumors. For example, SLC4A4 expression was shown to be higher in mucinous epithelial ovarian cancer and chronic myeloid leukemia than in adjacent normal tissue, suggesting that the biological processes in which SLC4A4 is involved are tumor specific [17]. Gene signatures, including SLC4A4, can indicate the prognosis of patients with colorectal cancer [18–20]. Although SLC4A4 participates in the pathologenesis of colon cancers, thes association between its expression and the tumor immune microenvironment and classical mutant genes remain incompletely understood.

This study was performed to assess the association between SLC4A4 expression and patient prognosis by high-throughput sequencing of colon adenocarcinoma and normal colon samples acquired from The Cancer Genome Atlas (TCGA). Relative differences in SLC4A4 expression tissue and their correlations with overall survival were verified in two independent cohorts, GSE41258 and GSE32323. The associations between SLC4A4 expression and major clinical features and prognosis in patients with colon adenocarcinoma were analyzed statistically. The expression of SLC4A4 in other tumors was also examined to determine whether its expression has a general or essential effect on tumorigenesis. To better understand the relationship between SLC4A4 expression and gene mutations, the correlations between SLC4A4 expression and the 20 most frequently mutated genes in colon adenocarcinoma were calculated. Gene set enrichment analysis (GSEA) was performed based on gene ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) databases. The associations between SLC4A4 expression and tumor-infiltrating immune cells (TIICs) were examined using of two web servers: tumor immune estimation resource (TIMER) and CIBERSORT. Finally, changes in SLC4A4 expression and their associations with patient prognosis were verified by gene expression profiling interactive analysis (GEPIA) and the Human Protein Atlas (HPA).

Material and Methods

CLINICAL INFORMATION AND EXPRESSION PROFILE:

Sequencing data on colon adenocarcinomas and corresponding clinical information in the TCGA were downloaded from the UCSC Xena database (https://xenabrowser.net/datapages/) [21], with org.Hs.eg.db (R package) used to annotate Ensembl ID. Data from some patients were excluded due to the absence of information on overall survival (OS), gender, age; tumor status, stage, or location; or lymph node or distant metastasis. Two independent cohorts (GSE41258 and GSE32323) were downloaded from the GEO database. GSE32323 contains 17 pairs of cancer and non-cancerous tissues, with gene expression measured by Affymetrix HG-U133 Plus 2.0 arrays. GSE41258 contains 390 samples, with genes measured by Affymetrix Human Genome U133A arrays. Primary tumor tissue from 153 patients with OS <120 months were selected and used in further analysis.

IMMUNE INFILTRATES ANALYSIS:

The association between SLC4A4 expression and the level of immune infiltration was assessed using the gene module in TIMER, a web server used for comprehensive assessment of levels of immune cell infiltration across different types of cancer [22]. TIMER uses a deconvolution method on high-throughput sequencing or microarray data to calculate the abundance of TIICs, including B cells, CD4+ T cells, CD8+ T cells, neutrophils, macrophages, and dendritic cells. High-throughput sequencing data were evaluated using the CIBERSORT database, which assessed the compositional differences of TIICs between groups with high and low SLC4A4 expression. CIBERSORT uses the LM22 gene signature to classify 22 hematopoietic cell phenotypes into 6 categories, including myeloid cells, B cells, T cells, dendritic cells, macrophages and natural killer cells [23]. CIBERSORT normalized the estimates of the 22 TIICs in each sample to a sum of 1, allowing comparisons across samples and across immune cell types. The Mann-Whitney U test was performed to determine whether TIICs belong to different groups. In addition, Pearson’s correlation analysis between TIICs was performed to determine their possible associations. P<0.05 was defined as statistically significant.

GENE SET ENRICHMENT ANALYSIS:

GSEA was performed on normalized mRNA data using clusterprofiler (R package) [24]. The samples were separated into two groups based on the expression of SLC4A4, and the expression of other genes in the two groups analyzed using the R package, edgeR [25]. Fold changes (FCs) were calculated to show the magnitude of change in expression between the two groups. The gene list imported in into clusterprofiler was sorted in descending order by log2 FC. Permutations were set at 1000, the minimum size of each gene set for analysis at 10 and the maximum number size of genes annotated for testing at 1000. All entries in the GO database and the KEGG pathway were analyzed to investigate the biological functions in which SLC4A4 was involved. Entries with P<0.05 were considered statistically significant.

STATISTICAL ANALYSES:

Logistic regression analysis was performed to determine whether SLC4A4 can distinguish among patients with different clinical characteristics (age, gender, tumor status, lymph node, metastasis, stage); odds ratio (OR) and 95% confidence intervals (CIs) were calculated. Univariate COX regression analysis was performed to evaluate the associations of SLC4A4 and clinical features with OS by calculating hazard ratios (HRs) and 95% CIs. All statistical analyses were performed using R software (version 3.6.2), with p<0.05 considered statistically significant.

COMPREHENSIVE ANALYSIS:

Gene Expression Profiling Interactive Analysis (GEPIA) is a web server for transcriptome analysis [26] that includes 8,587 normal samples and 9,736 tumor samples from the GTEx and TCGA databases, with output provided by a standard processing pipeline. The association between SLC4A4 expression and survival in patients with colon adenocarcinoma was assessed by GEPIA. In addition, a boxplot was drawn to show where the disease state can be considered a variable to show differences in SLC4A4 expression among different tissue. Results were verified in two cohorts. The expression profiles in normal and tumor tissues were compared in GSE32323, and the association between SLC4A4 expression and OS in GSE41258.

TIMER includes genes with mutation frequency >10% in colon adenocarcinoma. To clarify the correlation between SLC4A4 expression and genes with high mutation frequency, the relationships between the expression of SLC4A4 and the 20 genes with highest mutation frequency were calculated. These 20 genes included APC, TP53, TTN, PIK3CA, KRAS, PTEN, ATM, MUC4, SMAD4, SYNE1, MUC16, FBXW7, BRAF, KIT, CSMD3, PTCH1, NEFH, NF1, RB1, and FAT4. Because cytokines are important components of the tumor microenvironment, the correlations between SLC4A4 expression and the expression of cytokine encoding genes, including IFNG, IL6, IL10, TGFB1, TGFB2, and TGFB3, were calculated. Pearson’s correlation coefficients for all of these associatins were calculated using the R package, Hmisc.

The human protein atlas (HPA; www.proteinatlas.org) is a comprehensive database containing tissue, cell, pathology, brain, blood and metabolic atlases. The expression of SLC4A4 in normal and tumor tissue was assessed. The HPA contains 24028 antibodies to analyze 16975 protein-encoding genes in 44 major tissues and organs in the human body. The 13 million immunohistochemistry images obtained were each annotated by pathologists as being from normal or tumor tissue.

Results

EXPRESSION PROFILE OF SLC4A4 AND ITS CORRELATION WITH OTHER GENES:

Differences in SLC4A4 expression were evaluated using the sequencing data obtained from TCGA. A comparison of 453 colon adenocarcinoma and 42 normal colon tissue samples showed that SLC4A4 expression was lower in the adenocarcinomas (Figure 1A), and that the density distribution of SLC4A4 expression differed in tumor tissue and normal tissue (Figure 1B). Pan-cancer analysis revealed that SLC4A4 expression was lower in 13 tumor types, including bladder urothelial carcinoma, breast invasive carcinoma, colon adenocarcinoma, esophageal carcinoma, head & neck squamous cell carcinoma, kidney chromophobe cell carcinoma, renal cell clear cell carcinoma, renal papillary cell carcinoma, hepatocellular carcinoma, lung adenocarcinoma, lung squamous cell carcinoma, rectal adenocarcinoma and thyroid carcinoma, but was higher in cholangiocarcinoma and prostate adenocarcinoma. The correlations between SLC4A4 expression and the expression of genes encoding cytokines are presented in Figure 1C. The Pearson correlation coefficients of IFNG, IL6, IL10, TGFB1, TGFB2, and TGFB3 with SLC4A4 were 0.180 (P=2.01E-04), 0.117 (P=1.63E-02), 0.272 (P=1.21E-08), 0.219 (P=5.54E-06), 0.211 (P=1.19E-05) and 0.214 (P=8.77E-06), respectively. Figure 1D shows the correlations between SLC4A4 expression and the expression of genes with high mutation frequency in colon adenocarcinoma. The Pearson correlation coefficients of MUC4 and SMAD4 with SLC4A4 were 0.488 (P=7.95E-29) and 0.360 (P=1.81E-15), respectively.

RELATIONSHIP BETWEEN SLC4A4 EXPRESSION AND CLINICAL CHARACTERISTICS:

Figure 2A is a scatter plot showing the expression of SLC4A4, OS of patients with colon adenocarcinoma, and the average expression of SLC4A4 in patients who were dead or censored at the time of the last follow-up. The average expression of SLC4A4 was lower in patients who died than in the censored group, with the density peaks of these two groups being significantly different. Figure 2B shows that the expression of SLC4A4 was highly sensitive and specific in distinguishing colon adenocarcinoma from normal colon tissue (AUC=0.721). Univariate COX regression analysis showed that increased SLC4A4 expression was a protective factor (HR <1), being associated with increased OS. In addition, tumor status, lymph node disease, metastasis, and tumor stage wer associated with OS (Figure 3A). Logistic regression analysis comparing differences in SLC4A4 expression in patients with different clinical conditions found that the expression of SLC4A4 was significantly lower in groups with positive lymph node biopsy or distant metastasis than in groups without these conditions (Figure 3B). In addition, the average expression of SLC4A4 in stage IV colon adenocarcinoma was lower than in carcinoma tissue at other stages (Figure 3C). SLC4A4 expression was lower in stage IV tissue samples than in stage I, II, and III samples, with log2 FCs of −0.418 (P=5.30E-03), −0.585 (P=7.46E-06), and −0.348 (P=1.21E-02), respectively.

RESULTS OF FUNCTIONAL ENRICHMENT ANALYSIS:

The influence of SLC4A4 expression on biological processes was analyzed by GO and KEGG pathways (Table 1). GSEA, a knowledge-based approach for interpreting genome-wide expression profiles, revealed a significantly different (FDR <0.05) functional pathway induced by decreased SLC4A4 expression [27]. Five GO entries were found to be positively correlated with decreased SLC4A4 expression, including nuclear-transcribed mRNA catabolic process, positive regulation of metaphase/anaphase transition of the cell cycle, protein localization to nuclear bodies, positive regulation of establishment of protein localization to telomeres and protein localization to Cajal bodies (Figure 4A). The five most significant entries that negatively correlated with decreased SLC4A4 expression were chemical synaptic transmission, trans-synaptic signaling, action potential, regulation of postsynaptic membrane potential and regulation of calcium ion-dependent exocytosis. In the KEGG pathway analysis, the five most positively correlated pathways with decreased SLC4A4 expression were steroid biosynthesis, mismatch repair, base excision repair, DNA replication and proteasome, whereas the five most negatively correlated pathways are cholinergic synapse, cytokine-cytokine receptor interaction, cAMP signaling pathway, calcium signaling pathway and neuroactive ligand-receptor interaction. The chronological sequence of these biological processes and tumorigenesis cannot be determined, so causality should be carefully evaluated. Some results of GSEA have been demonstrated to be associated with colon adenocarcinoma, whereas others require further research.

DNA mismatch repair is one of the main causes of colon adenocarcinoma [28]. Checkpoint inhibitors such as pembrolizumab and nivolumab after progression on irinotecan- or oxaliplatin-based therapies may benefit patients who exhibit DNA mismatch repair [29]. This imbalance in base excision repair may determine treatment response of colon adenocarcinomas and may be prognostic for survival [30,31].

CALCULATION AND EVALUATION OF TIICS:

Because tumor-infiltrating lymphocytes may be independently predictive of OS in cancer patients [32–34], the association between SLC4A4 expression and immune cell components in the tumor microenvironment were analyzed by TIMER. Immune cell populations that positively correlated with SLC4A4 expression in colon adenocarcinomas included B cells (P=9.66E-08), CD8+ T cells (P=3.53E-06), neutrophils (P=4.48E-06) and dendritic cells (P=2.89E-07). These results indicate that SLC4A4 performs an important function in the immune microenvironment of colon adenocarcinoma (Figure 4B).

To further analyze the association between SLC4A4 expression and immune cell components, all colon adenocarcinoma samples were divided into two groups based on median SLC4A4 expression. Samples with higher than median SLC4A4 concentrations were assigned to the high expression group and those with lower than median SLC4A4 concentrations to the low expression group. The association between SLC4A4 expression and TIICs was further clarified using an immune infiltration evaluation tool, CIBERSORT, to analyze the immune cell composition in each group, followed by statistical analyses using the Mann-Whitney U test (Figure 5A). CIBERSORT uses the LM22 gene signature to annotate the transcriptome and determine the composition of immune cells in samples. Because SLC4A4 expression is low in tumor samples, this group was defined as the high-risk group. The numbers of resting dendritic cells, activated mast cells, monocytes and resting memory CD4+ T cells were significantly lower in the high-risk group. In contrast, the numbers of some eosinophils, macrophages M0, resting mast cells and resting NK cells were significantly higher in the high-risk group. The infiltration levels of resting dendritic cells, eosinophil, monocytes and resting NK cells were too low (Figure 5B). Because activated and resting mast cells are two states of the same cells, the comprehensive changes are consistent. Calculations of the correlations among 22 types of immune cells showed weak or moderate correlations between the proportions of components (Figure 5C).

VERIFICATION OF SLC4A4 EXPRESSION CHARACTERISTICS:

GEPIA includes high-throughput sequencing data from tumor and normal tissue in the TCGA and GTEx databases. Because the GTEx database contains only normal tissues, GEPIA survival analysis was based on tumor samples from TCGA. Levels of SLC4A4 expression were divided into quartiles, with samples in the highest and lowest quartiles defined as the high and low expression groups, respectively. Kaplan-Meier analysis showed that the survival ratio was lower in the high expression group (Figure 6A), and a box plot showed that SLC4A4 expression differed in tumor and normal tissue samples (Figure 6B). Similar results were observed in the GSE32323 cohort (Figure 6C) and similar survival probability in the GSE41258 cohort (Figure 6D), suggesting that SLC4A4 is a protective factor in colon adenocarcinoma. Immunohistochemical analysis of the HPA database showed that SLC4A4 expression was lower in glandular cells of colon adenocarcinoma tissue than in normal colon gland cells (Figure 6E), confirming that SLC4A4 is a protective factor and its expression lower in colon adenocarcinoma than in normal colon tissue.

Discussion

Colon cancer is a universal public health problem and a major reason cause of cancer-related mortality [35]. Worldwide estimates indicate d that more than 2,200,000 individuals will be newly diagnosed with colon cancer in 2030, and that about 1,100,000 will die of this disease [36]. Early diagnosis requires the identification of appropriate molecular markers. The serum concentration of TRIM72 has been reported lower in patients with than without colon cancer, with TRIM72 being a better diagnostic indicator of colorectal cancer than CEA or CA199 alone [32]. When these three markers are used in combination, they are highly diagnostic, with an area under the curve of 0.928 [37].

The epithelial-mesenchymal transition is a marker of heterogeneity and pathological progression of colon cancer. CYB5R1 was found to be highly expressed in tumor cells undergoing the epithelial-mesenchymal transition, with CYB5R1 expression necessary for the invasive phenotype of tumor cells [38]. These studies indicate the significant clinical significance of identifying phenotypes through molecular features. To our knowledge, the current study is the first systematic analysis of the association between SLC4A4 expression and various pathological features of colon adenocarcinoma, including phenotype, tumor immune microenvironment, high mutation frequency gene profile, possible signaling pathways and prognosis.

We found that the expression of SLC4A4 was lower in colon adenocarcinoma than in normal colon tissue, a finding verified in the GEPIA and HPA databases. Because the expression of SLC4A4 is also significantly reduced in 12 other tumor types, SLC4A4 likely inhibits tumorigenesis through a ubiquitous mechanism. When SLC4A4 expression decreases, its protective effect is weakened, which indirectly promotes tumorigenesis. Logistic regression analysis showed that patients with low SLC4A4 expression were at higher risk for lymph node invasion and distant metastasis. SLC4A4 expression was lower in stage IV tumors than in other stages, suggesting that SLC4A4 may be involved in the occurrence of advanced events in colon tumorigenesis. SLC4A4 expression showed moderate correlations with MUC4 and SMAD4 mutations. Genetic variations in MUC4 affect the survival of patients with colon adenocarcinoma and are often accompanied by microsatellite instability [39,40]. Clonal loss of SMAD4 expression accumulates with disease progression in adenomas and promotes inflammation-driven colon cancer [41,42]. GSEA showed that SLC4A4A expression was associated with steroid biosynthesis, mismatch repair, base excision repair, DNA replication and the proteasome. Steroid biosynthesis is important in patients with hormone-dependent tumors, as anti-tumor drugs such as CYP17A1 inhibitors designed to suppress steroid biosynthesis are used in the treatment of prostate cancer [43,44]. Few studies to date, however, have analyzed whether steroid biosynthesis promotes or suppresses the development of colon adenocarcinoma, suggesting the need for additional research. Systemic mismatch repair deficiencies are usually manifested as microsatellite instability, resulting in the unique clinicopathological characteristics of individual colon adenocarcinomas [45]. Microsatellite instability can also predict tumor resistance to certain chemotherapeutic agents, such as immune checkpoint inhibitors [46]. Base excision repair is critical for removing endogenously damaged DNA bases, with failure to remove these bases leading to chromosomal instability and tumor development [47].

Tumor cells have three unique attributes, genomic instability, evasion of apoptosis, and continued proliferation. The DNA replication stress model proposes that genes that drive cell proliferation can induce DNA replication, which leads to genome instability and the selection of cells that can escape apoptosis [48]. The proteasome functions significantly in the occurrence, growth and invasion of tumor cells, with the proteasome considered a target for anti-cancer drugs [49]. Several newly discovered natural proteasome inhibitors and drugs found to function similarly have been used in clinical treatment [50,51]. The GSEA indicates the importance of the mechanism by which SLC4A4 participates in and affects mismatch repair, base excision repair, DNA replication and the proteasome.

The immunosuppressive mechanism in the tumor microenvironment is a protective barrier for tumor cells. Treatment of tumor tissue with various immunomodulators makes the tumor itself a source of antigens, which activate the dormant immune system, a promising treatment method known as in situ vaccination [52]. Because changes in immune components in the tumor microenvironment are of great significance, TIMER and CIBERSORT were applied to explore the association between SLC4A4 expression and immune cells in colon adenocarcinomas. The numbers of B cells, dendritic cells, CD8+ T cells, and neutrophils were found to positively correlate with SCL4A4 expression, indicating that the proportion of these cells was lower in tumor than in normal tissue. B cells in the tumor microenvironment can not only trigger humoral immunity and secrete inflammatory factors, but also recognize, process and present antigens. Moreover, B cells maintain positive interactions with other immune cells, including T cells and all immune cells expressing Fc receptors [53]. The prognostic significance of CD3+ T cells and CD8+ T cells is usually consistent with that of tumor infiltrating B cells, with these T cell populations also having a stronger prognostic significance when coexisting with tumor infiltrating B cells [54]. Dendritic cells activate CD8+ T cells through cross-presentation of exogenous antigens, a step considered key to activate CD8+ T cells against tumors [55,56]. These findings suggest that reductions in the numbers of dendritic cells and CD8+ T cells lead to a poor prognosis in patients with colon adenocarcinoma. Based on the CIBERSORT algorithm, we compared immune cells in tumors expressing high and low levels of SLC4A4, find that eight of the 22 types of immune cells differed significantly in these two groups, suggesting that SLC4A4 widely affects TIICs.

In summary, this study systematically analyzed the significance of SLC4A4 expression in colon adenocarcinoma, especially its association with tumor infiltrating immune cells. SLC4A4 may be a marker for the diagnosis of colon adenocarcinoma and a target for treatment, as well as making possible combinations of targeted therapy and immunotherapy in colon cancer.

Conclusions

SLC4A4 extensively affected tumor-infiltrating immune cells in colon adenocarcinoma tissues and was associated with tumor characteristics. SLC4A4 may be a new biomarker for the diagnosis and prognosis of patients with colon cancer.