Integrated Bioinformatic Analysis of Competing Endogenous RNA Network of Choriocarcinoma.

BACKGROUND Numerous studies have demonstrated that noncoding RNAs are involved in choriocarcinoma (CC). The competing endogenous RNA (ceRNA) network plays an important role in the occurrence and development of carcinoma. However, the involvement of the ceRNA network in CC remains unclear. The current study aimed to investigate the regulatory mechanism of ceRNA in CC. MATERIAL AND METHODS We downloaded the messenger RNAs (mRNAs) expression profiles (GSE20510 and GSE65654) and microRNAs (miRNAs) expression profiles (GSE32346 and GSE130489) from GEO datasets. The limma package of R software was used to identify differentially expressed RNAs (DERNAs). Then, we performed functional annotation of the differentially expressed mRNAs (DEmRNAs). TargetScan, miRDB, miRWalk, and Starbase were used to construct a CC-specific ceRNA network and select key molecules. RESULTS The results identified a total of 177 DEmRNAs and 189 differentially expressed miRNAs (DEmiRNAs) between the trophoblast and CC cell line samples. Ten differentially expressed lncRNAs (DElncRNAs) were obtained based on experimental studies. The DEmRNAs were mainly enriched in cell proliferation, positive regulation of the apoptotic process, and cell death. A total of 10 genes were ascertained as hub genes. Based on DEmRNAs, DEmiRNAs, and DElncRNAs, a CC-specific ceRNA network was established. Five DElncRNAs, 15 DEmiRNAs, and 45 DEmRNAs were identified. In addition, LINC00261, MEG3, MALAT1, H19, and OGFRP1 were identified as 5 key lncRNAs in choriocarcinoma. CONCLUSIONS This study provides novel insights into CC mechanisms and identified potential therapeutic targets for CC.

Background

Choriocarcinoma (CC) is an aggressive, malignant trophoblastic neoplasm, either gestational or nongestational in origin [1]. Although CC is sensitive to chemotherapy, about 20% of patients with drug-resistant CC experience treatment failure and are at increased risk of metastatic lesions, hysterectomy, and decreased fertility [2,3]. Early metastasis and drug resistance in chemotherapy remain a clinical challenge in CC. Therefore, potential biomarkers and therapeutic targets in CC individualized treatment strategies need to be investigated.

Historically, most of the non-protein-coding parts of the human genome have been considered junk DNA. Although less than 2% of the human genome encodes proteins, the majority of all nucleotides can be detected under certain conditions [4]. MicroRNAs (miRNAs) are small endogenous RNAs with a regulatory function, consisting of 20–25 nucleotides. miRNA research in cancer is expanding because miRNAs are promising candidates for cancer biomarker development [5]. Functional studies have demonstrated that miRNA disorders are causal in many cancer cases, and miRNA mimics and miRNA inhibitors have shown promise as novel therapeutic agents [6]. Long noncoding RNAs (lncRNAs), a class of novel noncoding RNAs, are defined as transcripts of more than 200 nucleotides [7]. lncRNA-dysregulated expression is widespread in the pathogenesis of various human diseases, including cancer [8]. Studying lncRNA functions will provide researchers with tools and opportunities to develop lncRNA-based therapeutics for cancer [9].

Salmena and colleagues first proposed the concept of competitive endogenous RNA (ceRNA). Evidence has shown that pseudogene RNAs can act as a ‘sponge’ by competitively binding common miRNAs, releasing or weakening inhibitory effects through isolating miRNAs away from parental mRNA [10]. Accumulating evidence revealed the regulatory network of miRNAs and their target genes in carcinoma. Studies have shown that LINC01133 can competitively bind to miR-106a-3p and regulate adenomatous polyposis coli expression in gastric cancer [11]. Besides acting as a ceRNA sponge to absorb miR-570, LINC00612 can elevate the expression of PHD finger protein 14 in bladder cancer [12]. Previous studies suggested that lncRNA may act as ceRNA in CC. For example, lncRNA MALAT1 promoted the proliferation of CC by miR-218-mediated Fbxw8 regulation [13]. lncRNA RNA PCA3 contributes to the progression of CC by acting as a ceRNA against miR-106b [14]. A recent study reported targeted regulating of miR-515-5p by lncRNA LOXL1-AS1 in the proliferation and migration of CC [15].

How the ceRNA network is regulated in CC remains unclear. In the present study, we built a CC-specific ceRNA regulatory network using the RNA expression profiles collected from Gene Expression Omnibus (GEO) datasets. The research aimed to initially discover lncRNA-miRNA-mRNA ceRNA-mediated regulatory mechanisms in CC. This study may provide novel potential therapeutic targets for CC.

