Clinical characteristics of the glioma patients

The microarray data of GSE104291 was downloaded from The Gene Expression Omnibus (GEO) database involving 4 glioma samples and 2 normal samples, while detailed clinical information were not available.

The RNA sequencing data across 33 tumor types and normal tissues of 118,103 samples were downloaded from The Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression (GTEx) database, including glioma (contained 689 tumor samples and 1157 normal samples). A total of 298 female and 398 male patients were involved in the present study, of which 20.5% (n = 143) were over age 60. As for race, most of the cases 637 (93.3%) were white. The WHO grade included 224 (35.3%) G2, 243(38.3%) G3 and 168 (26.5%) G4. IDH (isocitrate dehydrogenase) status involved 246 (35.9%) with IDH-wide-type (IDH-WT) and 440 (64.1%) with IDH-mutant (IDH-MT). 1p/19q codeletion status included 171 (24.8%) codel and 518 (75.2%) non-codel. In terms of primary therapy outcome, 112 (24.2%) were progressive disease (PD), 147 (31.8%) patients were stable disease (SD), 64 (13.9%) patients were partial response (PR), and 139 (30.1%) patients were complete response (CR). Besides, as for histological type, 195(28%) of patients were astrocytoma, 134 (19.3%) of patients were oligoastrocytoma, 199 (28.6%) of patients were oligodendroglioma and 168 (24.1%) were glioblastoma.

From The Chinese Glioma Genome Atlas (CGGA), 325 tumors with both gene expression data and clinical futures were analyzed. The clinical characteristics of the glioma patients including age, gender, histology, grade, IDH status, 1p/19q codeletion status, chemo-status, radio-status, and MGMTp_methylation status and etc. Among the 325 participants, 203 (62%) were male and 122 (38%) were female, and the median age of all participants was 40.5 years.

Abnormally high expression of CDC6 in glioma

As shown in Fig. 1A, principal component analysis (PCA) plot showed that no batch effect was found in the GSE10429 dataset. After the R language limma package processed, a total of 3742 differentially expressed genes (DEGs) were identified between glioma tissues and normal tissues, including 2179 downregulated genes and 1563 upregulated genes. The volcano plots presented the expression of DEGs (Fig. 1B), and among them, the expression level of CDC6 was significantly upregulated in the GSE10429 dataset.

Fig. 1

The expression of CDC6 in different human cancers and its relationship with glioma prognosis. A PCA plot of the data showing no batch effect in the GSE104291 dataset. Red nodes represent the tumor cluster while blue nodes represent the normal cluster. B The volcano plots of DEGs in GSE104291. Red nodes represent upregulated genes, blue nodes represent downregulated genes, and gray indicates genes with no differential expression based on the criteria of p value < 0.05 and |log2 FC| > 1, respectively. C Wilcoxon rank sum test was used to analyze the differential expression of CDC6 in glioma tissues and normal tissues with the data of the GTEx database as controls. D CDC6 expression in TCGA tumors and normal tissues with the data of the GTEx database as controls. E Kaplan–Meier analysis of association between CDC6 expression and glioma prognosis in TCGA. F Kaplan–Meier analysis of association between CDC6 expression and glioma prognosis in CGGA. PCA, principal component analysis; DEGs, differentially expressed genes; log2 FC, log2 foldchange; TCGA, The Cancer Genome Atlas; GTEx, Genotype-Tissue Expression; CGGA, Chinese Glioma Genome Atlas. ns, p ≥ 0.05, *p < 0.05, **p < 0.01, ***p < 0.001

We analyzed the expression of CDC6 in 689 glioma samples of TCGA database and 1157 normal samples of TCGA database combined GTEx database, and found that CDC6 was significantly high expressed in glioma samples (Fig. 1C; p < 0.001). We further expanded the number of samples to evaluate the mRNA expression level of CDC6 across pan-cancer in TCGA tumors with the data of the GTEx database as controls. As shown in Fig. 1D, compared with normal tissues, CDC6 was significantly upregulated in 28 of 33 cancer types, including bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), glioblastoma multiforme (GBM), head and neck squamous cell carcinoma (HNSC), kidney renal clear cell carcinoma (KIRC), renal papillary cell carcinoma (KIRP), brain lower grade glioma (LGG), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), ovarian serous cystadenocarcinoma (OV), prostate adenocarcinoma (PRAD), stomach adenocarcinoma (STAD), oral squamous cell carcinoma (OSCC) and so on (all p < 0.05), while it was downregulated in acute myeloid leukemia (LAML) (p < 0.05), no difference was found in kidney chromophobe (KICH), pheochromocytoma and paraganglioma (PCPG), mesothelioma (MESO) and uveal melanoma (UVM). The results indicated that the mRNA expression of CDC6 was also highly expressed across multiple types of cancer.

