Exploration of prognostic biomarkers and therapeutic targets in the microenvironment of bladder cancer based on CXC chemokines

1.   Introduction

Bladder cancer (BLCA) is characterized by a high rate of morbidity and mortality, and is diagnosed in approximately 430,000 people per annum, of whom around 165,000 die of this cancer each year ^[1]. More than 90% of BLCA is considered to be urinary urothelial cancer according to its pathological histological subtype and smoking as a risk factor caused its formation ^[2]. Although the results of a prior clinical trial showed 5- and 10-year survival rates of approximately 50 and 36%, respectively, for muscle-invasive BLCA using neoadjuvant chemotherapy ^[3], other studies showed 5-year survival rates ranging from 15 to 20% and high recurrence rates ^[4,5]. Therefore, for patients with high-risk non-metastatic BLCA, the recommended therapeutic strategy is radical cystectomy with lymph node dissection because of the poor efficacy of drug treatment, including maintenance bacille Calmette-Guérin (BCG) ^[6,7]. To date, tumor stage and lymph node (LN) status have been identified as the most important prognostic factors after radical cystectomy ^[8], whereas in clinical practice, the value of other prognostic risk factors, including neutrophil-to-lymphocyte ratio, remains controversial ^[9,10]. Recently, molecular subtypes of BLCA derived from The Cancer Genome Atlas (TCGA) database attracted our attention. However, limited clinical value have restricted their application ^[11]. To date, although promising predictive biomarkers have been identified, they are not used routinely in clinical practice because of their finite values, and decision-making concerning therapy cannot be based only on these molecular markers, for example the tumor mutation burden, ribonucleic acid subtypes, and neutrophil-to-lymphocyte ratio ^[12]. Therefore, it is urgent to identify more therapeutic targets and prognostic biomarkers for use in patients with BLCA.

CXC chemokines represent the main component of approximately 50 chemokine family members, and their C-termini harbor a heparin-binding domain, which plays vital roles in the regulation of tumor-associated angiogenesis, tumor growth, and metastatic potential ^[13,14], in addition to mediating the migration of diversified leukocytes in non-tumor/tumor microenvironment (TME) ^[15,16,17]. The key role of CXC chemokines has been established as linking tumor development and metastasis; therefore, CXC chemokines and their receptors have been identified as therapeutic target in several cancers. Mounting experimental evidence shows that monoclonal antibodies against CXC chemokines or their receptors could inhibit tumor growth and/or metastasis, including breast cancer ^[18], non-Hodgkin lymphoma ^[19] brain tumors ^[20], pancreatic tumors, and bladder tumors ^[21]. Therefore, CXC chemokines are likely to become therapeutic targets or prognostic factors in many tumors including BLCA.

Previous studies have revealed the functional role of several CXC chemokines in BLCA ^{[22,23,24,25,26]}. However, the identification of optimal prognostic biomarkers and therapeutic targets originating from CXC chemokines in BLCA remains an urgent and unresolved problem. Therefore, with the aim of fully comprehending the function of CXC chemokines in BLCA, the present study investigated the expression of CXC chemokines and their potential as prognostic biomarkers and therapeutic targets using exhaustive bioinformatic analyses based on several established databases. In addition, we provided evidence that will allow clinicians to choose appropriate therapeutic targets and accurately predict clinical outcome in patients with BLCA.

2.   Materials and methods

2.1.   Oncomine

Oncomine (www.oncomine.org) integrates RNA and DNA-sequencing data from sources such as GEO, TCGA, and the published literature ^[27]. Screening of mRNA data allowed us to investigate CXC chemokine expression in BLCA. We set the significance thresholds of a p-value less than 1E-4, two-fold change, and overexpression or underexpression gene rank in top the 10%. The differential expression of CXC chemokine mRNAs in BLCA was analyzed using Student's t test.

2.2.   UALCAN

