Comprehensive Analysis of Differential Gene Expression to Identify Common Gene Signatures in Multiple Cancers.

BACKGROUND With the development of research on cancer genomics and microenvironment, a new era of oncology focusing on the complicated gene regulation of pan-cancer research and cancer immunotherapy is emerging. This study aimed to identify the common gene expression characteristics of multiple cancers - lung cancer, liver cancer, kidney cancer, cervical cancer, and breast cancer - and the potential therapeutic targets in public databases. MATERIAL AND METHODS Gene expression analysis of GSE42568, GSE19188, GSE121248, GSE63514, and GSE66272 in the GEO database of multitype cancers revealed differentially expressed genes (DEGs). Then, GO analysis, KEGG function, and path enrichment analyses were performed. Hub-genes were identified by using the degree of association of protein interaction networks. Moreover, the expression of hub-genes in cancers was verified, and hub-gene-related survival analysis was conducted. Finally, infiltration levels of tumor immune cells with related genes were explored. RESULTS We found 12 cross DEGs in the 5 databases (screening conditions: "adj p<0.05" and "logFC>2 or logFC<-2"). The biological processes of DEGs were mainly concentrated in cell division, regulation of chromosome segregation, nuclear division, cell cycle checkpoint, and mitotic nuclear division. Furthermore, 10 hub-genes were obtained using Cytoscape: TOP2A, ECT2, RRM2, ANLN, NEK2, ASPM, BUB1B, CDK1, DTL, and PRC1. The high expression levels of the 10 genes were associated with the poor survival of these multiple cancers, as well as ASPM, may be associated with immune cell infiltration. CONCLUSIONS Analysis of the common DEGs of multiple cancers showed that 10 hub-genes, especially ASPM and CDK1, can become potential therapeutic targets. This study can serve as a reference to understand the characteristics of different cancers, design basket clinical trials, and create personalized treatments.

Background

Pan-cancer research [1] and cancer immunotherapy [2] are entering a new era of anti-cancer research based on bioinformatics analysis tools. Pan-cancer research can help find new cancer-related molecular genetic traits and establish a molecular-based cancer classification for suitable treatments. For example, the application of trastuzumab in the treatment of gastric cancer, breast cancer [3], and colorectal cancer [4] should depend on the expression level of Her-2. Recently, the inhibition of PD-1/PD-L1 has led immunotherapy to become the focus of anti-cancer treatment. Immune checkpoint inhibitor (ICI) is a systemic treatment that blocks the negative regulatory signaling pathway to achieve tumor elimination. The efficacy of ICI is associated with PD-1/PD-L1 expression [5,6] for multitype cancer rather than being limited to a single type. Immunological expression level is regarded as a new criterion to guide treatment. Additional criteria of anti-tumor research might be discovered in the future. Different cancers may have the same genetic amplification changes that affect tumorigenesis and progression, but these genes are still partly unknown. Analysis of genomes of multiple cancers is necessary to identify common DEGs and thus deeply and macroscopically elucidate the characteristics of some common signaling pathways that are associated with oncogenesis and cancer treatment. Using the gene database GEO/ONCOMINE/TCGA, we analyzed the gene expression in various types of tumor tissues and related normal tissues to find hub-genes.

Material and Methods

Data sources

The research data in this paper were obtained from the Gene Expression Omnibus (GEO) database [7] (). Verification and further mining of data were performed using the Oncomine database [8] () and the TCGA database [9] ().

Dataset screening

Datasets involving differences in gene expression between tumors and normal tissues were screened from the GEO database. The following selection criteria were considered to control the heterogeneity and ensure the quality of research: samples were obtained from Homo sapiens, the research platform was the common large platform GPL570, “Expression profiling by array” was adopted, the sample size was controlled at 50 and above, the cancer species were included in the common clinical cancer, and the data were published for nearly 15 years. On the basis of the above selection criteria, the following data sets were selected: GSE42568 [10], GSE19188 [11], GSE121248 [12], GSE63514 [13], and GSE66272 [14].