Relationship between expression of CDC6 and survival in glioma

Kaplan-Meier analysis revealed that glioma patients from the TCGA database with high CDC6 expression was correlated with poor OS (Fig. 1E; p < 0.001). Then, similar analysis in glioma patients from the CGGA dataset demonstrated that OS was significantly decreased in patients with high CDC6 mRNA expression compared with those with low CDC6 mRNA expression (Fig. 1F; p < 0.001). Furthermore, we investigated the correlations between CDC6 expression and prognosis (OS) in different clinical subgroups of glioma from TCGA database. The results showed that the higher expression of CDC6 had a worse OS in following clinical subgroups, including subgroup of age > 60 years (Fig. 2A; p < 0.001), subgroup of female (Fig. 2B; p < 0.001), subgroup of primary therapy outcome (PD) (Fig. 2C; p < 0.001), subgroup of IDH status (WT) (Fig. 2E; p < 0.01) and subgroup of 1p/19q codeletion status (non-codel) (Fig. 2F; p < 0.001). However, no difference in OS was observed upon histological type (Glioblastoma) (Fig. 2D; p = 0.167).

Fig. 2

Transcriptional expression level and prognosis value of CDC6 in glioma from TCGA database. A-I The mRNA expression levels of CDC6 were analyzed in glioma patients according to (A) Age, (B) Gender, (C) Race, (D) IDH status, (E) 1p/19q codeletion, (F) WHO grade, (G-I) histological type, respectively. J-O Overall survival analysis towards the CDC6 expression was performed in subgroups of glioma patients: (J) Age: > 60, (K) Gender: Female, (L) primary therapy outcome: PD, (M) histological type: Glioblastoma, (N) IDH status: WT, (O) 1p/19q codeletion: non-codel. TCGA, The Cancer Genome Atlas; IDH, isocitrate dehydrogenase; WHO, World Health Organization; PD, progressive disease; WT, wild-type. ns, p ≥ 0.05, ***p < 0.001

Association between CDC6 expression and clinicopathological features of glioma patients

To evaluate the association between CDC6 expression and clinicopathological features of glioma patients, the clinical characteristics of 689 glioma patients including RNA sequencing data from TCGA database were analyzed. Based on the mean value of CDC6 expression, the patients were divided into high- and low- CDC6 expression groups, we performed Wilcoxon rank-sum test, Kruskal-Wallis test and logistic regression analysis. As shown in Fig. 2G-O, Our results showed that higher expression levels of CDC6 were observed in patients with age > 60 years (Fig. 2G; p < 0.001), IDH status (WT) (Fig. 2J; p < 0.001), 1p/19q codeletion status (non-codel) (Fig. 2K; p < 0.001), WHO grade (G3/G4) (Fig. 2L; p < 0.001), and worse histological type (Fig. 2M-O; p < 0.001). However, no statistically significant correlation was found between the expression levels of CDC6 and other clinical pathological characteristics, such as gender (Fig. 2H; p ≥ 0.05) and race (Fig. 2I; p ≥ 0.05). These data indicate that high expression of CDC6 is significantly associated with age, IDH status, 1p/19q codeletion status, WHO grade and histological type, respectively.

Univariate analysis using logistic regression demonstrated that CDC6 expression was correlated with poor prognostic features in glioma patients. High CDC6 expression was significantly correlated with age (> 60 vs. <= 60: OR = 4.212, 95%CI = 2.795–6.490, p < 0.001), WHO grade (G4&G3 vs. G2: OR = 9.413, 95%CI = 6.408–14.077, p < 0.001), IDH status (WT vs. MT: OR = 15.032, 95%CI = 10.032–23.106, p < 0.001), 1p/19q codeletion status (non-codel vs. codel: OR = 5.413, 95%CI = 3.638–8.230, p < 0.001), primary therapy outcome (PD vs. CR&PR&SD: OR = 2.656, 95%CI = 1.720–4.118, p < 0.001) and histological type (glioblastoma vs. astrocytoma&oligoastrocytoma&oligodendroglioma: OR = 31.688, 95%CI = 16.724–68.190, p < 0.001). However, there was no significant difference in overall survival (OS) upon gender and race (all p > 0.05). Taken together, these results indicate that the glioma patients with high CDC6 expression are associated with worse clinicopathological characteristics and may serve as a biomarker of poor prognosis.

Cox univariate and multivariate analysis of prognostic factors in glioma

