The SiHa (HTB-35), HeLa (CCL-2) and C-33A (HTB-31^TM) cell lines were purchased from American Type Culture Collection (ATCC, Manassas, VA, USA). SiHa and C-33A cells were cultured in high-glucose Dulbecco’s Modified Eagle’s Medium (DMEM) (11-965-092, Gibco^TM; Carlsbad, CA, USA). HeLa cells were cultured in a 1:1 mixture of DMEM/F12 medium (DMEM/F12) (Cat: DFP44, Caisson Labs, Inc, Smithfield, UT, USA). All culture media were supplemented with 10% fetal bovine serum (FBS) (Cat: A5256701, Gibco^TM, Grand Island, NE, USA), 100 U/mL penicillin, 100 µg/mL streptomycin and 0.25 µg/mL Amphotericin B. Cells were grown at 37 °C in a 5% CO₂ incubator. All cell lines were authenticated by STR profiling and tested negative for mycoplasma.

Total RNA was isolated using TRIzol (Cat: 15596026, Invitrogen, Waltham, MA, USA) according to the manufacturer’s instructions. The concentration and purity of total RNA were determined with a NanoDrop^TM 2000/2000c UV-Vis spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) based on the A260/A280 ratio. The integrity of total RNA was evaluated through 1% agarose gel electrophoresis. Total RNA was treated with RNase-free DNase I (Cat: EN0521, Thermo Scientific, Waltham, MA, USA). Finally, DNase-treated RNA was re-quantified using the NanoDrop^TM 2000/2000c UV-Vis spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) based on the A260/A280 ratio, followed by storage at –80 °C prior to cDNA synthesis.

cDNA synthesis for miR-19a-3p was performed with the TaqMan® MicroRNA Reverse Transcription kit (Cat: 4366596, Applied Biosystems, Vilnius, Lithuania) according to the manufacturer’s instructions. The reaction mixture contained 5 ng of total RNA, RT Primer Pool (10 µM), dNTPs with dTTP (10 µM), MultiScribe Reverse Transcriptase (50 U/µL), RT Buffer (10X), RNase Inhibitor (20 U/µL) and Nuclease-free water. Thermocycler settings were as follows: 30 min at 16 °C, 30 min at 42 °C, and 5 min at 85 °C.

cDNA synthesis for BDNF-AS was conducted using SuperScript III reverse transcriptase (Cat: 18080-044, Thermo Fisher Scientific, Waltham, MA, USA) following the manufacturer’s instructions.

cDNAs concentrations were quantified using the NanoDrop^TM 2000/2000c UV-Vis spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA), followed by storage at –20 °C prior to downstream use.

miR-19a-3p levels were quantified with the hsa-miR-19a-3p assay ID: 000395 (Cat: Assay ID: 000395, Applied Biosystems, Waltham, MA, USA) and the TaqMan microRNA assay kit (Cat: 4427975, Applied Biosystems, Waltham, MA, USA) using the ABI 7500 Real-Time PCR System (Applied Biosystems, Waltham, MA, USA). The reaction mixture contained 1 µL of cDNA, TaqMan™ Small RNA Assay (20$\times$), PCR Master Mix (2$\times$), and Nuclease-free water. The thermocycler settings were: 20 s at 95 °C, followed by 40 cycles of 95 °C for 3 s and 60 °C for 30 s.

BDNF-AS expression levels were analyzed by RT-qPCR in reactions containing 200 ng of cDNA, primers at a final concentration of 1 µM, Power SYBR® Green PCR Master Mix 2$\times$ (Cat: 4368577, Applied Biosystems, Waltham, MA, USA), and Nuclease-free water. Analyses were performed using a PCR CFX96TM Real-Time PCR Detection System (Bio-Rad, Foster City, CA, USA) with the following thermocycler settings: 10 min at 95 °C, followed by 40 cycles of 95 °C for 20 s, 60 °C for 20 s, and 72 °C for 50 s. Primers used in this study were BDNF-AS sense: 5^′-CCGTGAGAAGATCTCATTGGG-3^′ and BDNF-AS antisense: 5^′-GGGTCACAAGTCACGTAGCA-3^′; Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH) sense: 5^′-CCAGGTGGTCTCCTCTGACTT-3^′ and GAPDH antisense: 5^′-GTTGCTGTAGCCAAATTCGTTGT-3^′.