Material and Methods

Selection of the datasets

We selected datasets including trophoblast cells and CC cell lines from the GEO database (http://www.ncbi.nlm.nih.gov/geo/) of microarray data on gene expression obtained from mRNA and miRNA. Two expression profiling datasets, GSE20510 and GSE65654, which identify mRNA involved in CC, were obtained from the GEO (GPL96 platform, [HG-U133A] Affymetrix Human Genome U133A Array, and GPL 10558 platform, Illumina HumanHT-12 V4.0 expression beadchip). Trophoblast cells included first-trimester cytotrophoblasts and immortalized extravillous trophoblast cell line HTR-8, as well as CC cell lines including JEG-3 and BeWo. The GSE20510 dataset contained 3 trophoblast samples and 7 CC cell line samples. The GSE65654 dataset contained 2 trophoblast samples and 4 CC cell line samples. We also downloaded the expression profiles of miRNA GSE32346 (2 trophoblast samples and 4 CC cell line samples) and GSE130489 (3 trophoblast samples and 3 CC cell line samples).

Data Quality Control

Intra-group data repeatability in the was assessed by Pearson’s correlation test. We used the software and operating environment in the R programming language for drawing graphs and statistical analysis. The heat map in R software was used to visualize correlations between all samples. Principal component analysis (PCA) is an unsupervised statistical analysis technique and is a widely used sample clustering method. It is an efficient information compression method, which can reduce the dimensionality in the dataset and discover latent variables among samples. Sample clustering analysis was performed to assess the intra-group data repeatability.

Screening DEmRNAs and DEmiRNAs

We combined trophoblast samples and CC cell line samples data. The “Limma” package in R software was used to screen out DEmRNAs and DEmiRNAs. DEmRNAs were identified if | Log2 Fold Change (FC)| >1.5 and P value <0.05. Meanwhile, DEmiRNAs were defined if |Log2 FC| >2 and P value <0.05. Volcano maps were drawn using “ggplot2 package” operated by R software. The Circos was conducted by Omistudio, an online tool (http://www.Omistudio.cn/). The Circos delineated expression level and the overlaps between differently expressed RNA data visualized in the heatmap simultaneously. The Online Venn tool (http://bioinformatics.psb.ugent.be/webtools/Venn/) was used to create Venn diagrams. Venn diagrams exhibit common DEmRNAs and DEmiRNAs shared between the 2 datasets.

Functional Enrichment of DEmRNAs

The online tool DAVID (https://david.ncifcrf.gov/home.jsp/) (version 6.8) was used for enrichment analysis. DAVID was used to reveal Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis enriched by DEmRNAs. GO is an ontology database widely used in bioinformatics, including cellular components, molecular functions, and biological processes. KEGG is a database resource for understanding high-level functions and utilities of the biological system. The “ggplot2 package” in R software was used to visualize the enrichment analysis. Cytoscape plug-ins ClueGO and CluePedia were used to perform enrichment analysis for DEmRNAs again. The hypergeometric significance P value<0.05 was regarded as significant.

Protein–Protein Interaction (PPI) Network Construction

Search Tool for the Retrieval of Interacting Genes (STRING; http://string-db.org/) was used to map the common DEmRNAs. The PPI network was constructed by Cytoscape 3.7.1. The app Molecular Complex Detection (MCODE) from the Cytoscape software suite was used to perform the gene network clustering analysis. MCODE identified the most significant PPI network modules; the submodule of PPI network was a degree MCODE scores >5, cut-off=2, Maximum depth=100, node score cut-off=0.2, and k-score=2. In CytoHubba, a plug-in of Cytoscape, the key genes were explored under the degrees set (degrees ≥10). The expression level of hub genes was drawn using the “pheatmap package” in R.

Construction of a CC-Specific ceRNA Network

To improve the accuracy of DERNAs, we excluded the DEmRNAs and DEmiRNAs that were not paired in miRDB (http://www.mirdb.org/), TargetScan (http://www.targetscan.org/), and miRWalk (http://mirwalk.umm.uni-heidelberg.de/). The common DERNAs were collected to construct a choriocarcinoma-specific ceRNA network. The prediction of lncRNA-miRNA interaction was investigated by Starbase (http://starbase.sysu.edu.cn/). We searched PubMed, Web of science, and Medline for DElncRNAs in CC based on experimental research. Finally, Cytoscape software was used to visualize the CC-specific ceRNA network.

Statistical Analysis

The statistical significance of differences between the 2 groups was assessed using the t test. All analyses were conducted using R software 4.0.3. A P value <0.05 was considered statistically significant.

Results

Quality of the Datasets

Principal component analysis (PCA) and Pearson’s correlation test were used to test the intra-group data repeatability. As shown in Figure 1A1, 1A2, the samples in the trophoblast group had strong correlations in the GSE20510 and GSE65654 datasets, and the samples in the CC group were similar to the trophoblast group. Similarly, PCA results showed that the intra-group data repeatability for GSE20510 and GSE65654 was acceptable, and the trophoblast group and the CC group were separated from each other (Figure 1B1, 1B2). In addition, the samples in the trophoblast group had strong correlations in the GSE32346 and GSE130489, and the samples in the CC group were similar to the trophoblast group (Figure 1A3, 1A4). Based on the PCA, the intra-group data repeatability for GSE32346 and GSE130489 was acceptable, and the trophoblast group and the CC group were separated from each other (Figure 1B3, 1B4).

Figure 1

(A) A1–A4 represents Pearson’s correlation analysis of GSE20510, GSE65654, GSE32346, and GSE130489 datasets, respectively. The color reflects the intensity of the correlation. The greater the absolute value of a number, the higher the correlation. (B1–B4) PCA of samples from the GSE20510, GSE65654, GSE32346, and GSE130489 datasets. In the figure, the x-axis is principal component 1 (PC1) and the y-axis is principal component 2 (PC2). (C1–C4) The volcano plot of the GSE20510, GSE65654, GSE32346, and GSE130489 datasets. (D1) Circos of GSE20510 and GSE65654 was performed using the online tool Omistudio. The outer circle represents gene lists, the inner circle represents the gene expression level, and the red curve represents the overlap between differentially expressed gene lists. (D2) Circos of GSE32346 and GSE130489. (E1, E2) 177 sharing genes between the GSE20510 and GSE65654 datasets are shown in Venn diagram, and 189 differentially expressed miRNAs were contained in the GSE32346 and GSE130489 simultaneously.

DEmRNAs and DEmiRNAs Identified Between Trophoblast and Choriocarcinoma

After the analyses of GSE20510, GSE65654, GSE32346, and GSE130489 datasets with the “Limma package” in R software, the DEmRNA between trophoblast cells and CC cell lines were presented in volcano plots (Figure 1C1, 1C2). Meanwhile, the DEmiRNAs were shown in volcano plots (Figure 1C3, 1C4). The overlap between differentially expressed mRNA datasets of GSE20510 and GSE65654 and the expression level visualized in the heatmap simultaneously were presented in Circos (Figure 1D1). The DEmiRNA lists of GSE32346 and GSE130489 were also performed in Circos (Figure 1D2). A total of 1205 DEmRNAs in the GSE20510 dataset, 854 DEmRNAs in the GSE65654 dataset, 316 DEmiRNAs in the GSE32346, and 324 DEmiRNAs in the GSE130489 were identified. Venn diagrams were constructed and revealed that 177 DEmRNAs (Supplementary Table 1) and 189 DEmiRNAs (Supplementary Table 2) were contained within the 2 examined datasets simultaneously (Figure 1E1, 1E2).

GO and KEGG Enrichment Analysis

Functional enrichment analysis was conducted to explore the biological functions of DEmRNAs. The top 5 biological processes (BP) terms were enriched in positive regulation of transcription from RNA polymerase II promoter, cell proliferation, heart development, angiogenesis, and positive regulation of the apoptotic process (Figure 2A, Table 1). The top 5 cell components (CC) terms of DEmRNAs were plasma membrane, extracellular exosome, cytosol, extracellular region, and extracellular space (Figure 2B). The top 3 molecular function (MF) terms were protein binding, sequence-specific DNA binding, and protein homodimerization activity (Figure 2C). The KEGG pathway analysis disclosed that the main enriched pathways included pathways in cancer, signaling pathways regulating pluripotency of stem cells, Wnt signaling pathway, proteoglycans in cancer, and Rap1 signaling pathway (Figure 2D). We used Cytoscape plug-ins ClueGO and CluePedia to analyze DEmRNAs and found that positive regulation of stress-activated MAPK cascade, interleukin-1 beta production, response to hydrogen peroxide, and connective tissue development were involved (Figure 2E).

Figure 2

(A) The top 14 changes in DEmRNAs’ biological processes (BP) in trophoblast and choriocarcinoma samples through GO enrichment analyses. (B) Cellular components (CC) analysis of DEmRNAs. (C) Molecular functions (MF) analysis of DEmRNAs. (D) The KEGG pathway analysis of DEmRNAs. (E) Using Cytoscape plug-in ClueGO and CluePedia to display the network of enriched terms colored by cluster identity, where nodes that share the same cluster identity are typically close to each other. DEmRNAs – differentially expressed mRNAs; GO – Gene Ontology; KEGG – Kyoto Encyclopedia of Genes and Genomes.

Table 1

The enriched GO and KEGG term.

ID	Term	Count	P value
BP
GO: 0007507	Heart development	11	1.27E-05
GO: 0001525	Angiogenesis	10	3.49E-04
GO: 0030324	Lung development	6	8.83E-04
GO: 0008283	Cell proliferation	12	9.56E-04
GO: 0043065	Positive regulation of apoptotic process	10	0.00279
GO: 0060068	Vagina development	3	0.00326
GO: 0046330	Positive regulation of jnk cascade	5	0.00374
GO: 0048048	Embryonic eye morphogenesis	3	0.00405
GO: 0032940	Secretion by cell	3	0.00405
GO: 0009749	Response to glucose	5	0.00440
CC
GO: 0070062	Extracellular exosome	51	1.48E-06
GO: 0031012	Extracellular matrix	13	1.89E-05
GO: 0005578	Proteinaceous extracellular matrix	12	3.77E-05
GO: 0009986	Cell surface	15	4.79E-04
GO: 1903561	Extracellular vesicle	5	0.00113
MF
GO: 0005515	Protein binding	113	8.02E-05
GO: 0005520	Insulin-like growth factor binding	4	7.08E-04
GO: 0050840	Extracellular matrix binding	4	0.00213
GO: 0005201	Extracellular matrix structural constituent	5	0.00446
GO: 0004714	Transmembrane receptor protein tyrosine kinase activity	4	0.00633
KEGG
hsa04550	Signaling pathways regulating pluripotency of stem cells	7	0.00493
hsa04115	p53 signaling pathway	4	0.04003
hsa05169	Epstein-Barr virus infection	5	0.04925
hsa04350	TGF-beta signaling pathway	4	0.06958
hsa04310	Wnt signaling pathway	5	0.07099
hsa04960	Aldosterone-regulated sodium reabsorption	3	0.07184
hsa05200	Pathways in cancer	9	0.07585
hsa05205	Proteoglycans in cancer	6	0.07597
hsa04978	Mineral absorption	3	0.0885
hsa04015	Rap1 signaling pathway	6	0.08924

PPI Network Construction and Module Analysis

A PPI network with 127 nodes and 280 edges was acquired; the network had an interaction score >0.4 according to the STRING online database (Figure 3A). The nodes represent genes, and the edges correspond to the links between genes. Blue represents downregulated genes and red represents upregulated genes. We used MCODE and CytoHubba in Cytoscape to perform network gene clustering to identify the critical module and hub genes. As shown in Figures 3B, 4 key modules with 9 upregulated genes and 15 downregulated genes were identified. The list of key modules genes is shown in Table 2. Furthermore, 10 genes (CYR61, FSTL1, CTGF, FBN1, APOA2, MXRA8, CKAP4, MGAT4A, IGFBP3, and IGFSP7) were determined as hub genes with the greatest degree of network connection (Figure 4A). Hierarchical cluster analysis revealed that the hub genes could effectively differentiate the CC cells from the trophoblast cells (Figure 4B). Among these genes, MGAT4A, APOA2, and MXRA8 were highly expressed in the CC cell line samples, while FSTL1, IGFBP3, IGFBP7, CYR61, FBN1, CKAP4, and CTGF were downregulated in the CC cell line samples. The hub genes may play a vital role in the occurrence or development of CC.

Figure 3

(A) The protein–protein interaction (PPI) analysis in choriocarcinoma. The blue and red nodes represent downregulated and upregulated genes, respectively. The degree value is displayed by the size of the nodes. (B) The most significant modules in the PPI network.

Figure 4

(A) Cytohubba a plug-in of Cytoscape identified the hub genes from the PPI network. (B) Hierarchical cluster analysis suggested that the hub genes could effectively differentiate the choriocarcinoma samples from the trophoblast samples in the GSE20510 and GSE65654 datasets. Blue represents downregulated genes, red represents upregulated genes.

Table 2

The 9 upregulated and 15 downregulated key module genes identified by MCODE.

Genes	logFC	P value	adj. P value
MAGT4A	3.2192	1.23E-08	1.43E-07
MXRA8	3.5007	3.55E-08	3.07E-07
APOA2	2.6460	4.83E-10	2.12E-08
CCR7	2.4503	2.63E-05	7.39E-05
ACKR3	1.9622	7.69E-07	3.63E-06
GPC4	2.8810	5.55E-06	1.94E-05
MAGEA12	5.3826	1.69E-15	8.25E-12
MAGEA1	2.2676	5.76E-10	2.24E-08
MAGEA6	5.8546	2.66E-15	8.25E-12
FSTL1	−6.1075	5.19E-15	1.24E-11
IGFBP3	−2.8418	1.79E-09	4.10E-08
IGFBP7	−2.3101	2.83E-09	5.32E-08
CYR61	−2.3364	6.88E-06	2.32E-05
FBN1	−2.8740	5.25E-11	7.41E-09
CKAP4	−1.5676	0.001285	0.002333
NMU	−2.1097	4.14E-06	1.50E-05
GNAI2	−2.3328	1.01E-07	6.85E-07
GNG11	−3.0939	1.22E-13	1.46E-10
WNT5A	−1.9255	2.87E-06	1.10E-05
DKK1	−3.5953	2.28E-08	2.19E-07
JUN	−1.8796	4.73E-08	3.84E-07
CTGF	−5.4197	8.13E-12	2.17E-09
SRGN	−4.0831	1.46E-11	3.15E-09
BGN	−1.7001	1.02E-07	6.88E-07

CC-Specific ceRNA Network

Among the 15 obtained DEmiRNAs (Table 3), miR-10b, miR-155, miR-598, miR-411, miR-370, miR-134, miR-539, miR-10a, miR-376c, miR-495, and miR-152 served as the key upregulation miRNAs; miR-372, miR-520f, miR-934, and miR-1270 served as the key downregulation miRNAs. Due to the lack of high-throughput sequencing lncRNA in CC, we searched the database for DElncRNAs based on experimental research. Five DElncRNAs (LINC00261, MEG3, MALAT1, H19, and OGFRP1) were considered key lncRNAs regulating their target mRNAs’ expression by competing with miRNAs in CC. In total, 5 DElncRNAs (4 upregulated and 1 downregulated), 15 interacting DEmiRNAs (11 upregulated and 4 downregulated), and 45 DEmRNAs (18 upregulated and 27 downregulated) participated in the formation of CC-specific ceRNA regulatory network (Figure 5).

Table 3

DEmiRNAs interact with the DEmRNAs retrieved from miRDB, TargetScan and miRWalk database.

miRNAs	mRNAs
has-miR-10b	CTPS2
has-miR-155	DLX6, NR6A1, ANKRD6, TCF4, BDNF
has-miR-598	LIN28A, ID1, LBH
has-miR-411	PHKB, ACVR2B, MGAT4A
has-miR-370	NR6A1, RAB3B
has-miR-134	NR6A1, ACVR2B, FAT2, RFTN1, BDNF
has-miR-539	NR3C1, SERPINE2, DA2AP2, NR6A1, ICK, RAB3B, DNMT3B, WNT5A, FAM111A, NFIB
has-miR-10a	BDNF
has-miR-376c	NFIB, ATP10D
has-miR-495	ACVR2B, TCF4, SLC26A2
has-miR-152	TXNIP, CTPS2, ATP10D, NR3C1, CCNG1, MET, ELMO1, ATP2B4, DDAH1, TXNIP, ALCAM, TCF4, FBN1, MGAT4A
has-miR-372	TXNIP, EPHA2, RAB3B, NFIB, ANKRD6, TFPI2, MGAT4A, THAP10, MAGEA6, ALDH1A3
has-miR-520f	ASAP2, WNT5A, ATP10D, GNAI2, RND3, CCKBR, TSPAN12, NFIB
has-miR-934	RAB3B
has-miR-1270	ATP10D, NR3C1, TPD52L1, CCR7

Figure 5

The lncRNA-miRNA-mRNA ceRNA network constructed by Cytoscape. The diamond, square, and circle represent lncRNA, miRNA, and mRNA, respectively. Blue stands for upregulated RNAs, red stands for downregulated RNAs. lncRNAs – long noncoding RNAs; miRNAs – microRNAs; mRNA – messenger RNAs; ceRNA – competing endogenous RNA.

Discussion

Choriocarcinoma is a malignant tumor that occurs frequently in women of childbearing age. Emerging studies have indicated the regulatory roles of miRNAs and lncRNAs in the occurrence and development of CC. CeRNAs are transcripts that cross-regulate each other through competition with shared miRNAs. According to the ceRNA hypothesis, lncRNAs could serve as a sponge of target miRNAs [16]. CeRNA provides a novel way to explore cancer diagnosis, prognosis, and targeted treatment, but the ceRNA regulatory mechanism in CC has been unclear. Therefore, we constructed a CC-specific ceRNA network to investigate the regulatory mechanisms in CC.

In this research, we analyzed the differentially expressed RNAs to construct a ceRNA network. Herein, a train of bioinformatic analyses was performed on 4 independent high-throughput datasets (from trophoblast and CC cell line). By the analysis of GSE20510 and GSE65654 datasets, 177 common DEmRNAs were identified, of which 82 and 95 were upregulated and downregulated, respectively. It is well known that hub genes play critical roles in biological networks. Among the DEmRNAs, 21 hub nodes were identified by MCODE and 10 potential hub genes (CYR61, FSTL1, CTGF, FBN1, APOA2, MXRA8, CKAP4, MGAT4A, IGFBP3, and IGFSP7) were obtained using the CytoHubba plug-in of Cytoscape. Alpha-1,3-mannosyl-glycoprotein 4-beta-N-acetylglucosaminyltransferase A (MGAT4A) is a key enzyme in the pathway of synthesizing complex N-glycans. Several previous studies showed that MGAT4A is involved in the progression of choriocarcinoma, hepatocarcinoma, and colorectal carcinoma [17-19]. Chemokine receptors belong to the G-protein-coupled receptors (GPCRs) family, which appears to be involved in inflammatory diseases, tumor growth, and metastasis. A positive correlation between chemokine receptor expression and worse prognosis has been found in most but not all cancers. CC-chemokine receptor 7 (CCR7) is significantly overexpressed in colorectal tumor tissues compared with paired normal tissues. Significant differences were also found in overall survival and disease-free survival [20]. In breast cancer, it has been reported that the chemokine receptors CCR7 and their ligands are associated with metastasis [21]. This is consistent with our finding that CCR7 is a key node and has relatively high expression in CC.

The GO analyses showed that the DEmRNAs are mainly enriched in cell proliferation, positive regulation of the apoptotic process, and cell death. The KEGG pathway analysis also showed that the DEmRNAs mostly participate in pathways in cancer, signaling pathways regulating pluripotency of stem cells, Rap1 signaling pathway, and Wnt signaling pathway. The Wnt signaling pathway promotes decidualization, endometrial function, and trophoblast differentiation, and changes in Wnt signaling components were found in cancers of reproductive tissues, endometriosis, and gestational diseases [23]. Wnt/β-catenin signaling regulates StarD7 transcription in CC cell line JEG-3 cells through the canonical pathway [23].

Noncoding RNAs (ncRNAs) are of interest because of their versatile roles in biological and pathological processes, particularly the widely studied small miRNAs and long noncoding RNAs (lncRNAs) [24]. Through the analysis of GSE32346 and GSE130489 datasets, 189 common DEmiRNAs were obtained, of which 147 were downregulated and 42 were upregulated. Fifteen DEmiRNAs (11 upregulated and 4 downregulated) serve as key target miRNAs in the ceRNA network. It has been reported that miR-1270 is a negative regulator in liver carcinoma and is involved in post-transcriptional regulation of centromere protein M, and hepatitis B virus X protein can promote liver carcinoma by inhibiting miR1270 [25]. miR-1270 is aberrantly upregulated in papillary thyroid cancer. Reduced miR-1270 expression suppresses cancer cell proliferation and migration, which was reversed by suppression of cancer cell invasion downregulation [26]. miR-411 acts as an oncogene in osteosarcoma and esophageal squamous cell carcinoma [27,28]. On the contrary, miR-411 inhibits the growth and metastasis of bladder cancer, cervical cancer, and gastric cancer [29-31].

Five DElncRNAs, 15 DEmiRNAs, and 45 DEmRNAs were identified as abnormally expressed RNAs and we built a CC-specific network. This study found MALAT1 and LINC00261 can increase CCR7 via competing with miR-1270. The metastasis-associated lung adenocarcinoma transcript 1(MALAT1) was originally identified as a prognostic parameter for lung adenocarcinoma or squamous cell carcinoma [32]. A series of data has been collected linking MALAT1 with cancer types and provided new insights into its biogenesis, and molecular functions [33]. Several studies suggested that MALAT1 can be aberrantly expressed in choriocarcinoma, ovarian, cervical, colorectal cancers, and breast cancer [13,34-37]. MALAT1 may be an oncogene that promotes CC proliferation and tumor growth through miR-218-mediated F-box/WD repeat-containing protein 8 regulation [13]. Low expression of lncRNA LINC00261 was reported in multiple cancers, including breast, colon, and lung cancer. These results indicate LINC00261 is a tumor suppressor that halts cellular proliferation, invasion, and tumor growth [38-40]. LINC00261 was significantly underexpressed in CC tissues and CC cell lines. Overexpression of LINC00261 resulted in decreased cell proliferation and cell cycle arrest in G0/G1 phase, suppressed cell migration and invasion, and promoted cell apoptosis and the relative activities of caspase 3 and caspase 9 in CC cells [41]. We speculated that upregulating MALAT1 indirectly suppresses the expression of TCF4 via competitively binding to miR-155. miR-155 is involved in the syncytialization of trophoblast cells induced by lipopolysaccharides [42]. It has been reported that LIN28A has a functional role in regulating trophoblast differentiation and function, and that loss of LIN28A in human trophoblast is sufficient to induce differentiation [43]. In the CC ceRNA network, LIN28A was upregulated in CC and might be regulated by miR-598.

Conclusions

In the present study, we constructed a CC-specific ceRNA regulatory network by bioinformatics analysis. A CC-specific ceRNA regulatory network might provide novel insights into CC mechanisms and identify potential therapeutic targets for CC.

Supplementary Table 1

The list of DEmRNAs, including 82 upregulated DEmRNAs and 95 downregulated DEmRNAs.

Regulation	mRNAs
Upregulated (n=82)	C1orf106, PTN, NDUFAF4, AKR1B10, RGS16, YBX2, MXRA8, ANKRD6, JADE3, F2RL1, SLC2A3, SH2B2, DENND1A, DNMT3B, FYB, FAT2, CTSV, SARS, CLDN6, AGTPBP1, SNX10, ATP1A1, DLX6, ACTR6, CRYBA2, BNIP3, APOA2, IMPA2, CTH, NR6A1, SFN, CCKBR, TPD52L1, TYRO3, REEP1, CCR7, EPCAM, ZNF215, ICK, RFK, COL9A3, GPC4, SLC39A4, DNMT3L, TSPAN12, IRX4, KHDC1L, NLRP2, FZD5, HAND1, DDX3Y, CCNG1, ELOVL6, KHDRBS3, FBLN2, MAGEA6, TMEM38B, ELMO1, MYCN, ZNF415, LHFP, ATP2B4, NMRK2, C1orf105, MGAT4A, MAGEA1, CYB5R2, ACVR2B, CRIP1, ADAM8, SLC29A2, GATA2, MAGEA12, PCSK5, CBLC, LIN28A, HMGA1, TFF1, ZNF280A, PHKB, ACKR3, ANKZF1
Downregulated (n=95)	FSTL1, CD44, FAM111A, HSPB1, BID, BTBD3, ID3, MET, MARCKS, AXL, MSN, GNAI2, SERPINE2, PAPSS2, BMPR1A, ATP9A, EPHA2, SEL1L3, TRAF5, THAP10, NFIB, PLP2, RFTN1, DAZAP2, TSPO, RAB3B, BDNF, DDAH1, MEIS2, ZNF175, ALCAM, ID1, EMP3, ARMCX1, DPYD, IGFBP7, VASH1, SLC26A2, TCF4, GNG11, DKK1, ADORA2B, NR3C1, IL7R, IFI16, OXCT1, ASAP2, CAP2, ALDH1A3, SH3BGRL3, MEST, KLF2, SGK1, COL5A2, VIM, FAT1, GSTO1, HLA-A, RHOBTB3, ATP10D, NYNRIN, C11orf95, FRZB, CTPS2, WNT5A, NMU, INPP1, NT5E, SCCPDH, IGFBP3, RND3, SRGN, PMP22, CDR2L, TMSB4X, S100A6, TFPI2, JUN, GLCE, MT2A, FBN1, LBH, TGFBI, HOXB2, BGN, CYR61, TXNIP, ZAK, CTGF, ECI2, ZNF22, UCP2, EFEMP1, CKAP4, CBR1

Supplementary Table 2

The list of DEmiRNAs, including 128 upregulated DEmiRNAs and 60 downregulated DEmiRNAs.

Regulation	miRNAs
Upregulated (n=128)	hsa-miR-376b, hsa-miR-146a, hsa-miR-147, hsa-miR-642, hsa-miR-873, hsa-miR-380, hsa-miR-125b, hsa-miR-542-3p, hsa-miR-886-5p, hsa-miR-138, hsa-miR-127-3p, hsa-miR-654-3p, hsa-miR-126, hsa-let-7d, hsa-miR-29a, hsa-miR-488, hsa-miR-192, hsa-let-7g, hsa-miR-124, hsa-miR-133a, hsa-miR-616, has-miR-155, hsa-miR-222, hsa-miR-432, hsa-miR-184, hsa-miR-125a-5p, hsa-miR-125a-3p, hsa-miR-487a, hsa-miR-337-5p, hsa-miR-221, hsa-miR-195, hsa-miR-487b, hsa-miR-28-3p, hsa-miR-198, hsa-miR-504, hsa-miR-214, hsa-miR-323-3p, hsa-miR-200a, hsa-miR-376c, hsa-miR-27a, hsa-miR-335, hsa-miR-152, hsa-miR-127-5p, hsa-miR-539, hsa-miR-496, hsa-miR-146b-5p, hsa-miR-99a, hsa-miR-133b, hsa-let-7b, hsa-miR-544, hsa-miR-143, hsa-miR-328, hsa-miR-139-5p, hsa-miR-320, hsa-let-7a, hsa-miR-889, hsa-miR-28-5p, hsa-miR-194, hsa-miR-628-5p, hsa-miR-217, hsa-miR-199b-5p, hsa-miR-342-5p, hsa-miR-493, hsa-miR-886-3p, hsa-miR-199a-3p, hsa-miR-134, hsa-miR-574-3p, hsa-miR-342-3p, hsa-miR-409-5p, hsa-miR-542-5p, hsa-miR-99b, hsa-let-7e, hsa-miR-370, hsa-miR-199a-5p, hsa-miR-758, hsa-miR-10b, hsa-miR-431, hsa-miR-216b, hsa-let-7c, hsa-miR-103, hsa-miR-129-3p, hsa-miR-24, hsa-miR-100, hsa-miR-379, hsa-miR-299-5p, hsa-miR-876-5p, hsa-miR-369-5p, hsa-miR-382, hsa-miR-654-5p, hsa-miR-186, hsa-miR-10a, hsa-miR-21, hsa-miR-598, hsa-miR-655, hsa-miR-302c, hsa-miR-137, hsa-miR-331-3p, hsa-miR-410, hsa-miR-329, hsa-miR-153, hsa-miR-494, hsa-miR-411, hsa-miR-485-3p, hsa-miR-324-3p, hsa-miR-876-3p, hsa-miR-376a, hsa-miR-363, hsa-miR-196b, hsa-miR-145, hsa-miR-874, hsa-miR-338-3p, hsa-miR-148b, hsa-miR-340, hsa-miR-140-3p, hsa-miR-147b, hsa-miR-497, hsa-miR-369-3p, hsa-let-7f, hsa-miR-218, hsa-miR-551b, hsa-miR-433, hsa-miR-154, hsa-miR-541, hsa-miR-142-3p, hsa-miR-495, hsa-miR-29c, hsa-miR-132, hsa-miR-29b
Downregulated (n=60)	hsa-miR-518d-5p, hsa-miR-520b, hsa-miR-520c-3p, hsa-miR-518e, hsa-miR-517b, hsa-miR-1276, hsa-miR-520g, hsa-miR-1301, hsa-miR-519b-3p, hsa-miR-1270, hsa-miR-769-3p, hsa-miR-512-3p, hsa-miR-483-3p, hsa-miR-519d, hsa-miR-1179, hsa-miR-520e, hsa-miR-93, hsa-miR-934, hsa-miR-523, hsa-let-7i, hsa-miR-1292, hsa-miR-520a-3p, hsa-miR-375, hsa-miR-31, hsa-miR-526b, hsa-miR-518c, hsa-miR-519c-3p, hsa-miR-29b-2, hsa-miR-18a, hsa-miR-769-5p, hsa-miR-516-3p, hsa-miR-92a-1, hsa-miR-519a, hsa-miR-519e, hsa-miR-520h, hsa-miR-522, hsa-miR-518f, hsa-miR-148a, hsa-miR-24-2, hsa-miR-524-5p, hsa-miR-29b-1, hsa-miR-516a-5p, hsa-miR-213, hsa-miR-520a-5p, hsa-miR-517, hsa-miR-641, hsa-miR-181a-2, hsa-miR-520f, hsa-miR-151-3p, hsa-miR-373, hsa-miR-944, hsa-miR-371-3p, hsa-miR-518b, hsa-miR-524, hsa-miR-517c, hsa-miR-372, hsa-miR-15b, hsa-miR-151-5P, hsa-miR-520d-5p, hsa-miR-744, hsa-miR-512-5p