Univariate and multivariate Cox proportional hazard regression analyses were carried out with glioma patients from TCGA database. The univariate analysis indicated that high CDC6 expression was associated with the worse OS (HR = 4.608; p < 0.001). Other clinical parameters, such as gender (HR = 1.262; p = 0.062), age (HR = 4.668; p < 0.001), WHO grade (HR = 18.615; p < 0.001), IDH status (HR = 4.608; p < 0.001), 1p/19q codeletion status (HR = 4.428; p < 0.001), primary therapy outcome (HR = 3.542; p < 0.001) and histological type (HR = 9.114; p < 0.001) were also correlated with the worse OS time. Following multivariate Cox analysis, the results showed that gender (HR = 1.945; p = 0.004), age (HR = 4.689; p < 0.001), WHO grade (HR = 4.879; p = 0.006), IDH status (HR = 0.544; p = 0.026) and primary therapy outcome (HR = 3.643; p < 0.001) were independent prognostic factors in OS of glioma patients. However, we could not exhibit statistical significance of CDC6 expression in OS by multivariate Cox analysis (HR = 1.470; p = 0.078).

Functional enrichment analysis of DEGs

To elucidate the biological functions of CDC6, we analyzed the DEGs between high- and low- CDC6 expression groups based on the median CDC6 expression level. A total of 1357 DEGs were acquired with the threshold values of adjusted p value (p.adj) < 0.05 and |log2 FC| > 2, including 1321 upregulated genes and 36 downregulated genes, that were presented in volcano plots. Then, we conducted Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of DEGs, revealing that the primary biological process (BP) contained pattern specification process, embryonic organ development, regionalization, and anterior/posterior pattern specification. The cellular component (CC) was mainly enriched in collagen-containing extracellular matrix, kinetochore, condensed chromosome, centromeric region, and DNA packaging complex. The molecular function (MF) was primarily involved in receptor ligand activity, DNA-binding transcription activator activity, RNA polymerase II-specific, cytokine activity, and extracellular matrix structural constituent. Further functional enrichment analyses showed significantly enriched KEGG pathway in cytokine-cytokine receptor interaction, transcriptional misregulation in cancer, cell cycle, systemic lupus erythematosus and ECM-receptor interaction.

GSEA identifies CDC6-related signaling pathway

To identify CDC6-related signaling pathways in glioma, GSEA between high- and low- CDC6 expression datasets was conducted to reveal significant differences (p. adj < 0.05, false discovery rate (FDR) < 0.25) in enrichment of the Molecular Signatures Database (MSigDB) Collection (c5.bp.v7.2.symbols.gmt (Gene ontology) and c2.cp.v7.2.symbols.gmt (Curated)). In all, 3 GO items including single organism behavior, gated channel activity, cognition, transporter complex and ligand gated channel activity were showed significantly differential enrichment in CDC6 high expression phenotype. The results showed that the biological processes strongly associated with CDC6 were cell proliferation and immune-related pathways. 6 KEGG items that exhibited significantly differential enrichment in the CDC6-high expression phenotype were identified, including MAPK activation, DNA repair, P53 signaling pathway, focal adhesion, core matrisome and NF-κB activation. The results revealed that these pathways positively associated with CDC6 were responsible for regulation of cell proliferation/apoptosis and tumor invasion. Taking together, these findings suggest that these biological processes and signaling pathways of CDC6-high expression, which are critically important in development and metastasis of tumor, may be a potential target for the treatment of glioma.

The correlation between CDC6 expression and immune cell infiltration

We further explored the association between the expression level of CDC6 and immune cell infiltration level quantified by ssGSEA in glioma using Spearman correlation. The results showed that CDC6 expression was positively correlated with infiltration levels of Th2 cells, Macrophages, Eosinophils, etc., and negatively correlated with that of plasmacytoid dendritic cells (pDCs), natural killer (NK) CD56 bright cells, etc.

CDC6 promotes proliferation and inhibits apoptosis in glioma cells

To explore the effect of CDC6 on glioma progression, U87 and U251 cells were transfected with sh-CDC6–1004, and the transfection efficiency was confirmed by western blot analysis. The results showed that the expression levels of protein in the sh-CDC6–1004 group was significantly lower than that in the sh-NC group. CCK-8 was used to assess role of CDC6 on glioma cell proliferation, and the results showed that the proliferative viability of U87 and U251 cells in the sh-CDC6–1004 group was smarkedly lower than that in the sh-NC group, particularly at 72 h. Flow cytometry was further performed to analyze apoptosis of the transfected cells, the results indicated that, compared with the sh-NC group, the sh-CDC6–1004 group showed greater apoptosis of U87 and U251 cells. Our results suggest that CDC6 promotes proliferation and inhibits apoptosis of glioma cells.

CDC6 facilitates migration and invasion of glioma cells

To further clarify whether CDC6 affects the ability of migration and invasion of glioma cells, we performed the transwell assays. The results showed that the number of migrated and invaded U87 and U251 cells in the sh-CDC6–1004 group was smarkedly lower than that in the sh-NC group. These findings reversely suggest that CDC6 promotes migration and invasion of glioma cells.