The relative expression of miR-19a-3p and BDNF-AS was determined using the 2^{-Δ⁢Δ⁢Ct} method [18] with the respective normalization to U6 and GAPDH as housekeeping controls. Two biological replicates were included in these analyses, each with two technical replicates. miR-19a-3p and BDNF-AS expression data were Log₂-transformed. Significance was assessed using Student’s t-tests in GraphPad Prism v8.0.2 (GraphPad Software, La Jolla, CA, USA). A p-value 0.05 was considered statistically significant.

DEmiRNAs were identified when comparing CC and normal tissue samples from the GSE86100 dataset [19] using the GEO2R software (https://www.ncbi.nlm.nih.gov/geo/geo2r/) [20]. This dataset was selected based on the following criteria: (1) all cervical specimens were confirmed as squamous cell carcinoma through pathological examination, (2) patients had not received any chemotherapy treatment, and (3) control samples were confirmed to be normal using cytologic tests. The GSE86100 dataset (Platform: GPL19730, Agilent - 046064 Unrestricted_Human_miRNA_V19.0_Microarray, Probe Name version) consists of 6 HR-HPV-positive CC tissue samples (named as CC samples) and 6 uninfected normal mucosa tissue samples (named as controls samples). The microarray platform version used for this dataset was in format 8 $\times$ 60, including probes corresponding to 2006 human miRNAs. DEmiRNAs were identified based on the following criteria when comparing 6 CC and 6 controls samples: Fold Change (FC) $>$2 and p-value 0.001.

miR-19a-3p and BDNF-AS expression analyses were performed in CC and normal tissue samples from the Cancer Genome Atlas (TCGA) dataset using the miR-TV [21] and OncoDB [22] databases. Expression data were downloaded and transformed into the RPM Log₂ format. Significant differences were determined using Mann-Whitney U-tests, with p-value 0.05 as the threshold for significance.

miR-19a-3p expression was additionally analyzed in the GSE30656 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE30656) [23] dataset, while BDNF-AS expression was analyzed in the GSE67522 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE67522) [24, 25] and GSE63514 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE63514) [26] datasets from the GEO [20, 27] database. The GSE30656 dataset (Platform: GPL6955, Agilent-016436 Human miRNA Microarray 1.0, Feature Number version) contains 19 cervical adenocarcinoma/squamous cell carcinoma tissue samples and 10 cervical squamous epithelial samples with normal histology. The GSE67522 dataset (Platform: GPL10558, Illumina HumanHT-12 V4.0 expression beadchip) includes 20 HPV16 positive CC tissue samples and 22 HPV negative histologically normal controls/HPV16-positive non-malignant tissue samples. The GSE63514 dataset (Platform: GPL570, [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array) contains 28 CC tissue samples and 24 normal tissue samples. Expression analysis was done with the GEO2R [20] software. Expression data were Log₂-transformed, and differences between groups were analyzed with Student’s t-tests. A p-value 0.05 was considered significant.

Differentially expressed genes (DEGs) between CC patient samples with high and low levels of BDNF-AS expression from the TCGA-CESC dataset [28] were identified using the TCGA biolinks and DESeq2 packages [29, 30] using the R program, selecting those. mRNAs with an FC $>$1.5 and an adjusted p-value 0.05. These DEGs were visualized using the ggplot package of R program [31].

ROC curves analyses were performed based on CC and normal tissue samples from the GSE30656 dataset for miR-19a-3p, and the GSE67522 and GSE63514 datasets for BDNF-AS. These analyses were performed in GraphPad Prism v8.0 software (GraphPad Software, La Jolla, CA, USA). Area Under the Curve (AUC) values and 95% confidence intervals (CI) were calculated, and the Youden index method was used to calculate cut-off, sensitivity, specificity, and p-values [32].

Overall survival (OS) was evaluated using Kaplan-Meier curves based on CC tissue samples from the TCGA dataset included in the Kaplan-Meier Plotter database [33], comparing patient outcomes when stratified based on whether miR-19a-3p or BDNF-AS expression levels fell above or below the median value. Hazard ratios (HRs) and log-rank p-values were calculated, and p-value ˂ 0.05 was considered significant.

The target mRNAs of miR-19a-3p were identified with the miRDB [34] and TargetScan [35] databases using the default parameters. Common target mRNAs shared between these databases were identified with a Venn diagram. miR-19a-3p binding sites in BDNF-AS were identified using the Ensembl database [36].

Pathways and biological processes were analyzed using the Kyoto Encyclopedia of Genes and Genomes (KEGG) 2021 Human, Molecular Signatures Database (MSigDB) Hallmark 2020 and GO Biological Process 2021 libraries in the Enrichr database [37, 38, 39, 40]. Log₂-transformed p-values were calculated for each term, with p-value 0.05 being considered significant. The Gene Set Analysis (GSEA) was done with R software. The Normalized Enrichment Score (NES) was calculated, and a p-value 0.05 was considered as statistically significant.

The correlation between miR-19a-3p and BDNF-AS expression in CC tissue samples in the ENCORI/StarBase (Encyclopedia of RNA interactomes) v2.0 database [41] was analyzed. The expression of miR-19a-3p was transformed into the Log₂ (RPM+0.01) format, while BDNF-AS expression was transformed into the Log₂ (FPKM+0.01) format. Correlation coefficient (R) values were calculated, and p-value 0.05 was considered significant.

Genes correlated with BDNF-AS expression were identified in CC tissue samples from the TCGA dataset using the cBioportal database [42]. Genes with a Spearman’s Correlation $>$0.30 or –0.30, and q-value 0.001 were selected.

A protein-protein interaction (PPI) network was constructed using the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) v11.5 database [43], with the following parameters: (1) minimum required interaction score set to the highest confidence (0.900); (2) Active interaction sources: Experiments, Gene Fusion, and Co-expression; and (3) Hide disconnected nodes in the network: ON. A PPI enrichment p-value 0.05 was considered statistically significant.

