KIF20A as a potential biomarker of renal and bladder cancers based on bioinformatics and experimental verification

RCC and BC datasets

The search was conducted in GEO from December 31, 2000, to October 31, 2023, and the keywords for searching were “renal cancer”, “bladder cancer”, “clear cell renal cell carcinoma”, and “ccRCC”. The entry type was “series”. The study type was “expression profiling by array”. “Homo sapiens” was the selected item. We obtained RCC datasets GSE14762 [13] and GSE53757 [14] and BC dataset GSE121711 from GEO using configuration files generated from GPL4866, GPL570, and GPL17586. GSE14762 consisted of 10 RCC and 12 normal tissues, GSE53757 consists of 72 RCC and 72 normal tissues, and GSE121711 consisted of 8 BC and 10 normal tissues.

Batch effect removal and screening of DEGs

The inSilicoMerging package in R was used to merge the matrix. The removeBatchEffect function was used to remove batch effects. "Limma" was used to identify the DEGs. The DEGs were identified for tumor vs. control. The cutoff criteria for DEGs were p<0.05 and a fold change >3. The results were visualized using a volcano plot.

Weighted gene co-expression network analysis (WGCNA)

The first step was to preprocess the gene expression data to remove any technical variations. The next step was to construct a pairwise correlation matrix between all pairs of genes in the dataset. The correlation matrix was then transformed into an adjacency matrix using a soft-thresholding power function to emphasize strong correlations while downweighting weak correlations. Once the adjacency matrix was computed, modules or clusters of highly correlated genes were identified using hierarchical clustering or other clustering algorithms. These modules were assigned different colors for easy visualization. Module preservation analysis was conducted as follows: To ensure the robustness of the identified modules, preservation analysis was performed to evaluate the stability and reproducibility of the modules across different datasets or sample sizes. Module annotation was conducted as follows: After identifying the modules, the biological significance of each module was assessed by performing functional enrichment analysis to identify overrepresented biological pathways, Gene Ontology (GO) terms, or other annotation terms. The final step was to correlate the identified modules with the external traits of interest.

PPI network

