This article is mentioned in:

Introduction

Nasopharyngeal carcinoma (NPC) is a malignant tumor of the nasopharyngeal epithelium with an uneven geographical distribution, including 81% of novel cases occurring in Asia; however, only 9% in Africa (1). Regarding the number of newly diagnosed patients, the top 5 countries include China, Indonesia, Vietnam, India and Malaysia, accounting for 67% of patients with NPC worldwide, in 2012 (1). Despite the therapeutic enhancements in the treatment of NPC, late diagnosis, poor prognosis and metastasis influence survival (2). Therefore, it is crucial to investigate novel diagnostics and therapeutics intended for patients with NPC. Advances in radiotherapy and comprehensive chemotherapy strategies have enhanced outcomes in patients with primary NPC (3,4). Furthermore, NPC presents with highly malignant recurrence with local tissue invasion and distant metastasis, the main reason for radiotherapy failure (5).

In the course of tumorigenesis and development, various types of genetic and epigenetic abnormalities, such as methylation and allele imbalance contribute to the pathogenesis of NPC (5,6). Previous studies have reported that the molecular mechanism of NPC has provided an important understanding and revealed the etiology of NPC (7,8). Nasopharyngeal carcinoma is caused by multiple genes and different environmental factors, such as mutations in cyclin D1 and X-ray repair cross complementing 1, and an excess intake of salted fish and alcohol (8–10).

Microarray technology has been widely applied in the examination of general genetic abnormalities (11–13). Furthermore, the Gene Expression Omnibus (GEO) database is increasingly recognized by researchers with an increasing interest in data mining and bioinformatics and a number of tumor-associated genes, including for NPC, being discovered (14,15). However, only a limited number of studies integrate these microarray datasets in identifying essential genes and pathway regulatory networks in NPC (6,14).

In the present study, differentially expressed genes (DEGs) were screened based on gene expression profiling data obtained from GEO. Consequently, two datasets comparing gene expression in tumor and non-tumor tissues were used and DEGs considered to be involved in NPC development were identified. Gene Ontology (GO) term enrichment, Kyoto Encyclopedia of Genes and Genomes (KEGG), protein-protein interaction (PPI) network and gene co-expression network analyses were performed to uncover the underlying mechanisms of NPC occurrence and development.

Materials and methods

Microarray data

Two gene expression profiles, GSE12452 (16) and GSE34573 (17), were downloaded from the GEO database using the following criteria: i) Containing NPC and non-tumor samples; ii) based on the Illumina HiSeq 2000 (Homo sapiens; GPL11154) chip platform; and iii) sample size >30.

The array data for GSE12452 contains 41 samples, including 31 NPC and 10 normal nasopharyngeal tissues (accessed on December 23, 2017) (16). The array data for GSE34573 includes 64 samples, of which, 18 samples included expression profiles from chip scan (platform: GPL570) and used for further analysis in the present study, which included 15 NPC tumor and 3 normal control (NC) tissues (accessed on December 23, 2017) (17).

Data processing