To identify DEmiRNAs in CC, we analyzed the expression of miRNAs in CC tumor samples and normal tissue samples (n = 6 each) using the GSE86100 dataset with the GEO2R tool. Uniform manifold approximation and projection (UMAP) dimensionality revealed a clear separation between these tissue types, confirming that both groups were clearly distinguished (Fig. 1A). A total of 24 DEmiRNAs were identified between CC and control samples, of which 9 were downregulated, and 15 were upregulated in CC tissue samples (Table 1 and Fig. 1B).

Fig. 1.

Expression, diagnostic and prognostic value of miR-19a-3p in CC. (A) UMAP analysis of 6 CC and 6 normal tissue samples from the GSE86100 dataset. (B) Volcano plot of DEmiRNAs from the GSE86100 dataset. (C,D) Validation of miR-19a-3p expression in CC and normal tissue samples from the TCGA and GSE30656 datasets using the miR-TV and GEO databases. (E) Prognostic value of miR-19a-3p expression in patients with CC from the TCGA dataset through a Kaplan-Meier curve using the Kaplan-Meier Plotter database. Red line: high miR-19a-3p expression; Black line: low miR-19a-3p expression. (F) Diagnostic value of miR-19a-3p expression in CC and normal tissue samples from the GSE30656 dataset using a ROC curve in GraphPad Prism. *p-value 0.05; **p-value 0.01. UMAP, uniform manifold approximation and projection; TCGA, The Cancer Genome Atlas; RPM, Reads Per Million; GEO, Gene Expression Omnibus; ROC, receiver operating characteristic; AUC, Area Under the Curve.

While certain miRNAs identified among these DEmiRNAs (Table 1) have previously been characterized in CC, including miR-21-5p [44, 45] and miR-154-5p [46, 47], others remain poorly characterized in this cancer type. Given the relatively limited research on miR-19a-3p in CC to date, it was chosen as a target for further analysis.

To validate changes in miR-19a-3p expression in CC, we analyzed the levels of this miRNA in CC and normal tissue samples from the TCGA and GSE30656 datasets using the miR-TV and GEO databases. We observed elevated miR-19a-3p expression in CC in these two independent datasets (Fig. 1C,D), supporting our initial findings.

To determine the potential prognostic and diagnostic relevance of miR-19a-3p in CC, we evaluated patient OS based on Kaplan-Meier analysis and AUC values for ROC curves when stratifying patients from the TCGA and GSE30656 datasets according to miR-19a-3p expression levels using the KM-Plotter database and GraphPad Prism, respectively. This analysis revealed that miR-19a-3p expression is not associated with OS in patients with CC (Fig. 1E). However, we found a good AUC (0.8526, p-value = 0.0021) was observed, with excellent sensitivity and specificity of 94.74% and 70.00%, respectively (Fig. 1F), for a cut-off value of 6.602.

