Identification of hub genes and pathways in adrenocortical carcinoma by integrated bioinformatic analysis.

Abstract Adrenocortical carcinoma (ACC), a rare malignant neoplasm originating from adrenal cortical cells, has high malignancy and few treatments. Therefore, it is necessary to explore the molecular mechanism of tumorigenesis, screen and verify potential biomarkers, which will provide new clues for the treatment and diagnosis of ACC. In this paper, three gene expression profiles (GSE10927, GSE12368 and GSE90713) were downloaded from the Gene Expression Omnibus (GEO) database. Differentially expressed genes (DEGs) were obtained using the Limma package. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were enriched by DAVID. Protein‐protein interaction (PPI) network was evaluated by STRING database, and PPI network was constructed by Cytoscape. Finally, GEPIA was used to validate hub genes’ expression. Compared with normal adrenal tissues, 74 up‐regulated DEGs and 126 down‐regulated DEGs were found in ACC samples; GO analysis showed that up‐regulated DEGs were enriched in organelle fission, nuclear division, spindle, et al, while down‐regulated DEGs were enriched in angiogenesis, proteinaceous extracellular matrix and growth factor activity; KEGG pathway analysis showed that up‐regulated DEGs were significantly enriched in cell cycle, cellular senescence and progesterone‐mediated oocyte maturation; Nine hub genes (CCNB1, CDK1, TOP2A, CCNA2, CDKN3, MAD2L1, RACGAP1, BUB1 and CCNB2) were identified by PPI network; ACC patients with high expression of 9 hub genes were all associated with worse overall survival (OS). These hub genes and pathways might be involved in the tumorigenesis, which will offer the opportunities to develop the new therapeutic targets of ACC.


| INTRODUC TI ON
Adrenal gland is an important endocrine organ of the human body and is also one of the most common organs with high tumour metastasis rate. According to the latest edition of pathological and genetic classification criteria of adrenal tumours published by World Health Organization (WHO), adrenal tumours can be divided into five categories: adrenocortical tumour, adrenal medullary tumour, extra-adrenal paraganglioma, secondary tumour and other adrenal tumour. 1 Adrenocortical tumour mainly includes adrenocortical adenoma (ACA) and adrenocortical carcinoma (ACC). 2 ACC is a rare malignant tumour originating from adrenal cortical cells, 3 and its incidence ranges from 0.7/10 00 000 to 2.0/10 00 000, which is under a high degree of malignancy, aggressiveness, high recurrence rate and poor prognosis. Most of the patients are found to have metastases and relapse easily after treatment, and the overall 5-year survival rate was <35%. 4 Early and accurate diagnosis is particularly important for the treatment and prognosis of ACC. 5 At present, surgical resection is the only feasible method to cure ACC, but it is difficult to control its quality. 6 Therefore, identifying new therapeutic targets or biomarkers for prognosis, diagnosis or prediction of ACC is urgently needed.
In recently years, many microarray profiling studies have been performed in ACC, 7,8 and hundreds of differentially expressed genes (DEGs) have been obtained. However, the results are limited or inconsistent due to molecular heterogeneity, and the results are usually generated from a single cohort study. Until now, no reliable biomarkers have been used in ACC clinics. Hence, the bioinformatics methods integrating multi-cohorts analysed by gene microarray or RNAseq will be innovative and valuable for ACC research.
In this work, we downloaded three different Gene Expression Omnibus (GEO) datasets (GSE10927, GSE12368 and GSE90713) and screened differentially DEGs using the Limma package.
Then, the PPI network in STRING database was constructed to screen the hub genes and pathways, and GEPIA database was used to verify hub genes and potential pathways. This study will offer the opportunities to develop the new therapeutic targets of ACC.

| MATERIAL S AND ME THODS
The flow diagram of this study was shown in Figure 1. The raw expression data were operated through a series of databases and software. This study has been proved by the Henan University institutional committee.

F I G U R E 1
The flow diagram of this study

| Data collection
Gene expression profiles of GSE10927, 9 GSE12368 10 and GSE90713 11 were obtained from GEO database. The GSE10927 dataset includes 55 neoplastic samples and 10 non-neoplastic samples (33 cases of ACC, 22 cases of ACA and 10 cases of normal adrenal cortex). The GSE12368 dataset is composed of 28 neoplastic samples and 6 non-neoplastic samples (12 cases of ACC, 16 cases of ACA and 6 cases of normal adrenal cortex). The GSE90713 dataset includes 58 neoplastic samples and 5 non-neoplastic samples (58 cases of ACC and 5 cases of normal adrenal cortex).

| Screening differentially expressed genes
DEGs between ACC samples and non-neoplastic samples were screened by using Limma package based on R language. DEGs were defined as representing differences with |log 2 FC| > 1, P < .05. 12 to understand biological meaning behind large list of genes. 16 We performed GO terms and KEGG pathway analysis of DEGs by using DAVID database.