GEO2R (http://www.ncbi.nlm.nih.gov/geo/geo2r/), an interactive web tool, is used to compare two sets of data and is capable of analyzing GEO array data (18). GEO2R was applied in screening differentially expressed mRNAs in NPC and normal tissues. P<0.05 and |log2fold-change|≥1 were set as cut-off standards. A DEG heat map was subsequently generated using the online tool Morpheus (https://software.broadinstitute.org/morpheus/).

Analysis of functional and pathway enrichment

The GO database denotes a large number of gene annotation terms (19). KEGG is used to identify functional and metabolic pathways (20). Analyses were accessed from Database for Annotation, Visualization and Integrated Discovery (DAVID; v6.8; http://david.ncifcrf.gov/). The cut-off criterion was P<0.05.

Analysis of PPI network and hub proteins

Another online software, the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING; v11.0; http://string-db.org/), provides an important integration and assessment of the PPI (21). Cytoscape (v3.70) is a common open-source software tool that is useful in the visual evaluation of biomolecule interaction networks comprising of protein, gene and other forms of interactions (22). DEGs were mapped using STRING to assess the PPI information prior to visualization with Cytoscape. Nodes that are usually linked to the proteins are described as hub proteins and are crucial in the network (23).

Statistical analysis

The receiver operating characteristic (ROC) curve was constructed using SPSS 19.0 (IBM Corp., Armonk, NY, USA). P<0.05 was considered to indicate a statistically significant difference.

Results

Identification of DEGs

In the present study, two gene expression profiles were evaluated, including 46 NPC and 13 NC samples. On the basis of the GEO2R analysis, a total of 1,319 and 692 genes were identified for GSE12452 and GSE34573, respectively. Venn Diagram online software (http://bioinformatics. psb.ugent.be/webtools/Venn/) was used to determine the intersection of the two datasets as presented in Fig. 1. A total of 298 genes were identified as DEGs in the two sets, with 206 up- and 92 downregulated genes. Fig. 2 illustrates the DEGs expression heat map, including the top 40 genes.

Figure 1. Overlapping DEGs from the GSE12452 and GSE34573 microarray datasets. Venn diagram of (A) up- and (B) downregulated genes from the GSE12452 and GSE34573 datasets. DEGs, differentially expressed genes.

Figure 2. Heat map of the top 40 differentially expressed genes. Heat map for the (A) GSE12452 and (B) GSE34573 datasets. Each column represents a sample and each row a gene. Red indicates upregulated and blue downregulated genes.

Analysis of GO term enrichment

To functionally categorize the 298 DEGs, 82 enhanced GO terms and 7 KEGG pathways were obtained. Tables I and II summarize the enhanced top 10 GO terms by P-value. These genes may play a crucial role in the development of nasopharyngeal carcinoma. These GO terms were classified based on their association with biological processes, including response to epithelial cell differentiation, positive regulation of the intrinsic apoptotic signaling pathway, keratinocyte differentiation, cell division, mitotic nuclear division, G2/M transition of the mitotic cell cycle, mRNA export from the nucleus and DNA replication. In addition, there were numerous enriched GO terms associated with protein binding, activities of protein serine or threonine kinases, extracellular matrix (ECM) structural constituent and cell surface as presented in Fig. 3, divided by biological process, cellular component and molecular function. Fig. 4 illustrates the enriched KEGG pathways, including the phosphoinositide 3-kinase (PI3K)-protein kinase B (Akt) signaling pathway, ECM-receptor interactions, the cell cycle, arrhythmogenic right ventricular cardiomyopathy and focal adhesion pathways.

Figure 3. GO terms associated with the identified target genes. GO analysis of (A) up- and (B) downregulated target genes of the intersection of GSE12452 and GSE34573 according to biological process. GO analysis of (C) up- and (D) downregulated target genes of the intersection of GSE12452 and GSE34573 according to cellular component. GO analysis of (E) up- and (F) downregulated target genes of the intersection of GSE12452 and GSE34573 according to molecular function. GO, Gene Ontology; NADP, nicotinamide adenine dinucleotide phosphate.

Figure 4. Kyoto Encyclopedia of Genes and Genomes pathway analysis of the identified target genes. Pathways corresponding to the (A) up- and (B) downregulated genes. ECM, extracellular matrix; PI3K, phosphoinositide 3-kinase; Akt, protein kinase B.

Table I. Kyoto Encyclopedia of Genes and Genomes and GO functional annotation for upregulated genes (top 10).

Table I.

Kyoto Encyclopedia of Genes and Genomes and GO functional annotation for upregulated genes (top 10).

Category	Term	Count	P-value
GOTERM_CC_DIRECT	Extracellular exosome	76	<0.001
GOTERM_CC_DIRECT	Cytoplasm	74	<0.001
GOTERM_CC_DIRECT	Extracellular space	33	<0.001
GOTERM_CC_DIRECT	Vesicle	13	<0.001
GOTERM_BP_DIRECT	Oxidation-reduction process	12	0.0206
GOTERM_CC_DIRECT	Apical plasma membrane	11	<0.001
GOTERM_CC_DIRECT	Microtubule	10	0.0036
GOTERM_MF_DIRECT	Structural molecule activity	10	<0.001
GOTERM_CC_DIRECT	Cilium	9	<0.001
GOTERM_BP_DIRECT	O-glycan processing	7	<0.001

[i] GO, Gene Ontology.

Table II. Kyoto Encyclopedia of Genes and Genomes and GO functional annotation for downregulated genes (top 10).

Table II.

Kyoto Encyclopedia of Genes and Genomes and GO functional annotation for downregulated genes (top 10).

Category	Term	Count	P-value
GOTERM_MF_DIRECT	Protein binding	52	<0.001
GOTERM_CC_DIRECT	Nucleus	36	<0.001
GOTERM_CC_DIRECT	Cytosol	29	<0.001
GOTERM_CC_DIRECT	Nucleoplasm	24	0.0018
GOTERM_CC_DIRECT	Membrane	18	0.0154
GOTERM_MF_DIRECT	ATP binding	15	0.0056
GOTERM_BP_DIRECT	Cell division	11	<0.001
GOTERM_BP_DIRECT	Mitotic nuclear division	9	<0.001
GOTERM_CC_DIRECT	Cell surface	9	0.0028
GOTERM_CC_DIRECT	Centrosome	8	0.0029

[i] GO, Gene Ontology.

Module screening in the PPI network

Based on the STRING output, the top 10 hub nodes with higher degrees of interaction were screened, including cyclin-dependent kinase (CDK) 1, structural maintenance of chromosomes (SMC) 4, kinetochore-associated (KNTC) 1, kinesin family member (KIF) 23, aurora kinase A (AURKA), ATPase Family AAA domain containing 2 (ATAD2), NDC80 kinetochore complex component (NDC80), enhancer of zeste 2 polycomb repressive complex 2 subunit (EZH2), BUB1 mitotic checkpoint serine/threonine kinase (BUB1) and protein regulator of cytokinesis 1 (PRC1) (Fig. 5). The hub genes that were discovered in the PPI network analysis could play important roles in the aberrant signaling pathways and are potential targets for future research.

Figure 5. Protein-protein interaction network of the 298 differentially expressed genes. Yellow colored nodes represent hub genes.

Diagnostic values of CDK1, SMC4, KNTC1, KIF23, AURK and ATAD2 expression in NPC tissues

To verify the diagnostic value of the hub genes, expression levels in NPC tissues were evaluated using ROC curves. As presented in Fig. 6A and B, the area under the curve (AUC) for CDK1, SMC4, KNTC1, KIF23, AURKA and ATAD2 in NPC and non-cancerous nasopharyngeal tissue determined for the GSE12452 dataset were 0.948 [95% confidence interval (CI), 0.866–1.000; P<0.001], 0.961 (95% CI, 0.908–1.000; P<0.001), 0.955 (95% CI, 0.898–1.000; P<0.001), 0.968 (95% CI, 0.918–1.000; P<0.001), 0.958 (95% CI, 0.902–1.000; P<0.001) and 0.952 (95% CI, 0.888–1.000; P<0.001), respectively. In addition, as indicated in Fig. 6C and D, the AUC for CDK1, SMC4, KNTC1, KIF23, AURKA and ATAD2 in NPC and non-cancerous nasopharyngeal tissue for the GSE34573 dataset were 1.000 (95% CI, 1.000–1.000; P=0.008), 1.000 (95% CI, 1.000–1.000; P=0.008), 1.000 (95% CI, 1.000–1.000; P=0.008), 0.956 (95% CI, 0.858–1.000; P=0.015), 0.978 (95% CI, 0.913–1.000; P=0.011) and 0.822 (95% CI, 0.605–1.000; P=0.086), respectively. These results suggest that these genes may be diagnostic markers for nasopharyngeal carcinoma.

Figure 6. Diagnostic value of HUB genes. (A) ROC curves and (B) AUC values for the selected genes from the GSE12452 dataset. (C) ROC curves and (D) AUC values for the selected genes from the GSE34573 dataset. ROC, receiver operating characteristic; AUC, area under ROC curve; CI, confidence interval; DEG, differentially expressed gene; CDK1, cyclin-dependent kinase 1; SMC4, structural maintenance of chromosomes 4; KNTC1, kinetochore-associated 1; KIF23, kinesin family member 23; AURKA, aurora kinase A; ATAD2, ATPase Family AAA domain containing 2.

Discussion

In Southern China, NPC is a common malignant tumor (24,25). Sihui City in Southern China has the highest overall incidence rate, with 27.2 per 100,000 men and 11.3 per 100,000 women in 2003, while the average incidence rates in other parts of China are 4.31 per 100,000 men and 1.81 per 100,000 women, respectively (24,25). Advances in treatment, including radio- and chemotherapy, had an impact on clinical results, improving the long-term survival of patients with NPC (26). Despite the aforementioned, 70% of patients with NPC are diagnosed at mid to late stages and the 5-year survival rates for these patients are 40–70% (27). To gain insight into the pathogenesis of NPC, this study focused on DEGs in NPC compared with normal tissues based on two gene expression datasets obtained from the GEO database. A total of 206 up- and 92 downregulated DEGs were identified in the two datasets. It was indicated that DEGs, including CDK1, SMC4, KNTC1, KIF23, AURKA, ATAD2, NDC80, EZH2, BUB1 and PRC1, were hub genes in the PPI network analysis.

CDK1 was indicated as one of the most important genes in NPC compared with normal tissues. CDK1 is a well-preserved serine/threonine kinase, which forms multiplexes with its cyclin partners, contributing in the progression and regulation of the cell cycle (26). One previous study has reported that the CDK1-cyclin B1 signaling pathway is a target in NPC therapy (23). In previous research, overactivity of the CDKs were reported to be due to different genetic and epigenetic events in human types of cancer and inhibition causes cell cycle arrest and apoptosis (28). In a published study, CDK1 selectively targets human breast cancer cells, suggesting that further research into CDK1 suppression is required for a potential application in breast cancer treatment (29). In addition, CDK1 indicated to have a high AUC in the present study. Therefore, it was concluded that CDK1 may promote growth of NPC cells by regulating the cell cycle and it may serve as a biomarker.

SMC4 is a member of the SMC gene family (30). A study reported that SMC4 expression is substantially associated with tumor size, degree of tumor differentiation, progressive stages and vascular invasion of hepatocellular carcinoma (31). Clinical and experimental data propose that SMC4 stimulates the formation of colorectal cancer (32). In addition, SMC4 had a high AUC in the present study. Therefore, it is suggested that SMC4 may promote the development of tumors and may serve as a biomarker.

The ATAD2 gene is located on human chromosome 8q24.13 (33). ATAD2 is a member of the ATPase family. Caron et al (34) reported that ATAD2 is highly expressed in mouse testis and in many tumor cell lines, including leukemia, lymphoma, osteosarcoma, breast cancer, lung cancer and cervical cancer. Other studies have verified that high expression of ATAD2 is significantly associated with poor prognosis in patients with hepatocellular carcinoma, as well as ovarian cancer (35,36).

Furthermore, KNTC1, KIF23 and AURKA were reported to be the important hub proteins in the PPI network analysis, whereas KIF23 and AURKA had been implicated in tumor progression in the present study. The study by Takahashi et al (37) indicated that KIF23 may be a novel therapeutic target in malignant glioma. A previous study has examined and affirmed that KIF23 is a novel transcriptional target gene of the tumor suppressor protein p53 (38). It was reported that p53 decreases KIF23 mRNA and protein expression levels in various cell types. Previous studies have further reported that AURKA is a central kinase-encoding gene and overexpression causes amplification of centrosome, chromosomal instability, as well as transformation in mammalian cells (39,40). AURKA is associated with different types of cancer, including NPC (41), colorectal cancer (42), gastric cancer (43), esophageal cancer (44), multiple myeloma (45), breast cancer (46,47) and bladder cancer (48). These studies suggest that AURKA may be considered as a potential therapeutic target for NPC treatment. KNTC1 serves an important role in the mitotic checkpoint (49). It regulates the main mechanism that prevents chromosome instability usually observed in human tumors (50).

The present study indicated that a number of processes identified in the GO enrichment analysis were associated with NPC, including the physiological processes of cell growth and proliferation. These processes were associated with numerous cancers including lung cancer (51), breast cancer (52) and skin cancer (53), and any abnormality in these processes may cause cancer development. A total of 7 KEGG pathways were enriched and with P<0.05 in the current study. The ECM-receptor interaction pathway has been demonstrated to be involved in the tumorigenesis of various cancers (54–56) and in the current study, the KEGG pathway analysis indicated that these DEGs were also enriched in the ECM-receptor interaction pathway. The PI3K-Akt signaling pathway contributes to the regulatory function in cell survival, proliferation, metastasis, metabolism, angiogenesis and apoptosis (57–60). It is aberrantly expressed in many human tumors, including NPC, colon, breast and ovarian cancer (57–60). In the current study, it was established that the focal adhesion and the cell cycle were associated with the occurrence and development of tumors.

The development of bioinformatics tools has led to the development of various models and the establishment of databases (61–65). A number of them effectively predict potential miRNA-disease associations with high precision (61–65). These will provide novel directions for the diagnosis and treatment of patients with cancers, including NPC.

However, there are some limitations in the present study. First, different datasets in the GEO database microarray and platform were not unified. Second, it should be emphasized that the regulatory network and the potential mechanisms examined in the current study were solely based on bioinformatics analyses and were not verified by in vitro or in vivo experiments. Lastly, the datasets used in this study come from a single database and require further validation.

In summary, this study provided a bioinformatics analysis of DEGs, which may be involved in the development of NPC. Critical node proteins were identified, indicating an association with other proteins in the network, including CDK1, SMC4, KNTC1, KIF23, AURKA and ATAD2. These proteins may be used in NPC tumorigenesis through modulation of the cell cycle and nucleic acid metabolic processes, and may also be utilized as molecular biomarkers in the early diagnosis of NPC. There is a requirement for further experimental studies to validate the results of the present study.

Acknowledgements

Not applicable.

Funding

This work was supported by the National Natural Science Foundation of China (grant nos. 81460460, 81360405 and 81760542), the China Postdoctoral Science Foundation (grant no. 2016M602918XB) and the Research Foundation of the Science and Technology Department of Guangxi Province, China (grant nos. 2016GXNSFAA380252 and 2014GXNSFBA118114).

Availability of data and materials

The datasets generated and/or analyzed during the present study are available in the GEO repository, (https://www.ncbi.nlm.nih.gov/geo/).

Authors' contributions

RW designed the study. RW, KL and MK analyzed data and wrote the manuscript. ZZ and WQ analyzed data. All authors read and approved the final version of the manuscript.

Ethics approval and consent to participate

Not applicable.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

References