Identification of 11 Differentially Expressed Hub Genes and Their Upstream microRNAs in Nasal-type NK/T-cell Lymphoma Based on Clinical Sample Analysis

Nasal-type NK/T-cell lymphoma (nasal-type NKTL) is one of the most lethal cancers for. Our study aimed to identify hub differentially expressed genes (DEgenes) and their upstream microRNAs between nasal-type NK/T-cell lymphoma (NKTL) tumor samples and normal nasal tissues through integrated bioinformatics. The 503 DE-genes and 106 DE-miRNAs were identified between NKTL and human normal nasal samples. GO and KEGG analysis were significantly enriched in meiotic recombination, regulation of syncytium formation by plasma membrane fusion, deubiquitination, enriched in meiotic recombination, regulation of syncytium formation by plasma membrane fusion, and stem cell division. And 11 differential expression hub genes and their upstream microRNAs were identified between nasal-type NKTL and normal nasal samples. In summary, after a series of analyses, we found that 11 hub DE-genes and their upstream DE-miRNAs (CDC27miR-548c-3p, FREM2miR-373 * , ARHGAP29-miR-548c-3p, QSER1-miR-548c-3p, CD3EAP-miR-149 * , SF3A1miR-548c-3p, AQP4-miR-29b, ZFP36L2-miR-142-3p, SRP72-miR-16, TSC22D2-miR-16, TSC22D2-let-7f, DOCK5-miR-16) between nasal-type NKTL and normal nasal samples. They are highly likely to be serve as promising biomarkers in nasal-type NKTL.


Introduction
Non-Hodgkin lymphoma (NHL) is a highly heterogeneous proliferative disease that can further develop into B cell lymphoma (accounting for 80%-85% of NHL), and natural killer/T-cell lymphoma (accounting for 15%-20% of NHL). Nasal-type NKTL is a rare subtype of mature T-and NK-cell lymphoma leading to palatal deformity caused by cancer invasion through the hard palate, which accounts for approximately 10% of all peripheral NK/T-cell lymphomas in America [1]. Although the use of CCRT with platinum or L-asparaginase-containing chemotherapy and the avoidance of anthracyclines can improve the prognosis of NKTL, patients still have a poor median survival time of less than 2 years according to global NKTL patient statistics [2]. More unfortunately, nearly half of the patients with newly diagnosed NKTL continue to experience disease progression, and the prognosis in patients with relapsed or refractory NKTL remains unsatisfactory [3]. Previous studies have reported that genetic aberrations, the cellular microenvironment, EBV infection, and epigenetic dysregulation are closely related to NKTL occurrence, development, and metastasis [4]. A better understanding of the pathogenic mechanisms, especially those of nasal-type NKTL, will provide important insights into the biology of this disease. Therefore, it is urgently necessary to explore the molecular mechanisms underlying nasal-type NKTL occurrence, development, and metastasis to develop targeted therapies and improve prognosis.
MicroRNAs (miRNAs) are a class of noncoding single-stranded RNA molecules with hairpin structures of approximately 22 nucleotides that are involved in regulating the expression of protein-coding genes [5]. At the same time, miRNA regulates gene expression at the posttranscriptional level by pairing with the 3'-untranslated region (3'UTR) of the target gene mRNA. Each miRNA can have multiple target genes and several miRNAs can regulate the same gene. In recent decades, miRNAs have been considered to be involved in the regulation of all aspects of cancers, including proliferation, differentiation and metabolism [6]. However, few studies based on large-scale human tissue studies of differentially expressed miRNAs focused on nasal-type NKTL have been reported.
As gene-chip and RNA sequencing technology are widely used to identify genomic changes during tumour progression, the Gene Expression Omnibus (GEO) has played an important role in bioinformatics analysis. This provides an efficient tool to discover differentially expressed genes (DE-genes) and differentially expressed miRNAs (DE-miRNAs). Based on the GEO database and comprehensive bioinformatics analysis, this study focused on exploring the DE-genes in nasal-type NKTL and their upstream miRNAs, as well as the possible molecular mechanisms in the development of nasal type NKTL compared with normal human nasal tissues. Therefore, we mined and analysed the data from multiple datasets in the GEO database to highlight the characteristics of gene expression and their upstream miRNAs and provide a basis for further research in human nasal-type NKTL.