UALCAN (http://ualcan.path.uab.edu/analysis.html) is an effective online analysis and mining website of tumor data, mainly based on the TCGA database, and can be used to perform biomarker identification, expression profile analysis, and survival analysis ^[28]. In this study, CXC chemokine mRNA expression data was obtained using the "Expression Analysis" module of the "BLCA" dataset. The differential expression of CXC chemokine mRNAs in BLCA was analyzed using Student's t test.

2.3.   GEPIA

GEPIA (http://gepia.cancer-pku.cn/index.html) is an analysis website of RNA sequencing data based on TCGA and GTEx databases, which was developed by Tang et al. ^[29]. We investigated the relationship between pathological stage and CXC chemokines with differential mRNA expression and their prognostic significance by screening Oncomine and UALCAN databases using the "Single Gene Analysis" module of GEPIA in BLCA.

2.4.   GEO databases

GEO (Gene Expression Omnibus) databases collect and organize various expression chip data, and other chips, such as methylation chips, long noncoding RNA (lncRNA) chips, miRNA chips, chromosome number variation (CNV) chips, and even high-throughput sequencing data. CXC chemokine mRNA differential expression in BLCA was analyzed using Student's t test and the corresponding prognostic analysis was performed using the statistical software GraphPad Prism 8.0.1 (GraphPad Inc., La Jolla, CA, USA).

2.5.   cBioPortal

The data in cBioPortal (www.cbioportal.org) is derived from the International Cancer Genome Consortium (ICGC), TCGA, GEO, and other databases, which includes DNA methylation data, limited clinical data mRNA and microRNA expression data, non-synonymous mutations, protein level and phosphoprotein level (reverse-phase protein array (RPPA)) data, and DNA copy number data ^[30]. Based on TGCA data, cBioPortal provided genetic alterations and co-expression data for CXC chemokines. Bladder Urothelial Carcinoma (TCGA, Firehose Legacy) data from 413 cases were analyzed. Protein expression z scores (RPPA) and mRNA expression z scores (RNA Seq V2 (RNA-Seq by expectation-maximization (RSEM)) were entered using a z score threshold of ± 2.0.

2.6.   GeneMANIA

The GeneMANIA plug-in of Cytoscape software is used to study protein-protein interactions ^[31].

2.7.   STRING

STRING (https://string-db.org/) is an important software used to analyze protein-protein interactions ^[32]. To investigate the interactions of CXC chemokines, we constructed a network analysis of protein-protein interactions using STRING.

2.8.   DAVID 6.8

DAVID 6.8 (https://david.ncifcrf.gov/home.jsp) is used mainly for functional and pathway enrichment analysis for differentially expressed genes ^[33]. Differentially expressed CXC chemokines were analyzed using Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis and Gene Ontology (GO) enrichment analysis. The top 50 mutated genes of Bladder Urothelial Carcinoma (TCGA, Firehose Legacy; p < 0.05) were analyzed using DAVID 6.8. The "ggplot2" package in the R project visualized the data.

2.9.   Metascape

Metascape (http://metascape.org) is a powerful gene function annotation analysis tool, which can perform annotation enrichment for large quantities of genes or proteins, and construct protein-protein interaction networks ^[34]. Our study utilized the "Express Analysis" module to validate the enriched CXC chemokines and the top 50 mutated genes of Bladder Urothelial Carcinoma (TCGA, Firehose Legacy).

2.10.   TRRUST (version 2)

The TRRUST (version 2) database (https://www.grnpedia.org/trrust/) records transcription factor regulatory relationships, including the target genes corresponding to transcription factors, and the regulatory relationships among transcription factors ^[35].

2.11.   TIMER2.0

TIMER2.0 (http://timer.cistrome.org/) is a tumor immunity related database. "Gene module" refers to analyzing the association between CXC chemokine expression and immune cell types. "Clinical module" assesses the relationship among clinical outcome and tumor infiltrating immune cells and CXC chemokines ^[36].

2.12.   LinkedOmics

LinkedOmics (http://www.linkedomics.org/) was used to analyze tumor multi-omics data for 32 TCGA tumor types ^[37]. "LinkInterpreter module" was used to obtain biological insights into transcription factor, miRNA, and kinase targets of enriched CXC chemokines. A minimum number of genes of 3 and a simulation of 500 were performed in Gene Set Enrichment Analysis (GSEA) for the BLCA dataset.

2.13.   Statistical analysis

Heatmap of the correlations of differentially expressed CXC chemokines and bubble plot of enrichment analysis of the 50 top mutated genes and differentially expressed CXC chemokines in bladder cancer was conducted using R version 3.6.3. Visualization of gene-gene interaction network of the 50 top mutated genes and the differentially expressed CXC chemokines was performed with Cytoscape software. CXC chemokine mRNA differential expression in BLCA was analyzed using Student's t test and the corresponding survival analysis was performed using the statistical software GraphPad Prism 8.0.1 (GraphPad Inc., La Jolla, CA, USA). The plots of Figures 1 and 2 were achieved from online database UALCAN (http://ualcan.path.uab.edu/analysis.html).

3.   Results

3.1.   CXC chemokine expression status in patients with BLCA

ONCOMINE database analysis allowed us to analyze the mRNA expression levels of 16 CXC chemokines (CXCL1 to CXCL14, CXCL16 and CXCL17) in BLCA and normal bladder mucosa tissues. Table 1 and Figure 3 show the results. ONCOMINE data demonstrated the mRNA expression levels of CXCL13, CXCL10, CXCL9, and CXCL6 increased significantly and those of CXCL2, CXCL3, and CXCL12 decreased significantly in BLCA tissues compared with those in normal bladder mucosa tissues. This agreed with the results of Sanchez-Carbayo et al., who reported significantly increased mRNA levels of CXCL13, CXCL10, CXCL9, and CXCL6 in infiltrating bladder urothelial carcinoma and obviously decreased mRNA levels of CXCL2, CXCL3, and CXCL12 in superficial bladder cancer ^[38]. Lee et al. observed lower CXCL12 expression was in superficial bladder cancer (fold change = -4.495 and p = 2.52 ^e-21) and infiltrating bladder urothelial carcinoma (fold change = -3.371 and p = 1.41 ^e-12) compared with that in normal bladder mucosa ^[39].

Next, we used UALCAN to evaluate the mRNA expression levels of CXC chemokines in BLCA and normal bladder mucosa tissues, which showed markedly increased mRNA expression of CXCL1 (p = 8.50^e-03), CXCL5 (p = 1.12^e-06), CXCL7 (p = 1.42^e-02), CXCL9 (p = 4.06^e-08), CXCL10 (p = 3.91^e-09), CXCL11 (p = 6.29^e-08), CXCL16 (p = 1.44^e-03), and CXCL17 (p = 3.56^e-02) in BLCA tissues, and significantly decreased mRNA levels of CXCL2 (p = 2.30^e-02) and CXCL12 (p = 2.82^e-04) (Figure 1). We also revealed that CXC17 had the highest mRNA expression level in BLCA tissues based on the TCGA database (Figure 2). To investigate the role of the CXC chemokines in BLCA more comprehensively, including tumorigenesis and clinical prognosis, we assessed all the differentially expressed CXC chemokines that were acquired by the Oncomine and TCGA databases including CXCL1, CXCL2, CXCL3, CXCL5, CXCL6, CXCL7, CXCL9, CXCL10, CXCL11, CXCL12, CXCL13, CXCL16, and CXCL17.

3.2.   CXC prognostic significance in patients with BLCA

Next, we evaluated the relationship between the pathological stage of patients with BLCA and the differentially expressed CXC chemokines. We observed a significant association between CXCL2 (p = 0.0127) and CXCL12 (p = 8.68^e-07) expression and pathological stage (Figure 4). CXCL2 and CXCL12 mRNA expression levels were higher in more advanced pathological stages, indicating that they might associated with BLCA tumorigenesis. Meanwhile, GEPIA was used to assess the impact of differentially expressed CXC chemokines on clinical outcome. Patients with BLCA with high CXCL2 (p = 0.17), CXCL16 (p = 0.14), and CXCL13 (p = 0.12) mRNA expression had a statistical tendency towards longer disease-free survival (Figure 5). In addition, high mRNA levels of CXCL1 (p = 0.17) and CXCL12 (p = 0.0075) had a statistical tendency towards shorter overall survival (Figure 6).

GEO database analysis demonstrated that CXCL3 mRNA expression (p = 0.034) was higher in BLCA tissues than in normal tissues. CXCL12 (p < 0.001), CXCL9 (p = 0.002), and CXCL10 (p = 0.015) mRNA expression levels were higher in normal tissues than in BLCA tissues (Figure 7). Survival analysis showed that high levels of CXCL1, CXCL6, CXCL10, CXCL11, and CXCL13 mRNA expression were associated significantly with the poor overall survival (all p < 0.05). High levels of CXCL3 mRNA had the statistic tendency towards predicting shorter overall survival (p = 0.053) (Figure 8).

3.3.   CXC chemokine interaction analyses, neighbor gene network, co-expression, and genetic alteration in patients with BLCA

We analyzed comprehensively the molecular characteristics of the differentially expressed CXC chemokines. First, we used TCGA datasets to analyze their genetic alterations, which showed the rates of genetic alteration of CXCL1, CXCL2, CXCL3, CXCL5, CXCL6, CXCL7, CXCL9, CXCL10, CXCL11, CXCL12, CXCL13, CXCL16, and CXCL17 were 6, 5, 4, 5, 5, 4, 5, 6, 5, 2.7, 3, 5, and 11%, respectively in the BLCA samples queried (Figure 9A). In the BLCA samples, the most common alterations were amplification and high mRNA expression. Differentially expressed CXC chemokine co-expression analysis identified strong correlations among CXCL9, CXCL10, and CXCL11 expression; moderate to strong correlations among CXCL1, CXCL2, CXCL3, CXCL5, and CXCL6; and weak to moderate correlations among CXCL7, CXCL12, CXCL13, CXCL16, and CXCL17 expression (Figure 9B). Secondly, PPI network analysis with STRING obtained 13 nodes and 71 edges in the analysis of the interactions among the differentially expressed CXC chemokines (Figure 9D). Enrichment pathways identified in STRING for these CXC chemokines included the TNF signaling pathway, the IL-17 signaling pathway, cytokine-cytokine receptor interactions, and the chemokine signaling pathway. The results of analysis at GeneMANIA revealed that their functions were related to chemokine receptor binding, cell chemotaxis, and chemokine activity (Figure 9E).

In addition, the top 50 mutated genes of 200 most mutated cases with BLCA were identified using cBioPortal (Figure 9C), representing mutations in which lead to the formation of a tumor-associated antigen that is associated with the TME and is closely associated with the CXC chemokines. Next, we analyzed the associations of these mutated genes and the CXC chemokines. These top 50 mutated genes have close interactions with the differentially expressed CXC chemokines in BLCA (Figure 9F).

3.4.   Functional enrichment analysis of CXC chemokines in patients with BLCA

The functions of differentially expressed CXC chemokines and the top 50 mutated genes in BLCA were analyzed by DAVID6.8 and Metascape. The top 10 GO enriched terms in the biological process (BP), cellular component (CC), and molecular function (MF) categories were analyzed using DAVID 6.8. Chemokine-mediated signaling pathway, positive regulation of leukocyte chemotaxis, response to lipopolysaccharide, regulation of cell proliferation, inflammatory response, immune response, positive regulation of cAMP metabolic process, and cell-cell signaling were related to the development and progression of BLCA (Figure 10A). The top 10 GO terms in the CC category included, Z disc, sarcolemma, extracellular space, sarcoplasmic reticulum membrane, sarcoplasmic reticulum, extracellular region, T-tubule, MLL3/4 complex, costamere, and histone methyltransferase complex (Figure 10B). Chemokine activity and CXCR chemokine receptor binding were the two most significantly enriched GO terms in the MF category (Figure 10C). The top 10 KEGG pathways were the chemokine signaling pathway, cytokine-cytokine receptor interaction, cell cycle, the TNF signaling pathway, and microRNAs in cancer, HIF-1 signaling pathway, pathways in cancer, viral carcinogenesis, central carbon metabolism in cancer, which were intimately associated with BLCA tumorigenesis (Figure 10D).

Metascape analysis showed that the functions of the differentially expressed CXC chemokines and top 50 mutated genes in BLCA were predominantly enriched in chemokine receptor bind chemokines, T cell chemotaxis, and microRNAs in cancer (Figure 11A, B). Next, we extracted protein-protein interaction (PPI) and MCODE components to explore the correlations between BLCA and CXC chemokines in detail; the list of identified genes is shown in Figure 11C, D. Three most significant MCODE components were obtained and showed that the biological functions were involved in chemokine receptor binding to chemokines, chemokine-mediated signaling pathways, chemokine signaling pathways, activation of homeobox (HOX) genes during differentiation, and histone modification (Figure 11E).

3.5.   MiRNA targets, kinase targets, and transcription factor targets of CXC chemokines in patients with BLCA

Next, we investigated the transcription factor and kinase targets of the differentially expressed chemokines by analyses in the TRRUST and LinkedOmics databases. TRRUST analysis of CXCL1, CXCL2, CXCL5, CXCL7, CXCL10, and CXCL12 identified three transcription factors: Sp1 transcription factor (SP1), nuclear factor kappa B subunit 1 (NFKB1), and RELA proto-oncogene, NF-KB subunit (RELA) that were possibly related CXC chemokine regulation (Table 2). Critically, NFKB1 and RELA mediated the expression of CXCL1, CXCL2, CXCL5, CXCL10, and CXCL12, and SP1 mediated the expression of CXCL1 and CXCL5. Next, LinkedOmics database analysis was used to identify the top two kinase targets of these CXC chemokines (Table 3). LYN proto-oncogene, src family tyrosine kinase (LYN) was identified as the kinase target of CXCL1, CXCL2, CXCL3, CXCL5, CXCL6, CXCL9, CXCL10, CXCL11, CXCL12, and CXCL13 in the kinase-target network, and LCK proto-oncogene, src family tyrosine kinase (LCK) was identified as the kinase target of CXCL2, CXCL3, CXCL5, CXCL6, CXCL9, CXCL10, CXCL11, CXCL13, and CXCL16. Inhibitor of nuclear factor kappa B kinase subunit beta (IKBKB), protein kinase CGMP-dependent 1 (PRKG1), and Janus kinase 2 (JAK2) were identified as key kinase targets of CXCL1, CXCL12, and CXCL16, respectively, in the kinase-target network. Moreover, protein kinase N1 (PKN1), adrenergic, beta, receptor kinase 1 (ADRBK1) and checkpoint kinase 1 (CHEK1), polo like kinase 1 (PLK1) were identified as the top two target kinase for CXCL7 and CXCL17, respectively. Finally, analysis at the LinkedOmics database allowed us to identify miRNAs that targeted the mRNAs encoding CXC chemokines (Table 4). The miRNA miR200A interacted mainly with the CXCL1 and CXCL5 miRNA-target network. Chemokine mRNAs CXCL1 and CXCL6 might be targeted by miR-30E-5P, miR-30B, miR-30D, miR-30C, and miR-30A-5P. CXCL3 expression might be modulated by miR-409-3P and miR-526B. CXCL5 and CXCL6 might be targeted by miR-141 and miR-374, respectively. CXCL7 might be targeted by miR-154, miR-487, miR-525, and miR-524s. Finally, CXCL17 might be targeted by miR-490.

3.6.   Immune cell infiltration of CXC chemokines in patients with BLCA Patients

CXC chemokines influence immune cell infiltration and inflammatory responses in the TME, thus possibly affecting the clinical therapy and outcome in patients with BLCA. Using the TIMER database, we investigated comprehensively the correlation between infiltration of various immune cells and the differentially expressed CXC chemokines. CXCL1 expression correlated positively with the infiltration of dendritic cells (Cor = 0.427, p = 1.43^e-17, neutrophils (Cor = 0.444, p = 5.46^e-19), CD4+ T cells (Cor = 0.112, p = 3.26^e-02), and CD8+ T cells (Cor = 0.303, p = 3.39^e-09); Figure 12A). Similarly, CXCL2 expression correlated positively with the infiltration of dendritic cells (Cor = 0.345, p = 1.22^e-11, neutrophils (Cor = 0.355, p = 3.05^e-12), CD4+ T cells (Cor = 0.128, p = 1.46^e-02), and CD8+ T cells (Cor = 0.255, p = 7.61^e-07; Figure 12B). CXCL3 expression correlated negatively with the infiltration of dendritic cells (Cor = 0.4, p = 1.93^e-15, neutrophils (Cor = 0.424, p = 2.86^e-17), CD4+ T cells (Cor = 0.108, p = 3.89^e-02), CD8+ T cells (Cor = 0.302, p = 3.55^e-09), and B cells (Cor = -0.103, p = 4.88^e-02); Figure 12C). CXCL5 expression correlated positively with the infiltration of dendritic cells (Cor = 0.43, p = 7.83^e-18), neutrophils (Cor = 0.385, p = 2.70^e-14), CD4+ T cells (Cor = 0.116, p = 2.63^e-02), CD8+ T cells (Cor = 0.288, p = 2.10^e-08; Figure 12D). CXCL6 expression correlated positively with the infiltration of CD8+ T cells (Cor = 0.219, p = 2.42^e-05), neutrophils (Cor = 0.339, p = 3.50^e-11), and dendritic cells (Cor = 0.43, p = 7.83^e-18; Figure 12E). A similar correlation tendency was found between CXCL7 expression and dendritic cells (Cor = 0.239, p = 3.74^e-06), neutrophils (Cor = 0.174, p = 8.87^e-04), and CD8+ T cells (Cor = 0.164, p = 1.60^e-03; Figure 12F). CXCL9, CXCL10, and CXCL11 mRNA expression levels correlated positively with the infiltration of dendritic cells, neutrophils, CD4+ T cells, and CD8+ T cells (all p < 0.05; Figure 12G-I). CXCL12 expression correlated positively with the infiltration of dendritic cells (Cor = 0.156, p = 2.88^e-03), neutrophils (Cor = 0.208, p = 6.27^e-05), macrophages (Cor = 0.464, p = 7.25^e-21), CD4+ T cells (Cor = 0.146, p = 5.21^e-03), and CD8+ T cells (Cor = 0.134, p = 1.05^e-02; Figure 12J). With the exception of macrophages, CXCL13 expression correlated positively associated with the other immune cell types consisted of dendritic cells, CD4+ T cells, neutrophils, CD8+ T cells, and B cells (p < 0.05; Figure 12K). CXCL16 expression correlated positively with the infiltration of CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells (all p < 0.05; Figure 12L). CXCL17 mRNA expression correlated positively with B cells (Cor = 0.156, p = 2.83^e-03), and negatively with dendritic cells (Cor = -0.167, p = 1.35^e-03; Figure 12M). A Cox proportional hazard model was used identify factors associated with prognosis of patients with BLCA, which identified B cells (p = 0.032), macrophages (p = 0.001), and CXCL2 mRNA expression (p = 0.010). CXCL7 mRNA expression (p = 0.066) had significantly statistical tendency towards predicting the clinical outcome of BLCA patients (Table 5).

4.   Discussion

Preclinical data indicate that BLCA is one of the most immunogenic tumors ^[40], which is assumed to make it sensitive to favorable targeted therapies, not including PD-1/PDL-1 checkpoint inhibitors. Previous studies have unveiled a tumor microenvironment (TME)-related role of CXC chemokines ^[41], and functions in the immunotherapy of cancer ^[42]. Chemokines are involved in the process of cancer biology, including tumor proliferation and metastasis, and targeting of chemokines and chemokine receptors have been evaluated in several preclinical studies and clinical trials, such as in breast cancer ^[43], pancreatic cancer ^[44], renal tumor ^[45], and prostate cancer ^[46]. Nonetheless, in the prognosis and biology of BLCA, the significance of CXC chemokines are poorly understood.

The present study characterized the mRNA expression patterns CXC chemokines and investigated their pathological relationship with BLCA stages. Oncomine and TCGA database analysis identified 13 CXC chemokines that were differentially expressed in BLCA tissues compared with that in normal tissues, including upregulated expression (CXCL1, CXCL5, CXCL6, CXCL7, CXCL9, CXCL10, CXCL11, CXCL13, CXCL16, and CXCL17) and downregulated expression (CXCL2, CXCL3, and CXCL12). In addition, the TCGA analysis results showed that high expression of CXCL2 and CXCL12 was associated with advanced pathological stage, and analysis of the GEO and TCGA databases revealed that high levels of CXCL1, CXCL6, CXCL10, CXCL11, CXCL12, and CXCL13 mRNA expression were associated significantly with poor overall survival, indicating these differentially expressed CXC chemokines exert important roles in the proliferation and development of BLCA. This agreed partially with the results of previous studies, which reported that CXCL5 promotes chemoresistance in BLCA cells ^[47] and CXCL12 is associated with the genesis and development of BLCA ^[23].

We next investigated the differentially expressed CXC chemokines for their molecular characteristics in BLCA. Their most frequent genetic alterations were high mRNA expression and amplification, following by deep deletion, missense mutation, and truncating mutation, suggesting that genetic alterations might play a vital role in tumorigenesis and progression of BLCA. In particular, high mRNA expression and amplification of CXC chemokine genes could have the predominant impact on the BLCA processes. Previous studies reported that overexpression of several chemokines, including CXCL1, CXCL2, CXCL3, CXCL4, CXCL5, CXCL6, CXCL7, and CXCL8, was involved in breast cancer and melanoma metastasis, resulting from enhanced gene transcription ^[48] and copy number amplification of their encoding genes ^[49]. Moreover, mutants of CXCL8 had lower potency and efficacy in chemotaxis assays using neutrophils ^[50]. Together, our results and those of previous studies suggest the critical role of CXC chemokines in the tumorigenesis and cancer progression.

We next used GO and KEGG analysis to provide clues to the biological functions of these CXC chemokines. The results showed the CXC chemokines and top 50 mutated genes in BLCA are functionally associated with the TNF pathway, cell cycle, cytokine-cytokine receptor interactions, and the chemokine signaling pathway. Tumor metastasis, immune evasion, angiogenesis, and proliferation are reported to be associated critically with chemokine signaling pathways ^[51,52]. The TNF signaling pathway responds to cellular stress and inflammatory signals to activate cytokine cascades, which further affect the immune microenvironment. In addition, TNF interacts with different receptors to mediate cellular apoptosis or growth ^[53,54]. Moreover, the TNF signaling pathway contributes to tumor growth and metastasis in various cancers and acts as a therapeutic target ^{[55,56,57,58]}. These pieces of evidence showed that CXC chemokines might be developed as effective treatment targets for patients with BLCA in the future.

We explored the transcription factor targets of the differentially expressed CXC chemokines, their kinase targets, and miRNAs that targeted them. We observed that SP1, NFKB1, and RELA might regulate the expression of CXC chemokines. RELA regulates breast cancer progression via the CXCL13/CXCR5 axis ^[59]. In addition, via the IRF1-CDK4 signaling axis, RELA-induced Interferon-γ inhibits the proliferation of breast cancer ^[60]. RELA also regulates oncogene-induced senescence to decrease the proliferation in murine Kras-driven pancreatic cancer via the CXCL1/CXCR2 axis ^[61]. Additionally, RELA acts as an important transcription factor in the activation and stability of regulatory T cells ^[62] and regulates Th17 differentiation by miR-30a in experimental autoimmune encephalomyelitis ^[63]. NFKB1 is involved widely in carcinogenesis, including driving tumor progression or acting as a tumor-suppressor by phosphorylation and ubiquitination pathways ^[64]. Polymorphisms of the NFKB1 promoter are associated with a higher risk of bladder cancer ^[65], and aberrant NFKB1 activity leads to impaired M1 polarization of macrophages by repressing the expression of CXCL10 and IL-12, which forms an immunosuppressive microenvironment to induce the immune escape of tumors ^[64]. Our study provided evidence that an association is showed between CXC chemokines and the NF-κB signaling pathway during BLCA tumorigenesis and development. Our results also suggested Src family tyrosine kinases (LCK and LYN), PKN1, PRKG1, and CHEK1 might be the targets of the differential CXC chemokines, which play important roles in cell growth, division, migration, and survival signaling pathways ^{[66,67,68,69]}, and mediate the carcinogenesis of several tumors ^[70,71,72]. Therefore, these differentially expressed CXC chemokines might modulate BLCA development and progression dependent on the above-mentioned signaling pathways by regulating these kinases.

We found the miR-30 family is the main miRNA family targeting CXC chemokines such as CXCL1 and CXCL6 in BLCA, suggesting that this miRNA family plays crucial regulatory roles in BLCA carcinogenesis and progression by acting on mRNAs encoding CXC chemokines. This hypothesis was supported partially by previous reports that the miR-30 family is associated with the development of organs and the pathogenesis of diseases, including various cancers, such as breast cancer and colorectal cancer ^[73,74,75].

Chemokines function by mediating the migration of immune cells by interacting with their cognate ligands ^[41]. Immune cell infiltration, especially CD8+ T cells, could affect tumor progression, as well as predicting the sensitivity to immunotherapy and clinical prognosis ^[76,77]. The present study indicated that a distinct correlation is found between infiltrating immune cells (e.g., dendritic cells, neutrophils, macrophages, CD4+ T cells, CD8+ T cells, and B cells) and the differentially expressed CXC chemokines, suggesting that these chemokines acting as crucial components reflected the immune status in the TME and could provide evidence of drug targets in BLCA immunotherapy in the future.

5.   Conclusions

We analyzed the relationship between the CXC chemokines and the 50 top mutated genes in patients with BLCA and explored the CXC chemokines superfamily-based biomarkers associated with BLCA prognosis by public databases. The CXC chemokines superfamily-based biomarkers gene signature has prognostic value in BLCA patients and thus provides and guides the treatment of BLCA.

Acknowledgments

Xinke Zhang and Jiewei Chen designed this study. Xiaoqi Sun and Qunxi Chen analyzed the data and wrote the manuscript. Lihong Zhang analyzed the data. Xinke Zhang revised the manuscript. All authors approved the final version for submission.

This study was supported by the Medical Science Funding of Guangdong province (A2020139).

Conflict of interest

The authors declare that they have no conflict of interest.