Bioinformatics Analysis of Choriocarcinoma-Related MicroRNA-Transcription Factor-Target Gene Regulatory Networks and Validation of Key miRNAs

Introduction

Choriocarcinoma, which occurs secondary to hydatidiform moles, miscarriage, term pregnancy, preterm pregnancy, and ectopic pregnancy, is a highly aggressive and malignant trophoblastic neoplasm in the uterus.¹ If left untreated, choriocarcinoma tends to rapidly spread to other organs to establish distant metastases via blood vessels, causing distant organ failure.² Human chorionic gonadotropin (HCG) levels are the main diagnostic basis for gestational trophoblastic tumours. However, rising HCG concentrations can be observed in normal pregnancy, ectopic pregnancy, or gestational trophoblastic neoplasia, seriously affecting early screening and diagnosis of choriocarcinoma in the clinic.³ Therefore, it is necessary to elucidate the pathogenesis of choriocarcinoma and explore new markers for the early screening and diagnosis of choriocarcinoma.

MicroRNAs (miRNAs) are noncoding RNAs consisting of 20–25 nucleotides.⁴ miRNAs can regulate the posttranscriptional expression of target genes by binding to the 3ʹ-untranslated regions of messenger RNAs (mRNAs), leading to translational repression or mRNA cleavage.⁵ Increasing evidence has demonstrated that dysregulated miRNAs are involved in choriocarcinoma tumorigenesis, with altered expression patterns during tumour progression.^6,7 For example, miR-34a significantly inhibits the proliferation and invasiveness of choriocarcinoma cells via regulation of the Notch pathway,⁸ and miR-21 promotes the progression of choriocarcinoma.⁹ These findings showed that miRNAs regulate the occurrence and progression of choriocarcinoma. Unfortunately, the mechanism by which miRNAs regulate choriocarcinoma has not been systematically elucidated.

Based on high-throughput sequencing data from Gene Expression Omnibus (GEO) databases, comprehensive microarray analysis was utilized to explore the molecular mechanism to identify differentially expressed genes (DEGs) and miRNAs (DEMs) as well as transcription factors (TFs) involved in initiation and progression for further experimental validation.¹⁰ In our study, we identified DEMs and DEGs in choriocarcinoma cells by analysing one mRNA dataset (GSE20510) and two miRNA datasets (GSE32346 and GSE130489) downloaded from GEO databases. The top 400 upregulated and downregulated DEGs were selected, and protein-protein interaction (PPI) network construction and functional enrichment analysis were performed. Twenty-nine hub genes were obtained from 400 DEGs, and drug-hub gene networks were established. Moreover, we screened 70 differentially expressed target genes (DETGs) of DEMs derived from the intersection between predicted miRNA target genes and candidate DEGs. MiRNA-target gene regulatory networks were visualized by Cytoscape software. In addition, miRPath v.3.0 helped us perform functional enrichment analyses of selected hub DEMs (miR-29b-3p, miR-519c-3p and miR-520a-5p). TFs related to hub DEMs from the interaction network were identified by the iRegulon plugin. From the results of the database tools described above, we evaluated the effect of hub DEMs on the proliferation and migration abilities of choriocarcinoma cell lines through the transfection of miRNA mimics or inhibitors. We hope that our findings may provide potential therapeutic targets for inhibiting tumorigenesis and improving the prognosis of choriocarcinoma.

Materials and Methods

Data Acquisition

After careful review and comparison, mRNA expression datasets (GSE20510¹¹) along with miRNA expression from the GSE32346¹² and GSE130489¹³ datasets were downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/geo/) for subsequent analysis. The platform, variety of studies and proportions of choriocarcinoma cell samples and normal cell samples in each dataset are presented in Figure 1G. The above miRNA datasets were processed based on the GPL15497 and GPL26600 platforms. Hence, we normalized these data using R’s limma package.

Figure 1 Identification of differentially expressed miRNAs (DEMs) in choriocarcinoma. The before and after normalization and volcano plots of GSE20510 (A and B), GSE32346 (C and D) and GSE130489 (E and F). Characteristics of mRNA and miRNA expression profiles of choriocarcinoma (G). Venn diagram based on the miRNA expression profile of choriocarcinoma (H).

Screening of DEMs and DEGs

