Target genes in KIRC exhibit varied expressions between cancer and normal tissues.
The CXCL family genes in human kidney cancer and normal tissues were searched in GeneBank and 15 CXCL genes (CXCL1/2/3/4/5/6/8/9/10/11/12/13/14/16/17) were obtained. Oncomine and GEPIA analyses were performed to assess the expression levels of the 15 CXCL genes in KIRC vs. normal tissues. Eight CXCL genes (CXCL2/5/9/10/11/12/13/16) were differentially expressed (p <0.01, Fig. 1 A). In cancer tissues, seven genes (CXCL2/5/9/10/11/13/16) were up-regulated while CXCL12 was down-regulated. GEPIA was used to determine the correlation between gene expression levels and clinical stages. Four CXCL genes (CXCL5/9/10/11, Fig. 1 B) exhibited significant correlations with clinical stages (p <0.01). Oncomine was used to assess the expression levels of 8 CXCL genes in Pan-cancer (Fig. 1 C).
The PPI network and functional enrichments
The top 100 genes that were correlated with the 8 target CXCL genes were obtained and an intersection of all correlation genes established. Fifty-five genes were obtained. The STRING database was used to construct the PPI network and to perform GO/KEGG enrichment analysis (Fig. 2 A-B). From the PPI network, 13 genes were found to be closely associated with CXCL genes, such as, CCL5, CCL25, CCL27, and PF4. Enriched GO terms for the 8 CXCL genes were associated with: i. The chemokine-mediated signaling pathway (GO:0070098, GO:1990868, and GO:1990869) in biological processes; ii. The external side of the plasma membrane (GO:0009897), and secretory granule membrane (GO:0030667) in cellular components, and iii. Chemokine activities (GO:0008009) as well as chemokine receptor binding (GO:0042379) in molecular functions. KEGG analysis revealed that the CXCL genes were involved in viral protein interactions with cytokine and cytokine receptor (hsa04061) and Chemokine signaling pathway (hsa04062, Fig. 2 C).
TIMER correlation analysis revealed that several of the eight CXCL genes, including CXCL9 and CXCL10, CXCL9 and CXCL11, CXCL10 and CXCL11 were strongly correlated (Fig. 2 D-E).
Gene mutation analysis
CXCL2 and CXCL5 mutations and copy numbers were analyzed using the cBioPortal online tool. It was found that CXCL2 and CXCL5 were rarely mutated among KIRC patients (Fig. 3 A). The likelihood of two CXCL genes having the same copy number was low. Then, we searched for SNP mutations, somatic mutations, DNA methylation, and copy number variations for the two CXCL genes in UCSC Xena(Fig. 3 B). The findings were comparable to those obtained from the cBioPortal online tool.
Survival analysis
KM analysis revealed that only CXCL2 and CXCL5 had a survival value for KIRC patients (p<0.01, Fig. 4 A). ROC analysis showed that the one-year accuracy of CXCL2 and CXCL5 genes were 0.706 and 0.738, respectively (Fig. 4 B-C). To assess the association between clinical information and target genes, univariate and multivariate Cox proportional hazard regression analyses were performed (Table 1).
In the univariate Cox analysis, all seven indices were positively correlated (T stage, N stage, M stage, pathology stage, age and CXCL2 as well as CXCL5 levels; p<0.01). However, in multivariate Cox analysis, two indices (M stage and age) were positively correlated with target genes (p <0.05).
In vitro analyses
In caki-1 cells, CXCL2 genes were upregulated while CXCL5 did not exhibit significant variations when compared to 293T cells (Fig. 5 A). Moreover, WB analysis revealed differential CXCL2 protein levels between benign and malignant renal cells (Fig. 5 B). Analysis of the three CXCL2 siRNAs showed that siRNA-227 resulted in the lowest expressions of the CXCL2 gene. Therefore, caki-1 cells and caki-1 CXCL2 low cells (caki-1 cells mixed with siRNA-227) were used for subsequent experiments (Fig. 5 C).
The CCK assay showed that caki-1 CXCL2 low cells had weaker proliferative capacities, compared to caki-1 cells while the Transwell assay confirmed that caki-1 CXCL2 low cells had worse migration and invasive capacities, compared to caki-1 cells (Fig. 5 D-F). These findings imply that CXCL2 levels are upregulated in cancer cells and that CXCL2 expression levels affect cancer cell proliferation, migration as well as invasive capacities.
Upstream miRNA prediction and analysis
To identify the potential miRNAs upstream of CXCL2, 4 miRNA prediction databases, including miRDB, miRWalk, RNA22, and TargetScan were searched. Moreover, intersections of probable miRNAs were created using findings from the 4 databases. Six miRNAs upstream of CXCL2 were found.
Only one miRNA (hsa-miR-532-5p) was negatively correlated with CXCL2 (p = 4.87E-06, Fig. 6 A), and was found to be down-regulated in KIRC tissues (p <0.001, Fig. 6 B). Elevated hsa-miR-532-5p levels were associated with longer survival outcomes among KIRC patients (p = 0.003, Fig. 6 C).
Tumor-infiltrating immune cells
We investigated the association between CXCL2 levels and six immune cell components (B cells, CD4 T cells, CD8 T cells, Macrophages, Neutrophils, and DC cells) in KIRC patients using the TIMER database (Fig. 7 A-C). CXCL2 levels were positively correlated with CD4 T cells and Neutrophils (p<0.001), but were negatively correlated with B cells (p<0.01). However, COX analysis or KM analysis did not establish the survival value for CXCL2 levels combined with the three immune cell components in KIRC patients.
The risk score model
Cox proportional hazard regression was used to construct a CXCL gene signature with survival values. Six CXCL genes (CXCL1, CXCL2, CXCL3, CXCL5, CXCL6, and CXCL13) were obtained after univariate Cox regression analysis (p<0.01, Table 2), while 3 CXCL genes (CXCL1, CXCL5, and CXCL13) were obtained after multivariate Cox regression analysis (p<0.05, Table 3). The 3 CXCL genes (CXCL1, CXCL5, and CXCL13) were included in the survival signature and considered as the risk score (Fig. 7 B). The risk score = (0.115681*Exp CXCL1) + (0.060388*Exp CXCL5) + (0.118217*Exp CXCL13). KM analysis revealed that KIRC patients with lower risk scores had significantly longer survival times. At the same time, the association between the risk score and other clinical information was assessed. There was a strong relationship between the risk score and T stage, pathological stage, and pathological grade.