Differential expression gene extraction

Differentially expressed genes were extracted and analyzed using the online analysis tool GEO2R (), which is included in the GEO database. Samples in all the datasets ware defined as “normal” or “tumor” according to the actual situation, and the differential expression of each gene was obtained online. Differentially expressed genes (DEGs) were defined as “adj p<0.05” and “logFC>2 or logFC<–2”. A Venn diagram was drawn using the webtool () to obtain the cross DEGs.

Gene ontology

Gene ontology (GO) analysis is an important part of current functional genomics research, which refers to the high-throughput annotation of biological functions of all genes in the genome by using bioinformatics methods and tools. DAVID is a popular online programming feature () for GO analysis, including Biological Process analysis, Molecular Function analysis, and KEGG PATHWAY analysis [15]. The Kyoto Encyclopedia of Genes and Genomes (KEGG) () is a database that provides gene and genome functional significance at the molecular and pathway levels [16]. Multiple web tools such as DAVID can detect these contents of the KEGG database.

Protein–protein interaction network (PPI network) and hub-gene extraction

STRING () is an online searching tool for retrieving interacting genes that predicts the quality-controlled PPI networks [17]. We used STRING to construct a protein interaction network, and perform the gene co-expression analysis. Cytoscape is an open-source software platform for visualizing complex networks and integrating networks from the attributing data [18]. We visualized the results of the protein interaction network of STRING through Cytoscape. MCODE is an additional Cytoscape-based software that clusters a given network topologically to find areas of dense connectivity [19]. In the present study, we found more closely related genes by using MCODE (extraction conditions: degree cutoff: 10; K core: 2; max depth: 100), and named them the hub-genes.

Hub-gene verification

The Oncomine database () is one of the world’s largest oncogene chip databases and integrated data mining platforms [8]. It has comprehensive gene expression data and related clinical information. Some datasets were selected from the Oncomine database to test whether the selected hub-genes were differentially expressed in multitype tumors. The TCGA database is the world’s largest database of tumor information, and many online tools are available for extracting and analyzing TCGA data [9]. Whether such differences exist in the TCGA database was also investigated using online tools.

Hub-gene survival analysis

Gene Expression Profiling Interactive Analysis (GEPIA, ) is an interactive web resource for analyzing cancer transcriptome data, enabling researchers to collect valuable data on genes with interesting information [20]. We performed a survival analysis of hub-genes by using GEPIA.

Tumor immune infiltration levels