GEO2R (https://www.ncbi.nlm.nih.gov/geo/geo2r/) is an online analysis tool that allows for comparisons between choriocarcinoma cell samples and normal cell samples from normalized GEO datasets. The Benjamini and Hochberg false discovery rate (FDR) method was utilized to adjust the p-value to reduce false positivity. With the thresholds of adjusted p-value < 0.05 and |log2-fold change (FC)| > 1, the top 400 upregulated and downregulated given genes were selected as DEGs. Under the same conditions, DEMs were screened and analysed in overlapping sections of the GSE32346 and GSE130489 datasets through a Venn diagram to identify candidate DEMs.¹⁴ Heatmaps of DEMs and DEGs were visualized by TBtools, which is a toolkit that integrates various biological data handling tools.

Functional and Pathway Enrichment Analysis of DEGs

To further annotate the functions of DEGs, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were analysed using ClueGO V2.5.1,¹⁵ with thresholds of p-value <0.0 5 and gene count > 3.

Construction of PPI Network and Analysis of Modules

The STRING database (http://www.bork.emblheidelberg.de/STRING/) is an online tool that can provide predictive and visual analysis of interactive relationships among DEGs. We established PPI networks with an interaction score ≥0.4 used as the threshold. Cytoscape software was utilized to screen the top 400 upregulated and downregulated DEGs and construct the PPI network. Additionally, degree cut-off = 2, max. depth = 100, node score cut-off = 0.2 and k-Core = 2 were adopted as the cut-off criteria for identifying hub genes through Molecular Complex Detection (MCODE) v1.5.1.¹⁶

Drug-Hub Gene Networks Analysis

The drug-gene interaction database (DGIdb, www.dgidb.org) is a web resource that presents drug-gene interactions from databases, web resources, etc. The list of hub genes was uploaded to the DGIdb to predict targeted drugs. Finally, the drug-gene networks were constructed by Cytoscape.

Prediction of Target Genes for DEMs and miRNA–Target Gene Regulatory Network Construction

miRWalk2.0 (http://zmf.umm.uni-heidelberg.de/apps/zmf/mirwalk2/) was used to predict potential target genes of candidate DEMs based on six bioinformatic algorithms (miRWalk,¹⁷ Pictar2,¹⁸ PITA,¹⁹ RNA22,²⁰ RNAhybrid²¹ and TargetScan²²). Then, overlapping sections among these target genes and the top 400 upregulated and downregulated DEGs were obtained for further analysis. The miRNA-target gene regulatory networks between relative DEMs and the overlaps were visualized by Cytoscape. Among the upregulated miRNAs in this network, we selected the miRNA with the largest number of reverse regulation genes as the hub DEM. The same operation was performed for downregulated miRNAs in this network. Finally, the hub DEMs were selected from among the upregulated and downregulated DEMs based on the number of reverse-regulated genes.

Functional and Pathway Enrichment Analysis of Hub DEMs

GO and KEGG pathway analyses of hub DEMs were performed using miRPath v3.0 (http://snf-515788.vm.okeanos.grnet.gr/). GO functional analysis was mainly classified into biological process (BP), cellular component (CC), and molecular function (MF) categories. The heatmaps of GO and KEGG pathway analyses were also visualized with a p-value < 0.05 as the threshold.

Establishment of the miRNA-TFs-Target Genes Networks

Gene regulatory networks are controlled by miRNAs and TFs in a dependent or independent manner.²² Therefore, further research on the regulatory relationship between TFs and miRNAs can help elucidate the underlying mechanism of tumorigenesis. We selected hub DEMs and predicted potential TFs related to them via the iRegulon plug-in.²³ Relevant parameters were as follows: enrichment score threshold = 5.0, ROC threshold for AUC calculation = 0.05, rank threshold = 5000, minimum identity between orthologous genes = 0.05 and maximum FDR on motif similarity = 0.001. MiRNA-TF-target gene regulatory networks of hub DEMs were screened in Cytoscape.

Cell Culture and Transfection

The human trophoblast cell line HTR-8/SVneo and the choriocarcinoma cell lines JEG3 and JAR were purchased from American Type Culture Collection (ATCC) and cultured in RPMI 1640 medium (Gibco, USA) and MEM (Gibco, USA), respectively, containing 10% FEB, 100 U/mL penicillin, and 100 mg/mL streptomycin at 37°C in a humidified incubator with 5% CO₂. The cells were seeded in 6-well plates at the appropriate density and then transfected with negative control (NC) or corresponding hub DEM mimics or inhibitors (50 nM) using Lipofectamine 3000 (Invitrogen, USA) according to the manufacturer’s instructions. The sequences were synthesized by GenePharma (Shanghai, China).

Total RNA Extraction and RT-PCR Validation

Total cellular RNA was extracted by the TRIzol Total RNA method. Then, the RNA was reverse transcribed into cDNA using the RevertAid™ H Minus First Strand cDNA Synthesis Kit. RT-PCR was carried out by using SYBR^® Green MasterMix in an ABI PRISM^® 7500 Sequence Detection System (Applied Biosystems Inc: Foster City, CA). The reaction conditions were as follows: 45 seconds at 95°C, 40 cycles for 30 seconds at 95°C, 45 seconds at 60°C, 60 seconds at 72°C and a final cycle for 300 seconds at 72°C. The accumulation of PCR products was monitored in real time, and simultaneous determination of cycle threshold (Ct) values was also conducted. Gene expression values were determined using the ΔΔCT method with normalization to GAPDH. The primers were purchased from Sangon Biotech (Shanghai, China).

Cell Proliferation Assay

CCK-8 and 5-ethynyl-2́-deoxyuridine (EdU) assays were conducted to assess the proliferation of HTR-8/SVneo and JEG3/JAR cells following the manufacturer’s protocol. After 24 h of mimic or inhibitor transfection, the cells were seeded into 24-well (for EdU assays) and 96-well (for CCK-8 assays) plates at a suitable density. After 24 h of culture, CCK-8 solution (Beyotime, China) was added to each well (10 μL/well) and cultured for 1 h at 37 °C in an incubator. Finally, the absorbance at 450 nm, as the optical density (OD) of the samples, was detected by spectrophotometry at different time points. In addition, the cells in 24-well plates were incubated in EdU solution (100 μL/well) for 2 h, and then, after cell fixation with 4% paraformaldehyde for 30 min, permeabilization with 1% Triton X-100 for 10 min and staining with Hoechst 33,342 for 10 min, microphotographs were obtained to calculate the positive cell rate by fluorescence microscopy.

Cell Migration Assay

After 24 h of mimic or inhibitor transfection, the cells were harvested with 0.25% trypsin and resuspended to obtain a serum-free suspension. A 300 μL cell suspension at an appropriate concentration was placed into the upper layer of the Transwell chamber. Cell medium with 20% FEB was incubated in the lower layer of the Transwell chamber to promote cell migration. After 24 h of culture, the cells in the upper chamber were removed using cotton swabs, and those that migrated to the chamber membranes were fixed with 4% paraformaldehyde and stained with 0.1% crystal violet. Five randomly selected microscopic views were obtained through a microscope and analysed by ImageJ software to count the number of cell migration.

Statistical Analysis

All the experiments were repeated three times, and the data are presented as the mean ± SD. Analysis of variance was utilized to compare different sample means. SPSS for Windows, version 22.0 was used for analyses. Statistical analysis of variance (ANOVA) was used for statistical analysis between each group. P< 0.05 was considered statistically significant.

Results

Identification of DEGs and DEMs in Choriocarcinoma

As indicated in Figure 1, normalization of the expression profile datasets was carried out through the preprocess Core package in RStudio, and is the results are shown in boxplots. Two normalized miRNA expression profiles and one mRNA expression profile were downloaded from the GEO database. A total of 22 choriocarcinoma cell samples and 12 normal trophoblast samples were included in further analysis. A total of 2114 DEGs (1381 unregulated and 733 downregulated genes) were identified in GSE20510 with the mentioned threshold. Genes whose |log2 FC| ranked in the top 400 upregulated and downregulated genes were chosen as candidate DEGs (Supplementary Table 1). With the same threshold set, 149 and 251 DEMs were identified in GSE32346 and GSE130489, respectively. Furthermore, 140 candidate DEMs were identified by Venn diagram analysis overlapping these two datasets (Figure 1H). LogFC of 140 candidate DEMs from GSE32346 and GSE130489 were displayed in Supplementary Table 2.

Functional and Pathway Analysis of DEGs

GO and KEGG pathway analyses were performed on the 400 candidate DEGs. As shown in Figure 2, BP-associated results indicated that candidate DEGs were particularly enriched in cell migration, extracellular matrix organization and regulation of cell population proliferation; CC-associated results showed that candidate DEGs were mainly enriched in basement membrane, extracellular matrix and tight junction; MF-associated results indicated that candidate DEGs were significantly enriched in regulation of peptidase activity, regulation of transferase activity and positive regulation of catalytic activity. Furthermore, KEGG pathway results demonstrated that candidate DEGs were primarily enriched in the focal adhesion, Wnt and PI3K-Akt signalling pathways. Significant results of the GO enrichment analysis and KEGG pathway analysis of candidate DEGs in choriocarcinoma are shown in Supplementary Tables 3 and 4.

Figure 2 GO terms and KEGG pathway analysis of candidate DEGs in choriocarcinoma. (A) BP-associated category; (B) CC-associated category; (C) MF-associated category; (D) KEGG pathways. Abbreviations: BP, biological processes; CC, cellular component; MF, molecular function; KEGG, Kyoto Encyclopedia of Genes and Genomes.

PPI Network Construction and Module Analysis

To identify potential interactions among 400 candidate DEGs, a PPI network was constructed through the STRING database and visualized by Cytoscape (Figure 3A). Furthermore, hub genes were validated by MCODE, and we obtained 14 cluster modules. The highest-rated module included 29 nodes and 243 edges and was subjected to GO and KEGG pathway analysis (Figure 3B). The hub genes were mainly abundant in 167 GO terms, including positive regulation of developmental process, cell migration, basement membrane, extracellular matrix, tight junction and regulation of transferase activity. In addition, KEGG pathway analysis indicated that these genes were abundant in focal adhesion, proteoglycans in cancer, the PI3K-Akt signalling pathway and pathways in cancer (Supplementary Table 5).

Figure 3 PPI network construction and module analysis. A total of 400 DEGs were visualized in the DEG PPI network (A). Module analysis using MCODE (B). The nodes represent proteins, and the edges represent the interaction of the proteins. Red nodes represent upregulated genes, and blue nodes represent downregulated genes.

Drug-Gene Networks

A total of 99 potential drugs for treating choriocarcinoma patients were identified when the drug-gene interactions were explored using DGIdb (Supplementary Table 6). Additionally, drug-gene networks were constructed by Cytoscape (Figure 4). Most potential drugs might interact with hub genes either in unknown manners or as inhibitors.

Figure 4 Regulatory networks of drug-hub gene interactions. Red nodes and blue nodes represent upregulated DEGs and downregulated DEGs, respectively. Pink nodes represent target drugs.

Prediction of Target Genes for DEMs and miRNA–Target Gene Regulatory Network Establishment

Increasing evidence has demonstrated that miRNAs regulate tumorigenesis and progression mainly by binding to the 3ʹUTRs of mRNAs and thereby triggering degradation. MiRWalk2.0, a miRNA-target online tool, was implemented to predict target genes of candidate DEMs. After calculation by six bioinformatics algorithms, 2002 potential target genes were obtained and intersected with 400 candidate DEGs. Seventy overlapping genes, including 34 upregulated genes and 36 downregulated genes, were identified as DETGs for further analysis (Figure 5 and Supplementary Table 7). These DETGs were subjected to GO and KEGG pathway analyses and significantly enriched in blood vessel development, regulation of cell adhesion, cell migration, negative regulation of cell proliferation, and chemotaxis (Supplementary Table 8). The miRNA-target gene regulatory network between DEMs and DETGs was visualized using Cytoscape. As indicated in Figure 6, diamonds represent DEMs, red ellipses represent upregulated DETGs, and blue ellipses represent downregulated DETGs. The regulatory network included 34 upregulated DETGs, 36 downregulated DETGs, 79 DEMs and 253 interactions. According to the mechanism of miRNA negative regulation, DEMs were divided into two groups including upregulated and downregulated DEMs and ranked with the number of screened DETGs with an opposite trend to DEMs. Among the 79 DEMs, there were 64 upregulated DEMs and 15 downregulated DEMs. Among the upregulated DEMs, miR-29b-3p regulated 7 downregulated target genes, and it was the upregulated DEM with the largest number of opposite target genes (Supplementary Table 9). miR-519c-3p and miR-520a-5p were downregulated DEMs, both of which regulate three downregulated genes, and they are the downregulated DEMs with the largest number of downregulated genes (Supplementary Table 9). Herein, we identified 3 hub DEMs (upregulated miR-29b-3p and downregulated miR-519c-3p and miR-520a-5p). miR-29b-3p, miR-519c-3p and miR-520a-5p regulated 7, 3 and 3 screened DETGs, respectively.

Figure 5 Identification of DETGs between DEGs and DEMs predicted target genes.

Figure 6 Construction of miRNA-target gene regulatory network. Red nodes and blue nodes indicate upregulated and downregulated DETGs, respectively. Diamonds represent DEMs.

Functional and Pathway Analysis of Hub DEMs

The hub DEMs were subjected to GO and KEGG pathway analyses through miRPath v3.0. As indicated in Figure 7, these hub DEMs were enriched in cellular nitrogen compound metabolic process, protein modification and biosynthetic process in the BP-associated category, organelle and ion binding in the CC-associated category and nucleic acid binding transcription factor activity in the MF-associated category. The KEGG pathway analysis results demonstrated that hub DEMs were enriched in proteoglycans in cancer and lysine degradation.

Figure 7 GO and KEGG analysis of hsa-miR-29b-3p, hsa-miR-519c-3p and hsa-miR-520a-5p. BP (A), CC (B) and MF (C) terms and enriched (D) KEGG pathways.

Construction of the miRNA-TF-Target Gene Regulatory Network

To investigate the regulatory mechanism between hub DEMs and DETGs, miRNA-TF-target gene networks for these hub DEMs were established. The results showed that the regulatory networks for miR-29b-3p, miR-519c-3p and miR-520a-5p consisted of four TFs (STAT3, TFAP2C, TFAP2A and CHD1), three TFs (GATA2, HOXA11 and SMC3) and three TFs (FLI1, HNF4A and SOX11) (Figure 8).

Figure 8 Regulatory networks of the key miRNAs, target genes and transcription factors. Network of hsa-miR-29b-3p and four transcription factors (A); network of hsa-miR-519c-3p and three transcription factors (B); network of hsa-miR-520a-5p and three transcription factors (C). Green octagons indicate transcription factors. Purple circles indicate target genes regulated by hub DEMs and transcription factors, and blue circles indicate target genes regulated by hub DEMs.

Validation of Hub DEMs

As a rare disease, choriocarcinoma samples are difficult to obtain, and most tumour-associated databases do not contain their clinical data and tissue expression. Therefore, three hub DEMs were validated by RT-PCR and compared with normal trophoblast cells and choriocarcinoma cells. Consistent with the analysis of miRNA profile datasets, we found that miR-29b-3p was markedly upregulated while miR-519c-3p and miR-520a-5p were downregulated in JEG3 and JAR cells (Figure 9). The findings indicated that the analysis results were reproducible and reliable.

Figure 9 The expression of hsa-miR-29b-3, hsa-miR-519c-3p and hsa-miR-520a-5p in the comparison between the choriocarcinoma cell group and HRT-8 group. *p < 0.05 compared with HRT-8 group.

The Effect of Hub DEMs on Choriocarcinoma Cell Proliferation

MiR-29b-3p, miR-519c-3p and miR-520a-5p were overexpressed or silenced in JEG3 and JAR cells to explore the effect of hub DEMs on choriocarcinoma cell proliferation. A CCK-8 assay was used to examine proliferation. As demonstrated in Figure 10A–D, overexpression of miR-29b-3p enhanced cell proliferation in JEG3 and JAR cells compared with NC mimics, while silencing of miR-29b-3p reversed the above effects. Moreover, miR-519c-3p and miR-520a-5p mimics inhibited cell proliferation, whereas downregulation of miR-519c-3p and miR-520a-5p exerted opposite effects. Consistent with the results of the CCK-8 assay, the EdU assay results also verified our hypothesis through quantification of percentage of EdU positive cells (Figure 10E and F).

Figure 10 The effect of overexpression or silencing of hsa-miR-29b-3p, hsa-miR-519c-3p and hsa-miR-520a-5p on cell proliferation in JEG-3 and JAR cells. Viability of JEG-3 and JAR cells transfected with hub DEM mimics or inhibitors measured by CCK-8 (A–D) and EdU assays (E and F). The percentage of Edu positive JEG-3 and JAR. Proliferation ratio is assessed using EdU staining/Hoechst 33,342 (G and H). *p < 0.05 compared with NC-mimics in JEG-3 and JAR cells. **p < 0.05 compared with NC-inhibitor in JEG-3 and JAR cells.

The Effect of Hub DEMs on Choriocarcinoma Cell Migration

We next assessed the functional effects of miR-29b-3p, miR-519c-3p and miR-520a-5p on cell migration. After transfection with miR-29b-3p mimics, the migration abilities of JEG3 and JAR cells were greatly improved, with many cells crossing through the membranes compared with the NC mimic group. In contrast, the miR-29b-3p inhibitor notably decreased the migration abilities of the cells. In addition, cells transfected with miR-519c-3p and miR-520a-5p mimics migrated quickly through the membranes, and an increased number of cells was observed on the upper layer of the Transwell chamber. Simultaneously, the migration abilities of cells transfected with miR-519c-3p and miR-520a-5p mimic inhibitors dramatically decreased through calculating numbers of cell migration, as shown in Figure 11. In summary, upregulation of miR-29b-3p and downregulation of miR-519c-3p and miR-520a-5p contributed to the migration ability of choriocarcinoma cells.

Figure 11 The effect of overexpression or silencing of hsa-miR-29b-3p, hsa-miR-519c-3p and hsa-miR-520a-5p on cell migration in JEG-3 (A) and JAR (B) cells. Cell migration was measured by Transwell assays. *p < 0.05 compared with NC-mimics in JEG-3 and JAR cells. **p < 0.05 compared with NC-inhibitor in JEG-3 and JAR cells.

Discussion

Choriocarcinoma is a rare malignant trophoblastic neoplasm characterized by rapid proliferation and strong migration ability.² Unlike other types of cancer originating from somatic cells, choriocarcinoma is derived from chorionic villi and extravillous trophoblasts, indicating that other tumour-related targets may not be suitable to prevent and treat choriocarcinoma.¹ Studies have shown that the incidence of choriocarcinoma is 1:2416–1:2888 pregnancies in China,²⁴ 1:40,000 pregnancies in the United States and Europe,²⁴ and 1:193 pregnancies in the Philippines.²⁵ The mortality of low-risk and high-risk Chinese women is 0.3% and 15.1%, respectively.^3,26 Unfortunately, the mechanism underlying choriocarcinoma is still not fully elaborated, and early interventions are lacking. Therefore, identification of hub miRNAs and mRNAs involved in the tumorigenesis and progression of choriocarcinoma is urgently needed. Moreover, the synergistic effects of TFs and miRNAs in genetic regulatory networks are new hot spots for deepening the understanding of molecular mechanisms.²⁷ Unfortunately, to date, no associated research has been reported to construct miRNA-TF-target gene regulatory networks in choriocarcinoma. To the best of our knowledge, we are the first to attempt to elucidate interactions between DEGs and DEMs and construct miRNA-TF-target gene regulatory networks, providing novel therapeutic targets for early diagnosis, treatment and prevention.

Based on the gene expression profile, we identified the top 400 upregulated and downregulated DEGs from 2115 DEGs as candidate DEGs for further analysis. To explore DEG-associated functions and signalling pathways, we conducted GO category and pathway term enrichment analyses utilizing the ClueGO plugin. Functional enrichment analysis demonstrated that the DEGs were correlated with cell activity, including regulation of cell migration and cellular component movement, which are associated with tumour progression. Similarly, pathway enrichment analysis revealed that the DEGs were enriched in proteoglycans in cancer and were involved in cell proliferation, migration and neovascularization. SPARC is a matricellular glycoprotein that modulates the pathologies associated with tumorigenesis.²⁸ Jiang et al found that the mRNA and protein expression of SPARC was downregulated in JEG-3 cells compared with HTR-8/SVneo cells.²⁹ Furthermore, miRNAs in cancer were the primary pathway for DEG enrichment. This result suggested that miRNAs may play a vital role in mediating choriocarcinoma progression.

To explore the linkage between DEGs, we constructed a PPI network of candidate DEGs and obtained 29 hub genes, including upregulated DEGs and downregulated DEGs. Some of these hub genes were demonstrated to be involved in choriocarcinogenesis and progression in a previous study. Chen et al reported that TIMP3 was expressed less in JEG-3 and Bewo cells, consistent with the findings of the current study. Transferrin with adenovirus overexpressing TIMP3 decreased the possibility of tumour growth, adhesion and invasion.³⁰ Similarly, TIMP1 was markedly upregulated in JEG-3 cells compared with HTR-8/SVneo cells in a previous study.³¹ CD44, the receptor of hyaluronic acid, was demonstrated to be highly expressed in choriocarcinoma cells, facilitate tumour cell adhesion to the extracellular matrix and basement membrane and to provide sites for tumour invasion.³² Lan et al revealed that CDH1 expression was significantly increased in JEG-3 and JAR cells compared with HTR-8/SVneo cells, and DNA methylation of the CDH1 promoter contributed to choriocarcinoma cell motility and invasion.³³ TGF-β2 was dramatically downregulated in choriocarcinoma cell lines (JEG-3 and JAR) compared with a normal trophoblast cell line (HTR-8/SVneo).³⁴ Insulin-like growth factor binding protein (IGFBP) has become a popular gene family for oncology research and is involved in the proliferation, metastasis and invasion of tumour cells, thus mediating the development of tumours, including choriocarcinoma. Liu et al found that IGFBP7 was abundantly expressed in normal trophoblast tissue samples rather than hydatidiform moles and inhibited the proliferation and invasion capacities of JEG-3 cells.³⁵ IGFBP3 was dramatically downregulated in JEG-3 and JAR cells compared with normal trophoblastic cells.³⁶ Jiang et al showed that SPARC depletion upregulated the expression of IGFBP4, which promoted the invasion of JEG-3.²⁹ The above results were consistent with our study, strongly indicating that hub genes mostly participated in tumour proliferation and invasion.

Chemotherapy drug resistance remains one of the most important issues in the outcome of cancer treatment and affects the survival prognosis of patients. Thus, predicting and developing new drugs and designing a better treatment regimen to reduce drug resistance is crucial.^37,38 We used DGIdb to predict drug-gene interactions and identified 99 drugs targeting 29 hub genes, 47 of which were FDA-approved drugs. Many kinds of drugs are used to treat tumours. For example, bevacizumab is used to treat advanced ovarian cancer,³⁹ metastatic colorectal cancer⁴⁰ and non-small cell lung cancer.⁴¹ Bicalutamide is used to treat metastatic hormone-sensitive prostate cancer.⁴² Lapatinib is widely used in patients with non-small cell lung cancer,⁴³ breast cancer,⁴⁴ aggressive prolactinomas⁴⁵ and oesophageal cancers.⁴⁶ The nonconventional drugs predicted by the study may have a new use in choriocarcinoma.

Increasing evidence has demonstrated that miRNAs, as pivotal functional molecules, are regarded as biomarkers for early screening and potential targets to block tumorigenesis and tumour progression.¹² Various miRNAs are reportedly involved in the proliferation and migration of choriocarcinoma, suggesting that altered miRNA expression and miRNA-mediated downstream genes may be involved in the pathogenesis of choriocarcinoma.^12,47 In our present study, we identified 140 consistently altered DEMs in both the GSE32346 and GSE130489 datasets. Moreover, we established a miRNA-target gene regulatory network according to the DEMs and the overlaps among the DEGs and potential candidate target genes. DEMs were divided into two groups, upregulated and downregulated DEMs, and ranked with the number of screened DETGs with an opposite trend to DEMs. MiR-29b-3p, miR-519c-3p and miR-520a-5p were identified as hub DEMs.

As upregulated DEMs, miR-29b-3p regulated 7 downregulated DETGs (PXDN, SPARC, OXTR, LOX, EMP1 COL3A1 and FBN1) and was identified as a tumour promoter by increasing cell proliferation, migration and chemotherapy tolerance through various pathways in different cancers. A previous study demonstrated the upregulation of miR-29b-3p in breast cancer tissues.⁴⁸ Some studies found that miR-29b-3p improved the invasiveness of the breast cancer cell line MCF-7 by inhibiting the expression of SPARC and COL4A2.⁴⁹ Notably, Sugio A showed that in ovarian cancer, the expression of miR-29b-3p was positively correlated with the level of chemotherapy tolerance.⁵⁰ Okamoto K found that miR-29b-3p could improve the chemotherapy tolerance of cholangiocarcinoma to gemcitabine.⁵¹ Similarly, miR-29b-3p was highly expressed in a bladder cancer chemotherapy-sensitive cell line.⁵² The main treatment modality for choriocarcinoma is chemotherapy, but 20% of patients still fail to achieve a cure due to resistance after chemotherapy, so chemotherapy resistance has become a more difficult problem in clinical practice.² Based on the data described above, we speculate that miR-29b may be related to chemoresistance in choriocarcinoma and may become a clinical diagnostic and therapeutic target. Some of these miR-29b-3p-mediated DETGs should be considered. Overexpressed SPARC may inhibit gastric cancer cell invasion and migration by downregulating the expression of VEGF, MMP-9, and MMP-7.⁵³ In addition, methylation and low expression of FBN1 contributed to the progression of tumour cells in papillary thyroid carcinoma,⁵⁴ ovarian cancer⁵⁵ and endometrial cancer.⁵⁶

Numerous studies have shown that miR-519c-3p, a tumour suppressor gene, is downregulated in most cancers,⁵⁷ including cervical cancer, ovarian cancer, hepatocellular carcinoma and breast cancer. MiR-519c-3p overexpression suppresses tumour angiogenesis by targeting HIF-1α in lung adenocarcinoma cells.⁵⁸ In this study, miR-29b-3p regulated 3 upregulated DETGs (MAP7, RASL11B and SOX4). SOX4 overexpression activates Wnt/β-catenin and promotes the proliferation and migration of osteosarcoma cells.⁵⁹ The miR-519c-3p-regulated target gene SOX4 is regarded as a gene associated with embryonic development and an important transcription factor for T and B cells.⁶⁰ In addition, MAP7 interacts with RC3H1 to upregulate the cell cycle progression of cervical cancer cells by activating NF-κB signalling.⁶¹ Previous studies showed that low expression of miR-520a-5p promotes cancer cell migration and invasion in tissues, such as lung cancer⁶² and breast cancer.⁶³ A prospective pilot study showed that miR-520a-5p could be used as an effective foetal/placenta-specific marker for the diagnosis and continuous follow-up of GTD/GTN, especially for patients with HCG false-positive and HCG-negative results.⁶⁴ PGF, FZD5 and DTX4 are the gene targets of hsa-miR-520a-5p. PGF is a potent stimulant of cancer invasion, mediating the growth of vascular endothelial cells, activating angiogenesis, and promoting tumour cell invasion and metastasis.⁶⁵ Other studies have revealed that PGF regulates the expression of epithelial-mesenchymal transition-related proteins in cervical cancer.⁶⁶ FZD5 induced the proliferation of cancer cells by activating the Wnt signalling pathway in triple-negative breast cancer.⁶⁷

To further explore the function of hub DEMs, miRNA pathway enrichment analysis indicated that hub DEMs were enriched in the cellular nitrogen compound metabolic process, protein modification and biosynthetic process in the BP-associated category, organelle and ion binding in the CC-associated category and nucleic acid binding transcription factor activity and ion binding in the MF-associated category. The KEGG pathway analysis results demonstrated that hub DEMs were enriched in proteoglycans in cancer and lysine degradation. Moreover, we constructed a miRNA-TF-target gene regulatory network. MiR-29b-3p and four TFs, STAT3, TFAP2C, TFAP2A and CHD1, were verified to be coregulators of PXDN, ATP2B4, EMP1, FBN1, SPARC, BTG2, OXTR, LOX, and KDM5B. MiR-519c-3p and three TFs, GATA2, HOXA11 and SMC3, were verified to be coregulators of SOX4, DAZAP2, RASL11B, ZFPM2, NETO2 and MAP7. Additionally, miR-520a-5p and three TFs, FLI1, HNF4A and SOX11, were verified to be coregulators of PGF, DTX4 and FZD5. The networks provided us with a new way to elucidate the underlying mechanism and molecular targets of choriocarcinoma.

Since choriocarcinoma is a relatively rare disease worldwide, current data platforms and databases lack a summary and statistical analysis of clinicopathological characteristics and data related to the prognosis of choriocarcinoma. Moreover, the main treatment modality of choriocarcinoma is chemotherapy, so clinical samples and tissues are difficult to obtain.⁶⁸ Therefore, we selected two choriocarcinoma cell lines, JEG-3 and JAR cells, for further study and validation. Compared with HTR-8/SVneo,²⁹ upregulation of miR-29b-3p and downregulation of miR-519c-3p and miR-520a-5p were verified by RT-PCR, which was consistent with microarray analysis results. The results of overexpression and silencing these hub DEMs on cell proliferation and migration indicated that miR-29b-3p facilitated proliferation while miR-519c-3p and miR-520a-5p inhibited proliferation and migration in JEG-3 and JAR cells. Therefore, it is worth investigating the functions and underlying mechanisms of miR-29b-3p, miR-519c-3p and miR-520a-5p in choriocarcinoma.

Some limitations existed in the present study. Due to the lack of adequate clinical samples, the tissue expression and prognostic value of hub DEMs could not be evaluated. In addition, the expression and regulatory mechanism of TFs were not determined through experimental tests. Whether other noncoding RNAs, such as lncRNAs and circRNAs, play an essential role in the development of choriocarcinoma also deserves to be explored. Due to the lack of functional verification of the three key miRNAs in HTR-8 cells, there may be some deficiencies in tumour-specific verification of miRNAs. However, our experiment confirmed that these three key miRNAs can influence the regulation of choriocarcinoma cell proliferation and migration. In addition, we only performed in vitro experiments without an internal environment in vivo. We will perform animal experiments to strengthen the persuasiveness of the underlying mechanisms of these hub DEMs in the future.

Conclusion

In conclusion, we identified 400 DEGs and 140 DEMs based on bioinformatics analysis. Among these DEMs, three hub DEMs, miR-29b-3p, miR-519c-3p and miR-520a-5p, were confirmed in vitro as potential regulators of the proliferation and migration of choriocarcinoma cells. In addition, regulatory network analysis indicated the comodulating relationships among hub DEMs, TFs and corresponding DETGs in choriocarcinoma, which would provide novel diagnostic and therapeutic targets for choriocarcinoma.