Microarrays and Screening for Differentially Expressed Genes and Differentially Expressed miRNAs
We logged in to the National Center for Biotechnology Information (NCBI) GEO database (https://www.ncbi.nlm.nih.gov/geo) to search for the microarrays we needed. Then, to ensure the rigor and accuracy of our research, only datasets comparing mRNA/miRNA expression in human nasal-type NKTL samples with that in human normal nasal samples were included. Carefully read the complete information of the searched datasets and then selected the most suitable datasets for further analysis. Then, the GSE80631 [7] and GSE31377 [8] datasets were selected in our study. The dataset GSE80631 is a GPL6883 platform (Illumina HumanRef-8 v3.0 expression beadchip)-based dataset. We selected 13 human normal nasal samples and 19 human nasal-type NKTL samples from dataset GSE80631. GSE31377 is based on the platform of GPL8227 (Agilent-019118 Human miRNA Microarray 2.0 G4470B (miRNA ID version)). We selected 2 human normal nasal samples and 17 human nasal-type NKTL samples from dataset the GSE31377.
We used the "LIMMA" from the R package of the bioconductor project to normalize the data. To conduct the differential expression analysis of the miRNAs and genes in nasal-type NKTL samples compared to normal nasal tissues, the related codes were entered into R, and the DE-genes and DE-miRNAs were analysed through the "limma" package in the bioconductor package (http://www. bioconductor.org/). p-value<0.05 and |logFC| >1 were set as the cut-off criteria for identifying DE-genes and DE-miRNAs. At the same time, volcano plots of the DE-genes and DE-miRNAs were generated by R.

Prediction of Target Genes for DE-miRNAs
The potential target genes of the top 10 most upregulated and downregulated DE-miRNAs were obtained from the overlapping part of the predictive results derived from TargetScanHuman [9]  The TargetScanHuman Release 7.2, which was updated on March 2018, is a web server that predicts the biological targets of miRNAs by searching for the presence of conserved sites that match human 3'UTRs and the orthologues of each miRNA. The miRTarBase Release 7.0 is an experimentally validated microRNA-target interactions database, which is updated on September 15, 2017. miRDB Version 5.0 was released in August 2014, and is an online database for miRNA target prediction and functional annotations.

Enrichment Analysis
Metascape3.0 (http://metascape.org) [12] was introduced to perform functional and pathway enrichment analysis of the DE-genes and the predicted target genes of the selected top 10 upregulated DE-miRNAs and top 10 downregulated DE-miRNAs, including determination of the Gene Ontology (GO) terms for the biological process (BP), cellular component (CC) and molecular function (MF) categories, Reactome Gene Sets, Canonical Pathways, CORUM, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways.

Conjoint Analysis of DE-genes and DE-miRNAs
Due to inhibition by miRNAs of posttranscriptional gene expression via their specific binding to target messenger RNA (mRNA), 10 upregulated DE-genes were crossed with the predicted target genes of the top 10 downregulated DE-miRNAs, and 10 downregulated DE-genes were crossed with the predicted target genes of the top 10 upregulated Nasal-type NK/T-cell Lymphoma Based on Clinical Sample Analysis DE-miRNAs by calculating and drawing custom Venn diagrams in the web tool Bioinformatics & Systems Biology to obtain more convincing DE-miRNAs, DE-genes, and the one-to-one relationships.

Identification of Differentially Expressed Genes and miRNAs
To identify DE-miRNAs and DE-genes between nasal-type NKTL patient samples and normal nasal tissues, the GSE80631 and GSE31377 data were normalized and differential expression analysis was conducted using the "limma" software package. Then, visualization of normalized data is shown in Figure 1. Based on this differential expression analysis and our screening criteria (p-value<0.05, |logFC| >1), a total of 503 genes were found to be significantly differentially expressed in nasal-type NKTL samples when compared to normal nasal tissues, including 281 upregulated and 222 downregulated genes. At the same time, 106 significantly DE-miRNAs were discovered, of which 57 were upregulated and 49 were downregulated. Furthermore, the top 20 upregulated DE-genes and downregulated DE-genes (Appendix Table 3), downregulated DE-miRNAs and upregulated DE-miRNAs (Appendix Table  4) were ranked by |logFC| in Appendix Tables 3-4.

Pathway and Process Enrichment Analysis
Through Metascape3.0, enrichment analysis was carried out for DE-genes and the predicted target genes of DE-miRNAs. The networks of enriched terms were visualized using Cytoscape (Figures 2 and 3), where each node represents an enriched term and is coloured first by its cluster ID (Figures 2A1, B1, 3A1 and B1) and then by its p-value (Figures 2A2, B2, 3A2 and B2). More details on the pathway and process enrichment analyses are shown in Appendix Tables 7-8. The top 5 results for the pathway and process enrichment analyses for the upregulated DE-genes included meiotic recombination, regulation of syncytium formation by plasma membranes, cytokine-cytokine receptor interaction, and deubiquitination. Only the regulation of syncytium formation by plasma membranes belongs to GO biological processes, while the others belong to the KEGG pathway category. The top 20 clusters of significantly enriched terms for the upregulated DE-genes are shown in Figure 2A1 and A2. The top 5 pathway and process enrichment analyses for predicted target genes of the top 10 downregulated DE-miRNAs included covalent chromatin modification, pathways in cancer, intracellular signalling by second messengers, regulation of cellular amide metabolic process, and protein deacetylation. Except for pathways in cancer (KEGG pathway category) and intracellular signalling by second messengers (Reactome Gene Sets category), all above results belong to the GO biological processes category. The top 20 clusters of significantly enriched terms for predicted target genes of the top 10 downregulated DE-miRNAs are shown in Figure 2B1 and B2. The top 5 pathway and process enrichment analyses for the downregulated DE-genes included divalent metal ion transport, intracellular mRNA localization, hair follicle maturation, stem cell division, and the ERK1 and ERK1 cascades. They all belong to the GO biological processes category. There were only 17 clusters of significantly enriched terms for downregulated DE-genes; therefore, only 17 clusters are shown in Figure 3A1 and A2. The top 5 pathway and process enrichment analyses for 10 predicted target genes of upregulated DE-miRNAs included heart development, blood vessel development, muscle structure development, pathways in cancer, and response to growth factor. Only pathways in cancer belong to the KEGG pathway category, while others belong to GO biological processes. The top 20 clusters of significantly enriched terms for predicted target genes of the 10 upregulated DE-miRNAs are shown in Figure 3B1 and B2.

Discussion
The treatment of patients with NKTL has improved over the past decades. Gene therapy is increasingly becoming an effective way to treat cancer based on the identification of DE-genes in human disease. Abnormal expression of miRNAs is involved in regulating all aspects of the occurrence, invasion and metastasis of cancers, which play an important role in cancer diagnosis and treatment [6]. However, until now, DE-genes and DE-miRNAs focused on nasal-type NKTL compared with normal nasal NKTL have rarely been reported. Fortunately, the integration of bioinformatics and microarray technology has enabled more effective ways to identify DE-genes and DE-miRNAs, which can help to discover potential effective targets for the diagnosis and treatment of nasal-type NKTL patients. It is necessary to further study the functions of miRNAs, predict their target genes and understand the interaction between miRNAs and target genes in nasal-type NKTL at present.
In the present study, we screened an eligible gene dataset and miRNA dataset in the GEO database and obtained DE-genes and DE-miRNAs after standardized processing. A total of 503 significantly DE-genes were discovered, of which 281 were upregulated and 222 were downregulated. Furthermore, 106 miRNAs were screened out to be significantly differentially expressed in human nasal-type NKTL samples when compared to human normal nasal samples, including 57 upregulated and 49 downregulated miRNAs. A total of 1207 miRNA targets were predicted for the top 10 upregulated DE-miRNAs, and 564 miRNA targets were identified for the top 10 downregulated DE-miRNAs. At last, the one-to-one relationship between DE-genes and DE-miRNAs were obtained. Through enrichment analysis, we found that meiotic recombination, regulation of syncytium formation by plasma membrane fusion, deubiquitination, cytokine-cytokine receptor interaction, divalent metal ion transport, intracellular mRNA localization, hair follicle maturation, ERK1 and ERK2 cascade, and stem cell division were more likely to be involved in the development of nasal-type NKTL. Meiotic recombination events have been implicated in cancer development, including ovarian carcinomas, lung adenocarcinoma, head and neck squamous cell carcinoma, liver cancer, bladder urothelial carcinoma, and breast cancer, but meiotic recombination events in nasal-type NKTL have not yet been reported [13]. For nasal-type NKTL, regulation of syncytium formation by plasma membrane fusion is thought to be related to the Epstein-Barr virus (EBV) entry mechanism, which is useful for the development of infection inhibitors and developing EBV vaccine approaches [14]. Metalloprotease deubiquitinases can regulate the ubiquitin signalling cascade, which are thought to be new regulators in cancer progression, such as non-small cell lung cancer [15], but their involvement in nasal-type NKTL has not been previously reported.
In our study, six upregulated and five downregulated DE-genes and their upstream miRNAs were identified by analysing the one-to-one relationship between DE-genes and DE-miRNAs. CDC27-miR-548c-3p, FREM2-miR-373 * and ARHGAP29-miR-548c-3p were the top 3 upregulated genes and their upstream miRNAs in the nasal-type NKTL samples compared with normal nasal tissues. CDC27 (also known as ANAPC3) was recently reported to be related to tumour progression and poor survival in T-cell lymphoblastic lymphoma [16]. Additionally, FREM2 has been found to play an important role in the malignant progression of prostate adenocarcinoma [17], and it is also thought to be a biomarker of glioblastoma stem cells [18]. It has been demonstrated that miR-373 * can inhibit cell proliferation in colon cancer [19], and downregulation of miR-373 * was confirmed in childhood B-cell precursor acute lymphoblastic leukaemia [20]. ARHGAP29 is involved in the development and migration of tumors. Moreover, Margit Schraders et al. demonstrated that ARHGAP29 may be involved in mantle cell lymphoma development and maintenance through signal transduction [21]. There are few reports focused on CDC27, FREM2, ARHGAP29, QSER1, CD3EAP, SF3A1, miR-548c-3p, miR-373 * , miR-548c-3p, miR-548c-3p, miR-149 * and miR-548c-3p in nasal-type NKTL, and we believe these are worthy of further research.
In the present study, five downregulated DE-genes and their upstream miRNAs were identified (Table 2). AQP4 (miR-29b), ZFP36L2 (miR-142-3p) and SRP72 (miR-16) were the top 3 downregulated genes (upstream miRNAs) in nasal-type NKTL samples compared with normal nasal tissues. MicroRNA-29b is downregulated in cerebral ischaemia and participates in the occurrence and development of cerebral ischaemia by targeting AQP4 regulation [22]. Furthermore, AQP4 induces paraneoplastic syndrome by antibody production, which may be the mechanism of acute myelogenous leukaemia coexisting with optic neuritis after transplantation [23]. At the same time, AQP4 has been demonstrated to be involved in cell proliferation, invasiveness, migration, and apoptotic potential in glioma cells [24]. miR-29b has been reported to promote the occurrence and development of acute myeloid leukaemia [25] by participating in the regulation of tumour cell proliferation and apoptosis. Repression of miR-29b participates in regulation in cutaneous T-cell lymphoma through the miR-29b-BRD4 axis [26]. Moreover, it was demonstrated that miR-29s could alter the expression of targets involved in critical cancer pathways in Burkitt's lymphoma, such as cell cycle control, apoptosis inhibition and DNA methylation [27]. Moreover, T lymphoma invasion and metastasis 1 (Tiam1) was identified as a direct target of miR-29b [28]. However, so far, we have not found reports on miR-29b focused on nasal-type NKTL, which may become the innovation point and hot spot of future research. ZFP36L2 (also known as TIS11D, one of the NOTCH RNA repressors), dysfunction is associated with the pathogenesis of certain types of leukaemia, such as acute myeloid leukaemia and Burkitt's lymphoma [29]. In addition, ZFP36L2 plays a role in the development of acute lymphoblastic leukaemia [30]. However, there is little literature focused on the relationship between these genes and microRNAs (AQP4, ZFP36L2, SRP72, TSC22D2, DOCK5, miR-29b, miR-142-3p, miR-16, and let-7f) in nasal-type NKTL patients, which may deserve more in-depth research.
In our study, CDC27, ARHGAP29, QSER1, SF3A1 (miR-548c-3p), FREM2-miR-373 * , CD3EAP-miR-149 * , AQP4-miR-29b, ZFP36L2-miR-142-3p, SRP72, TSC22D2, DOCK5 (miR-16), and TSC22D2-let-7f were identified between nasal-type NKTL and normal nasal samples. Among them, eleven differentially expressed genes (CDC27, FREM2, ARHGAP29, QSER1, CD3EAP, SF3A1, AQP4, ZFP36L2, SRP72, TSC22D2, DOCK5) and six differentially expressed miRNAs (miR-548c-3p, miR-373 * , miR-149 * , miR-29b, miR-142-3p, let-7f) have been confirmed to be associated with other types of lymphoma, other cancers, or other non-neoplastic diseases. However, few DE-genes and DE-miRNAs focused on nasal-type NKTL have been explored, which may be an innovation of future research. The DE-genes and DE-miRNAs can serve as biomarkers in nasal-type NKTL, which is helpful to understand the pathogenesis of nasal-type NKTL patients. In addition, it may provide strong theoretical support for the early diagnosis of nasal-type NKTL patients and the development of novel therapeutic targets. However, it still has limitations: 1). Target gene prediction was only carried out for the top 10 DE-miRNAs; 2). Only the target genes of the top 10 DE-miRNAs were selected for further enrichment analysis; 3). Since there have been few studies focused on nasal-type NKTL, no dataset with a larger sample size can be found, so, 2 human normal nasal samples and 17 human nasal-type NKTL samples from the GSE313772 dataset were chosen in our study. 4). Due to a lack of experimental verification, more research should be carried out.