These results suggest that miR-19a-3p expression is elevated in CC and that this miRNA may serve as a potential diagnostic biomarker in this type of cancer, supporting its biological and clinical relevance in this malignancy.

To explore the potential role of miR-19a-3p in CC, we next identified potential targets of this miRNA using miRDB and TargetScan databases, identifying 829 common targets in both databases (Fig. 2A). Next, pathway and biological process analyses were performed for these 829 potential targets using the KEGG 2021 Human and GO Biological processes libraries from the Enrichr database. As shown in Fig. 2B, KEGG analysis revealed that the potential targets of miR-19a-3p participate in key pathways to CC, such as the mTOR signaling pathway, Forkhead Box O (FoxO) signaling pathway, and signaling pathways regulating pluripotency of stem cells. Interestingly, the identified biological processes were found to be involved in transcriptional regulation and protein phosphorylation, including the Regulation of transcription, DNA-template, Regulation of transcription by RNA polymerase II, and Negative regulation of gene expression processes (Fig. 2C). A Venn diagram analysis identified 13 targets of miR-19a-3p involved in the top 3 biological processes (Fig. 2D), including four transcription factors (ESR1, YY1, FOXP1 and REST).

Fig. 2.

Identification of potential pathways regulated by miR-19a-3p. (A) Identification of potential targets of miR-19a-3p using miRDB and TargetScan databases. (B) Pathways analysis of 829 potential targets of miR-19a-3p using KEGG Human 2021 library in Enrichr database. (C) BP analysis of 829 potential targets of miR-19a-3p using GO Biological Processes 2021 library in Enrichr database. (D) Identification of genes in common in the top three BP through a Venn diagram.

Together, our results suggest that miR-19a-3p overexpression may inhibit the expression of its target mRNAs, affecting key pathways and processes to promote CC progression.

Previous studies have reported that miR-19a-3p interacts with certain lncRNAs, such as XIST [48], LINC00094 [49], and KCNQ1OT1 [50]. Low levels of BDNF-AS expression in CC have also been reported, as has the tumor suppressor function of this lncRNA in this cancer type [51]. However, it remains unknown as to whether miR-19a-3p regulates BDNF-AS expression in CC. We therefore focused on this potential regulatory axis in CC.

To explore whether miR-19a-3p is involved in BDNF-AS regulation in CC, we searched potential binding sites for this miRNA in BDNF-AS using the Ensembl database. In total, three potential binding sites for this miRNA were located in intron 7, exon 8, and the 3 UTR^′ (Fig. 3A) of BDNF-AS. In addition, RT-qPCR revealed an inverse association between the expression patterns of these two ncRNAs in SiHa and HeLa CC cells (Fig. 3B). Moreover, we observed a negative correlation between miR-19a-3p and BDNF-AS expressions in CC tissue samples from the TCGA dataset using the ENCORI database (r = –0.1777, p-value = 1.93 $\times$ 10^-3) (Fig. 3C).

Fig. 3.

Identification of BDNF-AS as target potential of miR-19a-3p in CC. (A) Identification of binding sites for miR-19a-3p in intron, exon and 3^′ UTR regions in BDNF-AS using Ensembl database. (B) Expression analysis of miR-19a-3p in CC cell lines through RT-qPCR assays. Expression levels of miR-19a-3p and BDNF-AS were determined in HeLa and SiHa cells and they are relative to C-33A cells (dotted black line). *p-value 0.05; ***p-value 0.001. (C) Correlation analysis between miR-19a-3p and BDNF-AS expression in CC tissue samples from the TCGA dataset through ENCORI database.

These results suggest that miR-19a-3p may regulate the expression of BDNF-AS in CC, which could be due to their potentially direct interaction.

To establish the expression of BDNF-AS in CC, we analyzed CC and normal tissue samples from the TCGA, GSE67522, and GSE63514 datasets using the OncoDB and GEO databases. The results revealed that BDNF-AS expression is decreased in CC tissue samples (Fig. 4A–C).

Fig. 4.

