Identification of Novel Genes and Associated Drugs in Cervical Cancer by Bioinformatics Methods

Background

Cervical cancer is the most common malignant tumor of the female reproductive system. In recent years, cervical cytological screening, HPV testing, colposcopy, and histopathology have become widely available to facilitate early detection and treatment of cervical cancer [1]. However, the cure rate for cervical cancer remains low in many developing countries [2]. The majority of cervical cancers are linked to the human papillomavirus, which is a spherical DNA virus that causes the squamous epithelium of the skin’s mucous membranes to proliferate [3]. Despite the efforts of the human papillomavirus vaccination program, the incidence of cervical cancer remains very high [4]. Currently, cervical cancer screening, surgical techniques, radiotherapy equipment, chemotherapy, and other treatment methods have been relatively effective, but the clinical treatment effect for advanced and recurrent cervical cancer is still far from sufficient, resulting in poor prognosis of patients. With the development of modern molecular biology and genomics, targeted therapy has become one of the research hotspots in the treatment of advanced or recurrent cervical cancer, and has achieved good clinical efficacy in the treatment of cervical cancer. Current targeted therapies for cervical cancer mainly include antiangiogenic drugs, poly-ADP-ribose polymerase (PARP) inhibitors, immune checkpoint inhibitors, and epithelial growthfactor (EGF) receptor blockers [5]. Although these new prevention methods have made great strides in combination with traditional treatments, the 5-year disease-free survival rate for patients with advanced cervical cancer is only 45% [6]. A great deal of attention has been paid to the pathophysiology of cervical cancer, but the exact pathogenesis of cervical cancer is still not clear. Hence, an effective method is badly needed to explore the internal mechanism and targeted drugs affecting the occurrence of cervical cancer.

High-throughput sequencing has been widely used to study the pathological mechanism of diseases, especially in the prevention of cancer [7–9]. A variety of public databases are available for clinical and sequencing data download and analysis [10]. Reanalysis of the gene expression profiles associated with cervical cancer is important for understanding the pathogenesis and molecular mechanisms of cervical cancer [11]. Various studies have explored the pathogenic genes involved in cervical cancer. Dan et al identified PLOD2, ANLN, AURKA, and AR as possible therapeutic targets of cervical cancer through comprehensive analysis of cervical cancer-related gene expression profiles combined with bioinformatics [12]. Dan et al found that CDC45, GINS2, MCM2, and PCNA may be related to the occurrence and development of cervical cancer, and further verified that MCM2 knockdown can inhibit the proliferation of cervical cancer cells [13]. Although many studies have focused on certain genes in cervical cancer, it is still far from enough to reveal the underlying pathogenesis [14–16].Therefore, the present study used bioinformatics methods to find novel prognostic biomarkers of cervical cancer.

The DEGs in expression profiles GSE6791, GSE9750, GSE63514, and GSE67522 from the Gene Expression Omnibus were identified first in the present work. Then, the DEGs between cervical cancer and healthy groups were analyzed. GO annotation, KEGG pathway annotation, and PPI analysis were conducted among the DEGs using several bioinformatics methods. Further analysis of the clinical significance of these hub genes (including expression and survival analysis) via the Gene Expression Omnibus may provide new directions for understanding the role of these hub genes in the development and metastasis of cervical cancer. Discovery and development of new cancer drugs is expensive and time consuming. Our study aims to explore the targeted therapy of cervical cancer genes and related drugs, which will provide new ideas for drug use in the treatment of cervical cancer through the new use of old drugs.

Material and Methods

MICROARRAY DATA INFORMATION:

By searching for cervical cancer in the Gene Expression Omnibus (GEO; https://www.ncbi.nlm.nih.gov/geo/) database, 4 publicly available gene expression profiles – GSE6791 (GPL570 platform), GSE9750 (GPL96 platform), GSE63514 (GPL570 platform), and GSE67522 (GPL10558 platform) – were downloaded for this study. GSE6791 had 20 cancer samples and 8 normal samples. GSE9750 had 33 cancer samples and 24 normal samples, GSE63514 had 28 cancer samples and 24 normal samples, and GSE67522 had 20 cancer samples and 11 normal samples.

DATA IDENTIFICATION OF DEGS:

DEGs between cervical cancer specimens and normal specimens were identified via GEO2R online tools (https://www.ncbi.nlm.nih.gov/geo/geo2r/). The |logFC| ≥1.5 and adjusted P value <0.05 were set as the cut-off criteria. Then, the commonly occurring DEGs among the 4 datasets were screened in Venn software. The logFC ≤-1.5 represented downregulated genes, meanwhile logFC ≥1.5 represented upregulated genes among these DEGs.

FUNCTIONAL ENRICHMENT ANALYSIS:

Gene ontology analysis (GO) is a widely used method to identify biological properties of genes. It mainly includes the biological processes that genes participate in, the cell locations they are found in, and the molecular functions they play [17]. The Kyoto Encyclopedia of Genes and Genomes (KEGG) serves as a collection for assessing the signaling pathways of disease-related genes [18]. DAVID (https://david.ncifcrf.gov/) is a bioinformatics database that is used to analyze differential gene function and pathway enrichment analysis [19]. Therefore, the DAVID website was used for GO and KEGG analysis of differentially expressed genes (P<0.05).

PPI NETWORK OF DEGS:

The Retrieval of Interacting Genes (STRING: http://string-db.org) was used to determine the interaction among the DEGs [20]. The results with a minimum interaction score of 0.9 were analyzed and visualized in Cytoscape software (version 3.6.1; https://cytoscape.org). All parameters in cytoHubba were set by default. Nodes with disconnected PPI network are hidden. The top 10 hub genes from the PPI network were selected using the plugin cytoHubba [21] (score >15).

EXPRESSION ANALYSIS OF HUB GENES IN CERVICAL CANCER:

The expression of these hub genes in cervical tissue was examined using GEPIA2 (http://gepia2.cancerpku.cn/), a web-based tool that interactively analyzes gene expression profiles in cancer patients and normal patients from TCGA and GTEX data [22]. Kaplan-Meier curves were used to analyze the correlation between hub genes and overall survival in patients with cervical cancer [23]. The protein expression profile of these hub genes was obtained through the HPA database (https://www.proteinatlas.org/), which is a transcriptome- and proteomics-based technique that studies protein expression in different human tissues and organs at RNA and protein levels [24]. We further analyzed the expression of hub gene in cervical cancer tissues and normal tissues by using immunohistochemical data from the HPA database.

FUNCTION AND PROGNOSIS ANALYSIS OF HUB GENES IN CERVICAL CANCER:

cBioPortal (http://www.cbioportal.org) is a visual tool for studying and analyzing cancer genetic data [25]. The mutation rate and mutation form of these hub genes were analyzed using the cBioportal tool. The GEPIA2 database was again used to analyze the relationship between these hub genes and the clinicopathological staging of cervical cancer.

THE CONSTRUCTION OF DRUG-GENE INTERACTION:

QuartataWeb (http://quartata.csb.pitt.edu/) [26] is a comprehensive online platform for composite genomics, which can be used to search for and infer drug-target interactions, target cell association pathways, and drug-cell pathway associations. The top 10 hub genes were uploaded into the database to screen existing drugs or compounds.

STATISTICS ANALYSIS:

All analyses were performed using GraphPad Prism version 6.0 (GraphPad Software, USA). Data are presented as mean±standard deviation, and differences between groups were analyzed using the t test. P<0.05 was considered to indicate statistical significance.

Results

IDENTIFICATION OF DEGS:

DEGs (2103 in GSE6791, 815 in GSE9750, 1977 in GSE63514, and 427 in GSE67522) were identified by GEO2R (Figure 1A–1D). The overlap in the 4 datasets contained 78 genes, as showed in the Venn diagram (Figure 1E, 1F). As a result, 51 upregulated genes and 27 downregulated genes were obtained (Table 1).

FUNCTIONAL ENRICHMENT ANALYSIS OF DEGS:

The 78 overlapping DEGs were subjected to the GO enrichment and KEGG pathway analysis. In GO terms of Biological Process (BP) annotation (Figure 2A), it was mainly enriched in the “DNA replication initiation”, “cell division”, “cell proliferation”, “DNA replication”, “G1/S transition of mitotic cell cycle”, “mitotic nuclear divsion”, “mitotic cytokinesis”, and “G2/M transition of mitotic cell cycle”. In Cellular Component (CC) annotation (Figure 2B), it was enriched in the “nucleus”, “nucleoplasm”, “cytosol”, “midbody”, “centrosome”, “spindle”, “microtubule”, “spindle pole”, “nuclear chromosome, telomeric region”, “spindle microtubule”, “mitotic spindle”, and “MCM complex”. In Molecular Function (MF) annotation (Figure 2C), it was enriched in the “protein binding”, “single-stranded DNA binding”, “ATP binding”, “DNA helicase activity”, and “protein kinase binding”. Additionally, KEGG analysis showed that these DEGs were mainly enriched in the “cell cycle” and “DNA replication” (Figure 2D).

THE CONSTRUCTION OF PPI AND HUB GENES ANALYSIS:

To better explore the role of these DEGs in cervical cancer, we analyzed 78 DEGs using the STRING database (Figure 3A). There were 41 nodes and 218 edges enriched in the PPI network. The top 10 hub genes were selected using cytoHubba. As shown in Figure 3B, the top 10 hub genes were: CDK1, CDC20, AURKA, TOP2A, ASPM, NCAPG, KIF23, CENPF, KIF20A, and PRC1. All parameters in cytoHubba were set by default.

CLINICAL SIGNFICANCE OF 10 HUB GENES:

To further identify the reliability of the hub genes in cervical cancer, we searched the expression levels and prognostic information about these hub genes for cervical cancer via GEPIA2. All 10 hub genes (CDK1, CDC20, AURKA, TOP2A, ASPM, NCAPG, KIF23, CENPF, KIF20A, and PRC1) were reconfirmed via GEPIA2. Compared with normal cervical tissues, the expression level of these DEGs all was higher (Figure 4). Moreover, overall survival rate with higher CDK1 and PRC1 was significantly higher. However, lower AURKA was associated with higher overall survival rate (Figure 5A–5C), suggesting that CDK1, AURKA, and PRC1 play pivotal roles in the process of cervical cancer. Furthermore, the protein expressions of CDK1, AURKA, and PRC1 were significantly upregulated in cervical cancer compared to normal tissues. Among these highly expressed genes, CDK1 and AURKA spread to the cytoplasm in cervical cancer tissues, while PRC1 distribution showed no significant difference (Figure 6).

The occurrence of cancer depends on the accumulation of gene mutations. Due to the rapid cell division and abnormal DNA damage repair mechanism when cancer occurs, the genome of cancer cells becomes more unstable and prone to gene mutations. Therefore, we examined the mutation rate of these hub genes and analyzed their role in cervical cancer via the cBioPortal database. Our analysis of the mutation rate of the hub genes found that ASPM and CENPF had the highest genetic mutation rates, and CDK1, AURKA, and PRC1 were 0.7%, 2.9%, and 2.2%, respectively (Figure 7A). Among the hub gene mutations, CDK1 mutations were mainly shallow deletions, the AURKA and PRC1 mutations were mainly gene amplification (Figure 7B). Clinicopathological stage is an important index to evaluate the range of invasion and metastasis, as well as the degree of disease progression, metastasis, and prognosis of malignant tumors. Here, we examined 3 genes (CDK, AURKA, and PRC1) that were significantly associated with clinical prognosis. Notably, only AURKA expression increased with the expansion of cervical cancer invasion and metastasis (Figure 7C).

THE CONSTRUCTION OF DRUG-GENE INTERACTION:

Three genes associated with survival were uploaded into QuartataWeb for drug-gene interaction analysis. As shown in Table 2, AURKA was identified and matched with 21 predicted drugs that are FDA-approved: Fostamatinib, Regorafenib, Dasatinib, Brigatinib, Ponatinib, Pazopanib, Sorafenib, Imatinib, Nintedanib, Acetylsalicylic acid, Sunitinib, Bosutinib, Tamoxifen, Dabrafenib, Vandetanib, Lenvatinib, Cabozantinib, Baricitinib, Axitinib, Tofacitinib, and Midostaurin are shared by these genes. Thus, we speculate that these drugs could play a role in the treatment of cervical cancer.

Discussion

Cervical cancer is a malignant tumor originating from the cervix and is the most widespread malignant tumor of the female reproductive system. At present, cervical cancer ranks fourth among female malignant tumors in the world [27]. The development of cervical cancer is a complex process related to a variety of genetic and environmental factors. Although the development of the HPV vaccine has improved prevention and treatment of cervical cancer, the incidence of cervical cancer is still on the rise and its underlying pathogenesis remains unknown. Patients can be divided into 3 categories based on treatment: early, locally advanced, and advanced, which include stable, recurrent, or metastatic tumors. In addition to local or radical surgical treatment used to treat very early cervical cancer, radiotherapy and classical cytotoxic drugs are also widely used to treat patients with advanced cervical cancer [28]. Therefore, there is an urgent need to find reliable biomarkers and related drugs for diagnosis and prognosis of cervical cancer. In recent years, a large amount of microarray and sequencing data has been generated owing to the rapid development of bioinformatics, which provides comprehensive and convenient methods to screen genetic alterations and identify molecular mechanisms for the diagnosis, therapy, and prognosis of cancers.

In this study, GSE6791, GSE9750, GSE63514, and GSE67522 were analyzed for DEGs between cervical cancer tissues and normal tissues via GEO2R tools. A total of 78 DEGs (51 upregulated genes and 27 downregulated genes) were selected. Next, GO and KEGG analyses were performed using the DAVID website. In GO terms of Biological Process (BP) annotation, it was mainly enriched in the “mitotic nuclear divsion”, “cell proliferation”, “cell division”, and “DNA replication”. In Cellular Component (CC) annotation, it was enriched in the “nucleus”, “nucleoplasm”, “cytosol”, and “midbody”. In Molecular Function (MF) annotation, it was enriched in the “protein binding”, “protein kinase binding”, “single-stranded DNA binding”, and “ATP binding”. Additionally, in the KEGG pathway, it was enriched in “cell cycle” and “DNA replication”. Moreover, 78 DEGs were analyzed in the STRING database, and the top 10 hub genes (CDK1, CDC20, AURKA, TOP2A, ASPM, NCAPG, KIF23, CENPF, KIF20A, and PRC1) were selected using Cytoscape software.

A total 8 hub DEGs were enriched in the cell cycle-related pathways. The operation of the cell cycle is regulated by cyclin-dependent CDK kinases, and abnormal proliferation of the cell cycle can lead to carcinogenesis. Cyclin therapy is treatment of the cell cycle, and cyclin-dependent stimulated CDKs, Polo-like kinases, and WEE kinases are the main targets of cell cycle therapy drugs [29,30]. Cell cycle therapy can also activate tumor immunity, control metabolic function, and regulate transcription levels [31]. Patients’ responses to cell cycle therapy support the significance of cell cycle in the treatment of cervical cancer [32,33].

In the bioinformatics analysis, 10 hub genes were selected and verified based on TCGA. Finally, we identified 3 genes (CDK1, AURKA, and PRC1) as potential targets for predicting cervical cancer prognosis. High expression levels of CDK1 and PRC1 are associated with better clinical outcomes, while high expression of AURKA predicts poor survival in cervical cancer patients. The expression levels of these 3 genes were verified by immunostaining data from the HPA database. To confirm the important role of these 10 central genes in cervical cancer, mutation rates of these genes were further studied online in the cBioPortal. The results revealed that AURKA and PRC1 had higher mutation rates. Among the 3 central genes, only AURKA was positively correlated with clinical stage, while the other 2 had no statistically significant difference. Finally, based on the above analysis, AURKA is determined to be the key functional protein in these candidate genes.

AURKA is a serine/threonine kinase that plays an important regulatory role in cell mitosis. AURKA is overexpressed in several human malignancies and is the key to centrosome replication, maturation, and separation. Previous studies have shown that AURKA is elevated in most human cancers and is associated with adverse clinical outcomes in patients with lung cancer [34], liver cancer [35], and colorectal cancer [36]. A large number of studies have confirmed that AURKA, as an oncogene, affects many tumor biological processes, such as proliferation, genomic instability, aneuploid karyotype production, invasion, and metastasis. Xiaojiang et al [38] found that palmatine, an inhibitor of AURKA, can induce G2/M phase arrest and mitochondria-induced endogenous apoptosis of colon cancer cells by binding to AURKA. Silencing of AURKA resulted in the degradation of p27 and Bax in the cytoplasm, which decreased the viability of gastric cancer cells [39]. Shao et al reported that AURKA expression was significantly increased in tissues of cervical cancer patients. Further studies showed that hsa_circ_0075341 enhances the proliferation and invasion ability of cervical cancer cells by isolating miR-149-5p and upregulating AURKA expression [40]. Robert et al found that CDK1 activation initiated mitotic activation that AURKA was involved in. In addition, activated AURKA and PLK1 further react on CDK1 to promote mitosis and participate in spindle formation and centromere division in mitosis [41]. In the preset study, we found that the upregulation of AURKA was closely related to the survival prognosis and TNM clinical stage of cervical cancer, suggesting that AURKA is a potential target for clinical treatment of cervical cancer.

Using integrated bioinformatic analyses, AURKA was identified to play an important role in the development of cervical cancer, but few drugs designed to target this biomarker have been used to treat cervical cancer. In the present study, AURKA matched 21 FDA-approved drugs. Among the 21 drugs, the top 9 drugs – Fostamatinib, Regorafenib, Dasatinib, Brigatinib, Ponatinib, Pazopanib, Sorafenib, Imatinib, Nintedanib, and Acetylsalicylic acid (confidence score >0.3) – were analyzed. Fostamatinib is an oral tyrosine kinase inhibitor, which is involved in B cell receptor signaling and may be involved in regulating autoantibody production [42]. Regorafenib is a novel multikinase inhibitor that targets tumor genesis, tumor angiogenesis, and maintenance of signaling in the tumor microenvironment by inhibiting multiple protein kinases that promote tumor growth [43]. Dasatinib [44], Ponatinib [45], and ImatinibJournal [46] are tyrosine kinase inhibitors that effectively lock tyrosine kinase into inactive conformation by binding to highly conserved ATP binding sites, thereby inhibiting BCR-ABL activity. They can inhibit multiple oncogenic pathways activated by BCR-ABL tyrosine kinase and promote cell proliferation and survival, including the MAPK/ERK pathway, JAK-Stat pathway, and PI3K/Akt pathway. Brigatinib is an anaplastic lymphoma kinase inhibitor for the treatment of anaplastic lymphoma kinase-positive metastatic non-small cell lung cancer [47]. Pazopanib is a small-molecule inhibitor of a variety of protein tyrosine kinases, with potential antitumor activity [48]. Sorafenib works by stopping the progression of cancer and blocking the therapeutic replication of potentially malignant cells [49]. Nintedanib is a small-molecule kinase inhibitor that competitively inhibits both non-receptor tyrosine kinases and receptor tyrosine kinases, and these inhibitory effects increase fibroblast proliferation, migration, and transformation [50].

At present, the main drugs used for treatment for cervical cancer are Bevacizumab, Sunitinib, Olaparib, and cetuximab. In addition, Pazopanib and Nintedanib, obtained by our analysis, have also been widely used in the treatment of cervical cancer due to their anti-angiogenesis and tumor development inhibition [5,51]. However, with increased efficacy of these drugs, their adverse effects gradually increase. Therefore, we need to find more effective target genes and related drugs.

As mentioned above, the hub genes related to cervical cancer are mostly concentrated in the cell cycle pathway, and AURKA has been identified as a key gene affecting the occurrence and development of cervical cancer. Thus, we analyzed the mechanism of action of AURKA-related drugs in tumor therapy and found that these drugs are mostly targeted kinase inhibitors, which block tumor growth by affecting cell proliferation-related signaling pathways. Although the role of these drugs in the treatment of cervical cancer remains to be further determined, our results provide evidence of the important role of cell cycle signaling in the development of cervical cancer and once again support the possible efficacy of molecular-targeted therapies for cervical cancer.

Nevertheless, there exist some inevitable difficulties in this study that should be taken into consideration. For example, did not perform additional molecular and cellular biology experiments to verify their properties. However, our laboratory is currently carrying out relevant experiments to explore the role of AURKA in cervical cancer. The current analysis of samples collected from patients with cervical cancer is not statistically significant, so more samples need to be collected for verification. Although experimental data obtained from existing databases reduce confidence, the results of this study provide preliminary insights into the cell cycle signaling pathway involved in AURKA and its role in cervical cancer development, as well as AURKA gene-related targeting agents. This once again supports the use of cell cycle biologics and related potential drugs in the treatment of cervical cancer.

Conclusions

In summary, a total of 78 DEGs, including 51 upregulated genes and 27 downregulated genes, were identified in cervical cancer by integrated bioinformatic analysis. The further analysis identified that AURKA may account for the tumorigenesis and prognosis of cervical cancer. The potential drugs associated with AURKA gene were identified, which will provide new ideas for the treatment of cervical cancer. Although 4 data sets were been obtained and analyzed strictly and comprehensively, further experiments and molecular mechanism exploration are still needed to clarify how AURKA and its related hub genes affect cervical cancer proliferation. Our study provides a new direction for understanding the molecular mechanism and promoting clinical drug research for treatment of cervical cancer.