| PPI network construction and hub module selection
Search Tool for the Retrieval of Interacting Genes (STRING) database (http://www.strin g-db.org/) was used to evaluate protein-protein interaction (PPI). 17 In addition, the database was used to quantify the relationships among the DEGs. Then, we used Cytoscape software to construct PPI network. 18 The genes with the highest node score and the strongest connectivity were selected. P < .05 was considered to have statistical significance.

| Hub genes validation
GEPIA (http://gepia.cance r-pku.cn/) is a powerful interactive web server that can analyse the RNA sequencing expression data of 9736 tumours and 8587 normal samples from the TCGA and the GTEx projects. 19 GEPIA can be directly used for tumour/normal differential expression analysis according to cancer types, and the box plot will be shown to visualize the relationship. GEPIA was used to verify the hub genes and perform validation, P < .05 showed statistical significance. 19 TCGA-ACC RNA sequencing data with patient survival data were downloaded from the University of California, Santa Cruz (UCSC) Xena browser. 20 Clinicopathological parameters of ACC patients with primary tumours, including age at diagnosis, gender, pathologic stage, living status, and overall survival (OS), were used for survival-curve analysis.

| Statistical analyses
Clinicopathologic parameter association analysis, and univariate and multivariate Cox regression analysis of 9 hub genes were performed with SPSS 22.0. Statistical significance was set at probability values of P < .05.

| Identification of differentially expressed genes
Using P < .05 and |log 2 FC| > 1 as cut-off, we identified 1040 up-

| Gene ontology analysis of differentially expressed genes
Subsequently, GO analysis of up-regulated and down-regulated DEGs was carried out by using DAVID online analysis tool. 16 In terms of biological process (BP), up-regulated DEGs were significantly en-

| KEGG pathway enrichment analysis of differentially expressed genes
KEGG pathway analysis of all DEGs showed that most of upregulated DEGs were enriched in cell cycle, cell senescence, progesterone-mediated oocyte maturation, oocyte, p53 signalling pathway and folic acid resistance (Figure 4), while down-regulation of DEGs did not significantly enrich KEGG pathway.

| PPI network construction and hub gene selection
Through analysing STRING database 17 (Table S1). In terms of molecular function, these hub genes are significantly enriched in ATP binding (Table S1). KEGG pathway enrichment analysis showed that these hub genes are associated with progesterone-mediated oocyte maturation, cell cycle, oocyte meiosis, p53 signalling pathway and progesterone-mediated oocyte maturation (Table S1).

| Evaluate the prognostic value of hub genes
TCGA-ACC dataset was used to evaluate the prognostic value of nine hub genes by GEPIA. All patients with high hub gene expression F I G U R E 4 Significantly enriched pathway terms of up-regulated differentially expressed genes in ACC F I G U R E 5 Protein-protein interaction (PPI) network, module analysis and hub gene identification. Red nodes represent up-regulated genes. Green nodes represent down-regulated genes. A, PPI network of differentially expressed genes was constructed in STRING database. B, Top nine hub genes were selected by Cytoscape software based on the degree of each node were associated with worse OS (Figure 6). The additional univariate and multivariate Cox regression analysis showed that the hub gene BUB1 (budding uninhibited by benzimidazole 1) was an independent prognostic factor for ACC patients, and BUB1 was significantly associated with living status and clinical stage in TCGA data (Table 1 and Table 2). The analysis results of the remaining genes were shown in Tables S2-S17. In addition, all nine hub genes were validated to be significantly up-regulated in ACCs as they were in above three GEO datasets (Figure 7). RACGAP1 was found to be highly expressed in colorectal cancer 36 and breast cancer. 37 Finally, we evaluated the prognostic value of 9 hub genes by using GEPIA database and found ACC patients with high expression of 9 hub genes (CCNB1, CDK1, TOP2A, CCNA2, CDKN3 In conclusion, using multiple cohort profiling datasets and integrated bioinformatics analysis, we identified hub genes and potential pathways that may be involved in the progress of ACC.

E THI C S APPROVAL AND CON S ENT TO PARTI CIPATE
This study has been proved by the Henan University institutional committee.

CO N FLI C T O F I NTE R E S T
The authors declare that there is no conflict of interests. No animal or human studies were carried out by the authors for this article.

AUTH O R CO NTR I B UTI O N S
XG and LX involved in study concept and design; LX, JG, YG and XG involved in acquisition of data; JG, YG, XM, LZ, HL, ZY and YH involved in analysis and interpretation of data. LX, JG and XG drafted the manuscript; LX, JG and XG involved in critical revision of the manuscript for intellectual content.