Tumor Immune Estimation Resource (TIMER, () is a web tool for the comprehensive analysis of tumor-infiltrating immune cells [21]. It includes 7 analysis modules that help researchers obtain the immunological, clinical, and genomic features of a tumor. The “SCNA” module of TIMER can compare tumor immune infiltration levels among tumors based on different somatic copy number alterations using the two-sided Wilcoxon rank sum test.

Results

Datasets for research

After online screening, GSE42568, GSE19188, GSE121248, GSE63514, and GSE66272 met our research request. Five types of tumors were collected from the 5 datasets. Among them, GSE42568 analyzed the difference in gene expression between breast cancer and normal tissues, and 121 samples were included. GSE19188 used 156 samples to analyze early non-small cell lung cancer and normal tissue. GSE121248 analyzed HBV-related liver cancer and normal tissue. GSE63514 analyzed the difference between cervical cancer and normal tissues. GSE66272 analyzed the genetic differences in normal tissues adjacent to renal cell carcinoma with 54 samples. The basic characteristics of the included datasets are shown in Table 1.

Table 1

Basic characteristics of the included datasets in GEO database.

Dataset	Contributor(s)	Organism	Submission year	Tumor type	Samples
GSE42568	Clarke C etc.	Homo sapiens	2012	Breast Cancer	121
GSE19188	Philipsen S	Homo sapiens	2010	NSCLC	156
GSE121248	Hui KM	Homo sapiens	2018	Hepatocellular carcinoma	107
GSE63514	den Boon J	Homo sapiens	2014	Cervical cancer	128
GSE66272	Wotschofsky Z etc.	Homo sapiens	2016	Kidney cancer	54

The DEGs of the 5 datasets were analyzed by GEO2R online analysis software, and the statistically significant gene names were screened by EXCEL (OFFICE 2016). Among them, 1618 DEGs were screened by GSE42568, 635 by GSE1918, 176 by GSE121248, 1147 by GSE63514, and 2185 by GSE66272. We made a Venn diagram with a webtool and obtained 12 cross DEGs. Moreover, all the 12 genes were upregulated in cancers (ANLN, CDK1, CYP2B7P///CYP2B6, ECT2, PRC1, NEK2, ASPM, RRM2, TOP2A, BUB1B, DTL, and CTHRC1). The Venn diagram is shown in Figure 1.

Figure 1

Venn diagram, 5 data sets were GSE42568, GSE19188, GSE121248, GSE63514, and GSE66272. A total of 12 DEGS (ANLN, CDK1, CYP2B7P///CYP2B6, ECT2, PRC1, NEK2, ASPM, RRM2, TOP2A, BUB1B, DTL, CTHRC1) were obtained.

Some difficulties were encountered in using DAVID for gene annotation. DAVID software cannot be used for gene annotation because of the small number of DEGs. STRING () [17] is a website for building protein networks, which also include some simple gene annotations. Therefore, we chose STRING to add nodes through intelligence. The DEGs were annotated separately from the biological process, molecular function, and the KEGG PATHWAY. The biological processes of DEGs were mainly concentrated in cell division, regulation of chromosome segregation, nuclear division, cell cycle checkpoint, and mitotic nuclear division. The main molecular function is to regulate cell division-associated protease activity. The signaling pathway was enriched in the Cell cycle, Progesterone-mediated oocyte maturation, oocyte meiosis, HTLV-I infection, and p53 signaling pathway. The results are shown in Table 2.

Table 2

GO and KEGG analysis. FDR – false discovery rate.

Pathway ID	Pathway description	Gene count	FDR
Biological process
GO.0051301	Cell division	17	1.55E-20
GO.0051983	Regulation of chromosome segregation	9	5.30E-14
GO.0000280	Nuclear division	13	1.23E-13
GO.0000075	Cell cycle checkpoint	11	2.35E-13
GO.0007067	Mitotic nuclear division	12	2.93E-13
Molecular functions
GO.0004748	Ribonucleoside-diphosphate reductase activity, thioredoxin disulfide as acceptor	2	0.00454
GO.0004672	Protein kinase activity	6	0.0184
GO.0004674	Protein serine/threonine kinase activity	5	0.034
Cellular components
GO.0005819	Spindle	11	1.61E-12
GO.0000922	Spindle pole	7	1.86E-08
GO.0015630	Microtubule cytoskeleton	12	1.86E-08
GO.0044430	Cytoskeletal part	13	1.93E-08
GO.0031981	Nuclear lumen	17	6.10E-08
KEGG pathway
4110	Cell cycle	9	8.76E-13
4914	Progesterone-mediated oocyte maturation	6	2.79E-08
4114	Oocyte meiosis	6	1.06E-07
5166	HTLV-I infection	6	1.43E-05
4115	p53 signaling pathway	3	0.00271

PPI network and co-expression analysis

The PPI constructed by STRING was visualized by Cytoscape. The results are shown in Figure 2. The co-expression analysis of the genes is shown in Figure 3 (by STRING). Moreover, we extracted the most closely related genes by using MCODE (degree >10). These genes were the hub-genes we would focus on (Table 3). Gene information was sourced from the GeneCard online website (). The 2 pairs between each hub-gene were expressed in a proportional manner. In combination with Figure 3, we suspected that ASPM and CDK1 were closely related genes.

Figure 2

Protein interaction network, obtained 10 Hub-genes (degree >10) by MCODE in Cytoscape.

Figure 3

Co-expression analysis of DEGs. Deeper color means stronger association. CDK1/ASPM is more closely related to other genes.

Table 3

Hub-gene information.

Gene	Full name	Function
TOP2A	DNA topoisomerase II alpha	This gene encodes a DNA topoisomerase which act on during transcription
ECT2	Epithelial cell transforming 2	The protein encoded by this gene is a guanine nucleotide exchange factor and transforming protein
RRM2	Ribonucleotide reductase regulatory subunit M2	This gene encodes one of two non-identical subunits for ribonucleotide reductase
ANLN	Anillin actin binding protein	This gene encodes an actin-binding protein that plays a role in cell growth and migration, and in cytokinesis
NEK2	NIMA (Never In Mitosis Gene A)-related kinase 2	This gene encodes a serine/threonine-protein kinase that is involved in mitotic regulation
ASPM	Abnormal spindle microtubule assembly	This gene is essential for normal mitotic spindle function in embryonic neuroblasts
BUB1B	Budding uninhibited by benzimidazoles 1 homolog beta	This gene encodes a kinase involved in spindle checkpoint function
CDK1	Cyclin dependent kinase 1	The protein encoded by this gene is essential for G1/S and G2/M phase transitions of eukaryotic cell cycle
DTL	Denticleless E3 ubiquitin protein ligase homolog	Its related pathways are DNA Double-Strand Break Repair and Translesion synthesis by Y family DNA polymerases bypasses lesions on DNA template
PRC1	Protein regulator of cytokinesis 1	This gene encodes a protein that is involved in mitotic cycle

We conducted a meta-analysis of other cancers in the Oncomine database, and found that 10 hub-genes were also differentially expressed between other cancers and normal tissues. All of these were expressed at low levels in normal tissues and were highly expressed in tumor tissues, including sarcoma [22], esophageal cancer [23], gastric cancer [24], colon cancer [25], pancreatic cancer [26], bladder cancer [27], thyroid cancer [28], oral cancer [29], glioma [30], ovarian cancer [31], and lymphoma [32]. The results are shown in Figure 4. We obtained the same results in the TCGA database through the GEPIA online tool (Figure 5). To observe the results intuitively, we used the TIMER- Diff Exp module to explore the expression of CDK1 and ASPM between tumor and normal tissues. The results showed that CDK1 and ASPM were significantly higher in tumor tissues than in normal tissues (Figure 6).

Figure 4

Meta-analysis of the expression of 10 Hub-genes based on the Oncomine database. 1–12 represents a different study. They were related to gastric cancer, esophageal cancer, intestinal cancer, oral cancer, lymphoma, and glioma. Redder color means greater difference in the expression of a gene in tumor and normal tissues.

Figure 5

Expression of 10 Hub-genes based on the TCGA database. Green squares from light to dark indicate gene expression from low to high. The tumor type abbreviation follows the TCGA database abbreviation principle and also applies to Figures 6, 8, and 9.

Figure 6

ASPM and CDK1 expression in multitype tumors. P-value significant codes: 0 ≤ *** <0.001 ≤ ** <0.01 ≤ * <0.05 ≤. <0.1.

We used the GEPIA online analysis tool to analyze the survival data in the TCGA database. More than 3500 patients with lung cancer, liver cancer, renal cancer, cervical cancer, and breast cancer were included. The results showed that in lung cancer, liver cancer, renal cancer, cervical cancer, and breast cancer, the prognosis of patients with high expression of the hub-gene was worse than that of patients with low expression of hub-genes (p<0.05) (Figure 7).

Figure 7

Survival analysis of Hub-genes in patients with lung cancer, liver cancer, renal cancer, cervical cancer, and breast cancer patients.

We compared the immune cell infiltration levels of CDK1 and ASPM in 5 cases: deep deletion, arm-level deletion, diploid/normal, arm-level gain, and high amplification, including B cell, CD8+ T cell, CD4+T cell, macrophage, neutrophil, and dendritic cells. In breast cancer, liver cancer, and lung cancer with high amplification of ASPM, the above immune cells had obvious infiltration, suggesting a possible correlation between ASPM and immunity (Figure 8). Because of the limited data on CDK1 amplification, no firm conclusion could be drawn (Figure 9).

Figure 8

Tumor immune infiltration level of ASPM in multitype tumors. P-value significant codes: 0 ≤ *** <0.001 ≤ ** <0.01 ≤ * <0.05 ≤. <0.1.

Figure 9

Tumor immune infiltration level of CDK1 in multitype tumors. P-value significant codes: 0 ≤ *** <0.001 ≤ ** <0.01 ≤ * <0.05 ≤. <0.1.

Discussion

We established a method to find common genes of multiple cancers by referring to previous studies, and we verified the relationship between DEGs and the clinical factors. These genes might play an important role in the development of cancers and warrant further research.

We found 12 cross DEGs, all of which were cancer-upregulated genes. We performed GO and KEGG analyses on the DEGs to further clarify the reasons for the differential expression. The results of GO analysis showed that the biological processes of DEGs were mainly concentrated in cell division, cell cycle process, cell cycle, DNA replication, regulation of cell cycle, and cell cycle checkpoint. Previous studies have reported that cell cycle processes and dysregulation of the mitotic cell cycle can play important roles in the occurrence or progression of cancers [33,34]. In addition, KEGG signaling pathway enrichment analysis suggested that these differential genes were mainly involved in cell cycle, progesterone-mediated oocyte maturation, oocyte meiosis, HTLV-I infection, and p53 signaling pathway. Human T cell leukemia virus type 1 (HTLV-1) is the first retrovirus found to cause adult T cell leukemia. HTLV-1 persistently infects CD4(+) T lymphocytes. HTLV-1 encodes 2 oncoproteins – Tax and HBZ – which are required to initiate cell transformation and maintain cell proliferation, respectively. The development of HTLV-1 cancers is driven by clonal selection and amplification, during which host and viral factors synergistically disrupt genomic stability, immune surveillance, and other cancer-suppressor mechanisms [35]. P53 is a tumor-suppressor gene. In all malignant cancers, more than 50% of the mutations occur in this gene [36]. CDK1 protein interacts with iASPP protein and affects the proliferation and apoptosis of CRC cells through the p53 apoptosis pathway [37].

The results of co-expression analysis showed that hub-genes CDK1 and ASPM are closely related to the expression of other hub-genes. CDK1 modulates the centrosome and mitotic cells, which play important roles in the control of the eukaryotic cell cycle. The CDK1/cyclin A complex controls G2, whereas the CDK1/cyclin B complex governs orderly G2/M transition (i.e., entry into mitosis and maintenance of the mitotic state) [38]. As a key regulator of the cell cycle, CDK1 is a potent therapeutic target for inhibitors in cancer treatment. Leucine zipper cancer-suppressor 1 gene may play a role in cell cycle control by interacting with the Cdk1/cyclinB1 complex (). In addition, several studies tried to targeted CDK1. Flavopiridol (alvocidib) is the first potent inhibitor of cyclin-dependent kinases to reach clinical trial, and it is demonstrated to have sequence-dependent cytotoxic synergy with chemotherapy agents [39]. Recent studies have demonstrated the single-drug activity of dinaciclib in the treatment of recurrent myeloma [40]. Other compounds, such as SELICICLIB [41] and MILCICLIB [42], also exhibit anti-cancer activity. ASPM plays a role in the regulation of mitotic spindles. Most previous studies focused on the relationship between ASPM and neurodevelopment [43]. Mutations in the ASPM gene are common causes of autosomal recessive primary microcephaly [44]. Recent studies have reported an association between ASPM expression level and some cancers [45-47], but the specific mechanism remains unclear. Pai et al. found that ASPM interacts with disheveled-3, which is a cardinal upstream regulator of Wnt signaling, to increase protein stability by inhibiting proteasome degradation, thereby enabling the Wnt-induced beta-catenin transcriptional activity in prostatic cancer cells [48]. In the present study, ASPM high amplification was associated with cancer immune cell infiltration. However, this study used a small number of samples, and literature to support the mechanism of ASPM and immune infiltration is lacking. Therefore, future anti-cancer research should explore ASPM functions.

Other hub-genes have also been linked to cancers. TOP2A encodes a DNA topoisomerase that acts during transcription. In prostate cancer, the high stage and high Gleason scores are always followed by TOP2A amplification or overexpression. In advanced prostate cancer, TOP2A amplification is associated with androgen resistance and poor survival [49]. ECT2 promotes the exchange of guanine nucleotides on small GTPases of Rho family members, such as RHOA, RHOC, RAC1, and CDC42 [50]. The expression of the genes increases at the beginning of the cell cycle and remains high in the G2 and M phases. Cancer-specific ECT2 gene amplification results in the overexpression of ECT2 in human cancers [50]. Some studies suggest that P53 is a novel molecule that upregulates ECT2 in gastric cancer cells, and its modification after translation plays a key role in the tumor cycle [51]. One of 2 distinct subunits of RNA reductase is encoded by RRM2. A number of in vivo and in vitro studies aimed at downregulating RRM2 in the treatment of malignant melanoma [52] and fibrosarcoma [53]. Silencing of ANLN and HSPA4L inhibits the proliferation, migration, and apoptosis of nasopharyngeal carcinoma cells. miR-497 is a potent cancer suppressor that inhibits cancer phenotype by targeting ANLN and HSPA4L in nasopharyngeal carcinoma [54]. Inhibition of phosphoinositide 3-kinase/AKT activity in NSCLC cells reduces the stability of ANLN, resulting in decreased nuclear ANLN levels [55]. Recently, cholangiocarcinoma with fewer anti-cancer targets has been found to increase the expression of NEK2 in a cancer-specific manner compared with normal fibroblasts. Expression of exogenous NEK2 does not affect the growth of cholangiocarcinoma cells, whereas inhibition of NEK2 expression by siRNA inhibits cell proliferation and induces cell death [56]. The protein encoded by BUB1B has been localized to the centromere and plays a role in inhibiting the late-promoting complex/loop, delaying the onset, and ensuring proper chromosome segregation. Spinal checkpoint function is impaired in many cancers. GO annotations related to DTL include ubiquitin-protein transferase activity [57]. Protein ubiquitination and degradation represent druggable vulnerabilities in cancer cells [58]. Solving this problem is of great significance in combating cancer. Eliminating the regulation of PRC1 leads to cell division defects, which promote chromosomal instability (CIN), leading to cancer heterogeneity and cancer evolution [59]. Targeting PRC1 in aneuploidy cancer may induce apoptosis by normalizing CIN or creating chromosomal disorders in genomically stable cancers.

Conclusions

This study analyzed multiple genomes of lung cancer, liver cancer, kidney cancer, cervical cancer, and breast cancer and found the gene co-expression characteristics of cancers of various tissues and organs. Results confirmed that some genes (TOP2A, ECT2, RRM2, ANLN, NEK2, ASPM, BUB1B, CDK1, DTL, and PRC1) were related to patient prognosis, and ASPM might be associated with immune cell infiltration. We suspect that these genes may be potential therapeutic targets. Exploration of genomes of multiple cancers provides ideas for a comprehensive understanding of cancer characteristics, and is helpful for basket clinical trials and personalized treatments. In addition, further research on genes and immunity is needed.