Expression, diagnostic value, and prognostic value of BDNF-AS in CC. (A–C) BDNF-AS expression in CC and normal tissue samples from the TCGA, GSE67522, and GSE63514 datasets was assessed using the OncoDB and GEO databases. (D,E) The diagnostic value of BDNF-AS expression in CC and normal tissue samples from the GSE67522 and GSE63514 datasets was assessed using ROC curves in GraphPad Prism. (F) The prognostic value of BDNF-AS expression in CC tissue samples from the TCGA dataset was assessed using a Kaplan-Meier curve using the Kaplan-Meier Plotter database. Red line: high BDNF-AS expression; Black line: low BDNF-AS expression. *p-value 0.05; ****p-value 0.0001.

To investigate the diagnostic value of BDNF-AS expression in CC, a ROC curve analysis was performed based on the expression of this lncRNA in CC and normal samples from the GSE67522 and GSE63514 datasets using GraphPad Prism. The results revealed an AUC of 0.7170 (p-value: 0.0162) and 0.6920 (p-value: 0.0179), with sensitivities of 65.00% and 60.71%, respectively, and corresponding specificities of 68.18% and 66.67% (cutoff: 4.931 and 1.881), respectively (Fig. 4D,E). In addition, the prognostic value of BDNF-AS expression in CC was analyzed using a Kaplan-Meier curve in CC tissue samples from the TCGA dataset using the KM-Plotter database. Notably, low BDNF-AS expression was found to be associated with poor OS (Fig. 4F).

These results suggest that BDNF-AS expression is decreased in CC and may serve as a potential diagnostic and prognostic biomarker.

To identify the potential pathways regulated by BDNF-AS in CC, we first identified the genes correlated with BDNF-AS expression in CC tissue samples from the TCGA dataset using the cBioportal database. BDNF-AS expression was positively correlated with the expression of 826 genes (Supplementary Table 1), while it was negatively correlated with the expression of 306 genes (Supplementary Table 2). These 1132 genes were further analyzed using the KEGG 2021 Human library in the Enrichr database. Interestingly this approach identified key pathways in CC, such as the Glycolysis/Gluconeogenesis, Cellular senescence, Tumor Necrosis Factor (TNF) signaling, Hedgehog signaling, Proteoglycans in cancer, and Central carbon metabolism in cancer pathways (Fig. 5A). In addition, biological process analysis revealed that BDNF-AS may be involved in processes including the following: Positive regulation of protein localization to cajal body, Ribonucleoprotein complex biogenesis, “De novo” Post-translational protein folding and Ribosome biogenesis (Fig. 5B). Finally, the 1132 genes were used to construct a PPI network (Fig. 5C) with 1054 nodes, 204 edges, an average node degree of 0.387, an average local clustering coefficient of 0.102, an expected number of edges of 110, and a PPI enrichment p-value of 1.11 $\times$ 10^-15 (Fig. 5D).

Fig. 5.

Pathways and biological processes of genes that correlate with BDNF-AS in CC. (A) Pathways associated with the 1132 genes that were negatively or positively correlated with BDNF-AS expression in CC tissue samples, as assessed using the KEGG 2021 Human library in the Enrichr database. (B) The top biological process terms associated with the 1132 genes that were negatively or positively correlated with BDNF-AS expression, as assessed using the Biological Process 2021 library in the Enrichr database. (C,D) Construction and analysis of a PPI network incorporating the 1132 genes correlated with BDNF-AS expression in CC tissue samples.

In addition, we identified DEGs when comparing the top 30 CC patients with the lowest levels of BDNF-AS expression and the top 30 patients with the highest levels of BDNF-AS expression (Fig. 6A). For this comparison (low vs high BDNF-AS expression), 1796 DEGs were identified, including 1125 upregulated genes and 671 downregulated genes (Fig. 6B). GSEA analysis revealed that these DEGs were associated with key CC-related pathways, such as the apoptosis, mTORC1 signaling, and IL6/JAK/STAT3 signaling pathways. In patients with CC and high BDNF-AS expression, these pathways were found to be suppressed as compared with CC patients and exhibiting low levels of BDNF-AS expression (Fig. 6C).

Fig. 6.

Enriched pathways in CC patients with high BDNF-AS expression. (A) BDNF-AS expression in CC patients with high BDNF-AS expression (top 30) and low BDNF-AS expression (top 30) from the TCGA dataset. (B) DEGs were identified when comparing those CC patients with low and high levels of BDNF-AS expression from the TCGA dataset. (C) Enriched pathways derived from GSEA analysis of the identified DEGs. ****p-value 0.0001.

These results suggest that BDNF-AS may be involved in the negative regulation of genes involved in key pathways related to CC pathogenesis.