STRING was performed to build a PPI network with the DEGs. Cytoscape software was used for visualization, and MCODE (https://cytoscape.org/) identified the hub module.

Functional enrichment analysis

We used the R package clusterProfiler (version 3.14.3) for enrichment analysis. Enrichment analysis was performed using the Metascape database. Functional enrichment analysis is a computational method used to determine whether a set of genes or proteins is enriched for particular biological functions, pathways, or processes. This type of analysis is commonly used in genomics and proteomics research to gain insights into the biological mechanisms underlying a set of experimental observations.

GSEA

The gene expression data were normalized, and genes with low expression were filtered out to reduce noise. A set of gene sets was defined, usually based on prior knowledge of biological pathways or molecular functions. Enrichment scores were calculated. Statistical significance testing was conducted as follows: The significance of the enrichment scores was assessed by permutation testing, comparing the observed scores to the distribution of scores obtained by randomly permuting the sample labels. This study aimed to determine the differential gene list input into the KEGG REST API.

Heatmap of gene expression

The R package heatmap was used to make a heatmap of the degree of expression of core genes in GSE14762 and GSE53757, visually displaying expression differences of core genes among RCC, BC and normal tissue.

Survival analysis

The GEPIA database (http://gepia.cancer-pku.cn/) derives its data from the TCGA (https://www.cancer.gov/ccg/research/genome-sequencing/tcga) and GTEx databases. We used GEPIA for survival analysis.

Immunoinfiltration analysis

The LM22 gene signature was used to define 22 immune cell subtypes. We applied an integrated bioinformatics approach and used the CIBERSORT (http://CIBERSORT.stanford.edu/) software package to analyze the batch-corrected merged matrices from GSE14762 and GSE53757. We selected samples with sufficient confidence levels based on a cutoff of p < 0.05.

Comparative toxicogenomics database analysis

Key genes were input into CTD (http://ctdbase.org/), the diseases associated with key genes were found, and a bar graph of gene expression differences was created using Excel.

Verification of the role of KIF20A in apoptosis by Western blotting

Normal renal cells (293 [HEK-293]) and renal cancer cells (Caki-1) were obtained from the National Biomedical Experimental Cell Resource Bank (Beijing, China). Cells were divided into 4 groups: normal renal cancer cells (Con), renal cancer cells (RC), RC/KIF20A-knockdown group (RC/KIF20A-KO), and RC/KIF20A-overexpression group (RC/KIF20A-OE). Sample preparation was conducted as follows: First, the protein of interest was extracted and purified from the sample. This involved homogenization, cell lysis, and/or chromatography, depending on the sample type. Protein separation was conducted as follows: Extracted protein was separated based on molecular weight using SDS-PAGE (sodium dodecyl sulfate-polyacrylamide gel electrophoresis). Protein transfer was conducted as follows: After separation, the proteins were transferred from the gel to a solid support, either nitrocellulose or a PVDF (polyvinylidene fluoride) membrane, using electroblotting. During this step, proteins were immobilized on the membrane. Blocking was conducted as follows: To prevent nonspecific binding of the primary antibody, the membrane was incubated with blocking solution, either 5% nonfat dry milk or 1% BSA in Tris-buffered saline (TBS) with 0.1% Tween-20, for approximately one hour at room temperature or overnight at 4° C. Primary antibody incubation was conducted as follows: The membrane was incubated with a primary antibody specific to the protein of interest diluted in blocking solution for approximately one hour at room temperature or overnight at 4° C. Secondary antibody incubation was conducted as follows: After washing, the membrane was incubated with a secondary antibody conjugated to an enzyme such as horseradish peroxidase (HRP) or alkaline phosphatase (AP), which recognizes and binds to the primary antibody. The membrane was then washed again to remove any unbound secondary antibody. Signal detection was conducted as follows: The membrane was treated with substrate solution specific to the enzyme conjugated to the secondary antibody. This produced a chemiluminescent or chromogenic signal that could be detected using X-ray film or a specialized imaging system. Data analysis was conducted as follows: Since the intensity of the signal is proportional to the amount of protein present in the sample, the signal intensity was quantified using image analysis software, and the results were compared between different samples or treatments.

Availability of data and materials

The datasets generated and/or analyzed during the current study are available from the corresponding author on reasonable request.

Differential gene expression analysis

A total of 1720 DEGs were found from the batch-corrected merged matrices of GSE14762 and GSE53757, according to preset cutoff values (Figure 1A). A total of 764 DEGs were identified from GSE121711 (Figure 1B).

Figure 1. Differential gene analysis. (A) DEGs in GSE14762 and GSE53757. (B) DEGs in GSE121711.

Functional enrichment analysis

GSEA

We used GSEA to perform enrichment analysis of the entire genome. The major enriched pathways were identified as organic acid metabolism, drug metabolism, mitochondria, and metabolism of cysteine and methionine (Figure 2A–2D).

Figure 2. Functional enrichment analysis. (A) BP, (B) CC, (C) MF, (D) KEGG analysis by GSEA.

Metascape enrichment analysis

In the Metascape enrichment analysis, the GO enrichment analysis results showed processes such as small molecule catabolic processes, carbon metabolism, and organic hydroxy compound metabolic processes (Figure 3A). Additionally, we also output the enriched network colored by enrichment items and p values (Figures 3B, 3C, 4A, 4B) to visually represent the associations and confidence levels of each enrichment analysis.

Figure 3. Metascape enrichment analysis. (A) Heatmap of enriched terms across input differently expressed gene lists, colored by p values, via Metascape; (B) Network of enriched terms colored by cluster identity, where nodes that share the same cluster identity are typically close to each other; (C) Network of enriched terms colored by p value, where terms containing more genes tend to have a more significant p value.

Figure 4. Metascape enrichment analysis. (A) PPI network drawn by Metascape. (B) The screened core gene group.

WGCNA

In the WGCNA, the soft threshold power was set to 6 (Figure 5A, 5B), and interactions were found between these modules (Figure 5C, 5D). A module-phenotype correlation heatmap (Figure 6A) and scatterplots of the GS versus MM correlation of relevant hub genes (Figures 6B–6E, 7, 8A) were generated.

Figure 5. WGCNA. (A) β=6, 0.85. (B) β=6, 73.85. (C, D) A hierarchical clustering tree of all genes was constructed, important modules were generated, and the interactions between these modules were then analyzed.

Figure 6. WGCNA. (A) Module-to-phenotypic correlation heatmap. (B–E) GS vs. MM correlation scatterplot.

Figure 7. The correlation scatterplot in WGCNA.

Figure 8. The conjoint analysis for the WGCNA and the DEGs. (A) Correlation analysis by WGCNA. (B) A Venn diagram was generated by intersecting the differentially expressed genes (DEGs) selected by WGCNA and DEGs and was then used to create and analyze a protein-protein interaction network.

We also generated a Venn diagram for the WGCNA with differential genes screened as DEGs and assessed the intersection (Figure 8B).

PPI network

The merged matrix results of GSE14762 and GSE53757 are shown in Figure 9A. The interaction network of GSE121711 shows KIF20A as the core gene (Figure 9B).

Figure 9. Construction and analysis of the protein-protein interaction (PPI) network. (A) GSEGSE14762, GSE53757. (B) GSE121711.

Survival analysis

Compared with normal kidney tissue, KIF20A in RCC had high expression (Figure 10A). Patients with lower KIF20A expression had good overall survival (Figure 10B). Compared with normal bladder tissue, in BC tissue KIF20A was highly expressed (Figure 10C). High KIF20A gene expression is associated with a worse prognosis.

Figure 10. (A) Box plot of KIRC. (B) Survival curve. (C) Box plot of BLCA. (D) Histogram of CTD.

CTD analysis

We found that the KIF20A gene is linked to liver cell carcinoma, necrosis, chemical- and drug-induced liver damage, weight loss, delayed effects of prenatal exposure, liver enlargement, inflammation, and proliferation (Figure 10D).

Gene expression heatmap

We visualized the gene expression heatmap for core genes in the samples, with Figure 11A showing the merged matrix results of GSE14762 and GSE53757. The results showed that KIF20A was highly expressed in RCC (Figure 11A). The core gene heatmap in GSE121711 showed that KIF20A was highly expressed in BC tissue (Figure 11B). The KIF20A gene was highly expressed in the RCC and BC samples, while it was expressed at low levels in normal tissue samples, indicating its potential regulatory role in these cancers.

Figure 11. Heatmap. (A) GSEGSE14762, GSE53757. (B) GSE121711.

Immune infiltration analysis

We conducted immune infiltration analysis of the batch-corrected merged matrix of GSE14762 and GSE53757 using the CIBERSORT package and obtained the proportion of immune cells in the whole gene expression matrix at 95% confidence (Figure 12A), as well as the heatmap of immune cell expression from the dataset (Figure 12B). We also performed coexpression correlation analysis on the infiltrating immune cells and obtained a coexpression pattern map between immune cell components (Figure 12C).

Figure 12. Immunoinfiltration analysis. (A) Whole-gene expression matrix results of proportions of immune cells. (B) Immune cell expression calorigram. (C) Plot of coexpression patterns between immune cell components.

Role of KIF20A in the pRB Ser 780/CyclinA signaling pathway

KIF20A expression in the RC group was upregulated, as shown by western blotting. The core proteins (including pRB Ser 780, CyclinA, E2F1, CCNE1, and CCNE2) in the pRB Ser 780/CyclinA signaling pathway were also upregulated in the RC group. Bcl2 expression was high in the RC group, and Bax expression in the RC group was downregulated, which showed that apoptosis was inhibited in the RC group. When KIF20A was knocked down, the core proteins (including pRB Ser 780, CyclinA, E2F1, CCNE1, and CCNE2) in the pRB Ser 780/CyclinA signaling pathway were also downregulated in the RC/KIF20A-KO group, and apoptosis was enhanced. However, when KIF20A was overexpressed in the RC tissues, the core proteins (including pRB Ser 780, CyclinA, E2F1, CCNE1, and CCNE2) in the pRB Ser 780/CyclinA signaling pathway were also upregulated in the RC/KIF20A-OE group, and apoptosis was inhibited. The results showed that KIF20A might inhibit apoptosis in RC via the pRB Ser 780/CyclinA signaling pathway (Figure 13).

Figure 13. KIF20A might inhibit apoptosis in RCs via the pRB Ser 780/CyclinA signaling pathway.