Blood Plasma Small Non-Coding RNAs as Diagnostic Molecules of Progesterone Receptor-Negative Phenotype of Serous Ovarian Tumors

1. Introduction

According to world statistics, ovarian cancer ranks seventh in cancer mortality among women [1]. Despite ovarian cancer mortality has declined by more than 30% over the past 50 years due to improvements in treatment, survival rate is still fewer than 50% at 5 years after diagnosis [2]. The main methods of treatment of advanced ovarian cancer are a combination of surgery and chemotherapy (paclitaxel and platinum drugs in combination with bevacizumab or PARP inhibitors), the optimal order of which has not yet been determined. Cancer antigen 125 (CA125) is commonly used as a marker of ovarian cancer, but elevated levels are detected in only 50% of disease stage I and 80% of disease stage III-IV, and should be interpreted in conjunction with clinical signs and ultrasound findings due to nonspecificity and occurrence in other diseases [3,4,5,6]. It is used to monitor the patient’s response to neo-adjuvant and adjuvant treatments [7,8], or to predict the overall survival probability at 3 months after the end of primary treatment: risk of death is 51% at 24 months and up to 79% at 60 months if a CA125 value is above 35 [9]. Significantly improvement diagnosis of early stage invasive serous ovarian/tubal carcinoma was the result of the two fully completed largest population clinical trials, using multimodal screening strategy (measurement serum level CA 125 and transvaginal ultrasound) [10,11]. But there was no significant reduction in ovarian cancer mortality compared to no screening cohort.

Different researcher teams are still elaborating new strategies ovarian cancer screening selecting biomarkers alone or in combination with CA 125, such as ROMA, CPH1, OVA1, Overa [12,13]. These approaches were developed to improve detecting early-stage disease aimed to reduce ovarian cancer mortality. But the aggressive behavior of a tumor is caused by its biological properties (extensive stromal reaction and increased invasiveness), which determine the failure of cytoreductive surgery, as well as сhemoresistance [14]. Recent studies revealed that hormone receptor status defines tumor invasive properties and longevity, in particular, low levels of progesterone receptor (PGR) expression are associated with a more aggressive disease course and worse outcome in LGSOC and HGSOC [15,16], endometrioid carcinoma [17]. On the contrary, improved survival rate in patients with poorly differentiated epithelial ovarian tumors was associated with high serum progesterone in combination with expression of PGR [18]. The protective effect of progesterone may be due in part to PGR mediated suppression of progesterone receptor membrane component-1 providing the sensitivity of ovarian cancer cells to platinum-based chemotherapy [19].

The main regulators of signalling pathways in a cell through gene expression level control at the transcriptional and post-transcriptional levels are small non-coding RNAs, including microRNAs (miRNAs) and piwiRNA (piRNA). Changed profiles of many miRNAs, implicated in pathogenesis of the gynecological diseases, have been identified [20,21] and showed a clear histotype specific pattern [22]. Role of piRNAs in diverse kinds of cancers have been also demonstrated [23]. Recently, PIWI proteins and piRNAs play a prometastatic role in ovarian carcinoma, contributing to disease progression, and are considered as potential diagnostic and prognostic biomarkers for ovarian cancer [24,25].

Therefore, in the present study, it was of interest to find the relationship between the level of progesterone receptor expression in serous tumors and small non-coding RNAs, which are potential regulators of CA125 and epithelial-mesenchymal transition, forming an aggressive tumor phenotype and its chemoresistance. It was important to develop a noninvasive method of diagnosing progesterone receptor-negative tumor for the proper management of patients with this type of tumor.

2. Results

2.1. Analysis of tumor-specific miRNAs that regulate the level of CA125 in blood serum

In order to identify potential mechanisms of changes in CA125 secretion levels in patients with serous ovarian tumors, we compared miRNA sequencing data from tumor tissue (Table S1) and miRWalk data (http://mirwalk.umm.uni-heidelberg.de/). Comparison of the miRNA expression pattern in HGSOC tissue samples from patients P25, P26, P28 (Table 1) with benign cystadenoma tissue samples from patients P4, P6, P7 (Table 1) revealed 144 differentially expressed miRNAs (Table S1, Sheet 1), among which 64 miRNAs were significantly downregulated (Table S1, Sheet 2). According to miRWalk, 31 of 64 miRNAs (Table S1, Sheet 3) are potential regulators of expression levels of mucin 16 (MUC16), matrix metalloproteinase 7 (MMP7) and matrix proteinase 9 (MMP9). At that, MMP7 and MMP9 presumably detach from MUC16 its extracellular domain, CA125 [26]. Moreover, increased expression/activity of MMP7 and MMP9 provides conditions for metastasis of ovarian cancer cells by remodeling the extracellular matrix and enhancing their migration and attachment to secondary sites [27]. In addition, MMP-7 upregulation can occur due to the interaction of the carboxy-terminal portion of the MUC16/CA125 protein and mesothelin causing an increase of invasive tumor properties [28,29]. The levels of the following miRNAs were analyzed in FFPE sections of tumor tissue from 38 patients (Table 1) by real-time quantitative PCR (Table S5): hsa-miR-199a-5p, hsa-miR-424-3p, hsa-miR-424-5p, hsa-miR-134-5p, hsa-miR-214-3p, hsa-miR-125b-5p, hsa-miR-139-5p.

The relative miRNA expression level was calculated from the difference between the threshold cDNA amplification cycles (Сt) of the analyzed miRNA and the reference endogenous SNORD68 (Figure 1). Analysis of the significance of differences in miRNA expression level (-∆Ct) in the compared groups was performed using a two-sided Wilcoxon-Mann-Whitney test (Table 2).

It follows from Figure 1 and Table 1 that in malignant ovarian neoplasms (SBT, LGSOC and HGSOC) the median expression level of analyzed miRNAs was lower than that in benign serous ovarian tumor (BSC), which is consistent with sequencing data (Table S1). A significant downregulation of miR-214-3p and miR-424-5p, potentially regulating MMP7 and/or MUC16, in the groups of serous ovarian carcinomas (LGSOC and HGSOC) was observed. Potential additional regulators of the expression level of MUC16 and MMP7 in the HGSOC group are miR-125b-5p and miR-199a-5p, whose expression level was significantly reduced relative to the BSC group. It should be noted that only the LGSOC group compared to the BSC group showed a significant decrease of miR-134-5p and miR-424-3p expression level, which potentially regulate the expression level of target genes MUC16 and MMP9. There were no significant differences between the SBT and the BSC groups in the expression levels of all analyzed miRNAs.

2.2. Analysis of MMP7, MMP9, MUC16, and PGR gene expression levels in the serous ovarian tumors

The level of expression of the MMP7, MMP9, MUC16 and PGR genes was analyzed in 38 samples of FFPE sections of serous tumors BSC, SBT, LGSOC, HGSOC (Table 1) by quantitative real-time PCR, using GAPDH, TUBA1B, ACTNB as reference genes (Figure 2, Table 3and Table S5). For MMP7, MMP9, and MUC16 mRNAs, there is a trend toward increased expression level, and for PGR mRNA – a trend toward decreased expression levels in malignant serous tumors (SBT, LGSOC, HGSOC) compared to benign serous tumors (BSC). A significant increase of MMP7, MMP9, and MUC16 gene expression levels was observed in the serous ovarian carcinoma groups (LGSOC and HGSOC) relative to the BSC group (Table 3). A statistically significant increase in the expression level of MMP7 and MUC16 was detected in the SBT group relative to the BSC group (Table 3). Multidirectional significant changes in the expression level of the MMP7 and MUC16 genes and their potential regulators miR-214-3p, miR-424-5p, miR-125b-5p and miR-199a-5p were revealed in the HGSOC group (Figure 1 and Figure 2). In the LGSOC group, multidirectional significant changes in the expression level of the MMP7, MMP9, and MUC16 genes and their potential regulators miR-214-3p, miR-424-5p, miR-134-5p, and miR-424-3p were revealed. It is important to note that along with a significant increase in the expression level of MMP7, MMP9 and MUC16, there is a significant decrease in the expression level of PGR in the LGSOC and HGSOC samples (Figure 2), that may provide the invasive properties of the tumor and the ability to metastasize.

2.3. Analysis of microRNA-regulators of EMT in the serous ovarian tumors tissue of

In a previous study [16], we detected differential expression of EMT-associated miR-16-5p, miR-17-5p, miR-20a-5p, miR-93-5p in the blood plasma of PGR-negative HGSOC patients. The present study aimed to analyze the expression level of these miRNAs in FFPE sections of serous tumors as a function of PGR expression level. The relative level of miRNA expression was calculated from the difference between the threshold cDNA amplification cycles of the analyzed miRNA and the reference endogenous SNORD68 (Figure 3).

Figure 3 and Table 4 show that in malignant ovarian tumors (SBT, LGSOC and HGSOC) the median expression level of miRNAs that potentially regulate EMT was higher than that in benign serous ovarian tumors (BSC), which is consistent with miRNA sequencing data in tissues of serous tumors (Table S1, Sheet 1 and Sheet 4). A significant increase in the expression level of miR-16-5p, miR-17-5p, miR-20a-5p, and miR-93-5p was found in the HGSOC group relative to the BSC group, which, along with a significant decrease in the PGR expression level and an increase in the expression level MMP7 and MUC16, potentially regulated by miR-214-3p, miR-424-5p, miR-125b-5p and miR-199a-5p with decreased levels, indicates a more aggressive behavior of the HGSOC in terms of the ability to metastasize compared to SBT and LGSOC. In the LGSOC group, a significant increase in miR-16-5p expression was found compared to BSC. In the SBT group, there were no significant changes in the analyzed miRNAs responsible for EMT.

2.4. Analysis of piRNA expression in serous ovarian tumors

One of the functions of piRNA is to regulate the stability of the cell genome by interacting in the nucleus with the retrotransposon transcript as part of the RISC complex, which contacts with histone deacetylase and histone methyltransferase and DNA methyltransferase, which block further transcription of the retrotransposon, thereby preventing its activity and integration into various genome regions [30]. In addition to the suppressive activity of piRNA against transposons, their regulatory effect on various signaling pathways in the cell is also known through both target mRNA destabilization/inhibition of translation and target mRNA stabilization/activation of translation [30] [31]. Due to proven genomic instability and changes in the activity of numerous signaling pathways in cancer cells, in particular, in serous ovarian carcinomas [15,32,33,34,35,36], we analyzed piRNA expression profiles in SBT tissue from patients P13, P15, P16 (Table 1) and in HGSOC tissue from patients P25, P26, P28 (Table 1) relative to that in BSC tissue from patients P4, P6, P7 (Table 1) by deep sequencing method (Table S2, Sheet 1 and Sheet 2). Ninety-seven piRNAs and 77 piRNAs differentially expressed in the HGSOC and SBT groups, respectively, were identified compared with BSC (p < 0.1), with the resulting lists of 37 piRNAs overlapping. For further analysis of 38 samples of FFPE tissue sections of serous ovarian tumors, 19 piRNAs were selected for quantitative real-time PCR (Table S5), six of which had an altered expression level in HGSOC and SBT, and the remaining 13 piRNAs – only in HGSOC. The relative level of piRNA expression was calculated from the difference between the threshold cycles of cDNA amplification of the analyzed piRNA and reference endogenous SNORD68 (Figure 4).

Out of 19 piRNAs, a Ct value of less than 35 cycles was observed only in hsa_piR_004307, hsa_piR_009295, hsa_piR_019914, hsa_piR_020813, hsa_piR_022437, among which significant differences were found for hsa_piR_004307, hsa_piR_009295 and hsa_piR_019914 in the HGSOC group relative to BSC.

2.5. Correlation analysis of the expression level of tumor-associated miRNA, piRNA, mRNA, progesterone receptor and the level of CA125 in the blood serum of patients.

A Spearman correlation matrix was constructed to understand the potential relationships of the ovarian tumor process characteristics which form the molecular-biological portrait of serous tumors of different types (Figure 5, Table S3). The analyzed samples were ranged, according to the diagnosis (type of the serous tumor) in the following way: “BSC” < “SBT” < “LGSOC” < “HGSOC”. The expression level of hsa-miR-199a-5p, hsa-miR-214-3p, hsa-miR-424-5p, and hsa-miR-125b-5p were significantly inversely correlated with the expression level of MUC16 mRNA in the tumor tissue and with the CA125 level in the blood serum, but directly correlated with the level of PGR mRNA and PGR protein (according to Allred score). On the contrary, the expression level of miR-16-5p, miR-17-5p, miR-20a-5p, and miR-93-5p is significantly inversely correlated with the level of PGR, and directly correlated with the expression level of hsa_piR_004307, hsa_piR_009295, hsa_piR_019914, hsa_piR_020813, hsa_piR_022437 in the tumor tissue and with the CA125 level in the blood serum of patients. Among piRNAs, the hsa_piR_004307 expression level was significantly inversely correlated with the PGR protein level in the tumor tissue. MMP7 mRNA is significantly directly correlated with the level of MUC16 mRNA and MMP9 mRNA, which, in turn, is inversely correlated with the level of PGR. Note that the level of PGR mRNA is significantly correlated with the level of PGR protein, both of which are significantly inversely correlated with the level of CA125 in the blood serum of patients with serous ovarian tumors.

2.6. Partial least squares discriminant analysis (PLS-DA) of the molecular biological parameters determining certain type of tumor.

Since significant correlations were found between various molecular biological parameters depending on the type of serous ovarian tumor, it seemed interesting to evaluate the contribution of each parameter to the formation of a specific tumor type, namely: BSC, SBT, LGSOC, and HGSOC. All data obtained in sections 2.1 - 2.4 were used for partial least squares (PLS) analysis with the plotting the graph presented in Figure 6, which clearly shows the formation of clusters of samples depending on the type of serous tumor. The greatest contribution to the separation of BSC, SBT, LGSOC, and HGSOC groups is made by molecules with VIP > 1: PGR expression level, miRNAs responsible for EMT (miR-16-5p, miR-17-5p, miR- 20a-5p, and miR-93-5p), hsa_piR_004307 (potential regulator of genome stability and signaling pathways in the cell), miR-214-3p (potential regulator of MMP7 and MUC16 expression levels), CA125 level in patients' blood serum, where the PGR expression level plays a primary role.

2.7. Logistic regression models for diagnosing a progesterone receptor-negative serous ovarian tumor based on miRNA and piRNA expression levels in tumor tissue

The relationship between PGR expression level in tumor tissue and response to adjuvant chemotherapy was found, as shown in Table 1: in PGR-negative LGSOC and HGSOC, there was disease stabilization/progression in 75% and 100% cases, respectively. Due to the presence of correlations between the expression level of PGR and the expression level of small non-coding RNAs (miRNA and piRNA) as well as protein-coding mRNAs responsible for biological properties of a certain type of serous ovarian tumor, it was decided to develop a logistic regression model for tumor identification by PGR-phenotype. Based on the "-ΔCt" values obtained in Sections 2.1–2.4, logistic regression models were developed in RStudio program by finding the optimal combination of predictor variables with their stepwise inclusion and exclusion according to the contribution to the model and their significance, where the PGR expression level was dependent variable (0 - PGR-positive tumor, 1 - PGR-negative tumor in accordance with Allred score).

Table 6. Parameters of Figure 7 logistic regression models.

Coefficient	coefficient value (95% CI)	Wald test	p_value	Odds ratio (95% CI)	sensitivity	specificity
model 1					0.7857	0.9167
(Intercept)	6.549(2.224;12.346)	2.618	0.008	698.816(9.245;230163.148)
piR_020813	-1.504(-3.151;-0.316)	-2.144	0.032	0.222(0.042;0.728)
piR_004307	1.142(0.198;2.475)	2.045	0.04	3.136(1.219;11.885)
miR-93-5p	1.037(0.34;1.935)	2.614	0.008	2.823(1.405;6.928)
model 2					0.9286	0.7083
(Intercept)	6.701(2.428;13.488)	2.429	0.015	813.662(11.341;720717.037)
piR_020813	-1.807(-3.818;-0.415)	-2.127	0.033	0.164(0.021;0.66)
miR-16-5p	1.803(0.678;3.517)	2.56	0.01	6.073(1.97;33.695)
piR_022437	1.248(0.259;2.713)	2.063	0.039	3.483(1.296;15.081)
model 3					0.9286	0.75
(Intercept)	6.803(2.115;13.645)	2.397	0.016	900.756(8.292;844027.632)
piR_004307	2.748(0.756;5.757)	2.23	0.025	15.613(2.13;316.434)
piR_019914	-2.168(-4.702;-0.446)	-2.046	0.04	0.114(0.009;0.639)
miR-93-5p	0.881(0.209;1.811)	2.248	0.024	2.415(1.233;6.12)
model 4					0.8571	0.875
(Intercept)	4.319(1.219;8.825)	2.29	0.021	75.133(3.384;6806.095)
piR_020813	-1.541(-3.308;-0.29)	-2.051	0.04	0.214(0.036;0.747)
piR_004307	1.129(0.129;2.466)	1.951	0.051	3.093(1.138;11.782)
miR-16-5p	1.327(0.374;2.65)	2.364	0.018	3.773(1.454;14.161)
model 5					0.8571	0.75
(Intercept)	7.585(2.542;14.208)	2.603	0.009	1969.531(12.709;1481859.036)
piR_009295	-1.37(-2.869;-0.194)	-2.061	0.039	0.254(0.056;0.823)
piR_004307	1.901(0.452;3.84)	2.276	0.022	6.698(1.572;46.542)
miR-93-5p	0.772(0.174;1.53)	2.3	0.021	2.164(1.191;4.619)
model 6					0.7857	0.875
(Intercept)	6.903(2.61;12.621)	2.769	0.005	995.67(13.604;302894.002)
piR_020813	-1.077(-2.305;-0.093)	-1.965	0.049	0.34(0.099;0.91)
miR-93-5p	1.071(0.382;1.977)	2.697	0.006	2.92(1.466;7.221)
piR_022437	0.72(0.006;1.611)	1.819	0.068	2.054(1.006;5.008)

Table 6. Parameters of Figure 7 logistic regression models.

Coefficient	coefficient value (95% CI)	Wald test	p_value	Odds ratio (95% CI)	sensitivity	specificity
model 1					0.7857	0.9167
(Intercept)	6.549(2.224;12.346)	2.618	0.008	698.816(9.245;230163.148)
piR_020813	-1.504(-3.151;-0.316)	-2.144	0.032	0.222(0.042;0.728)
piR_004307	1.142(0.198;2.475)	2.045	0.04	3.136(1.219;11.885)
miR-93-5p	1.037(0.34;1.935)	2.614	0.008	2.823(1.405;6.928)
model 2					0.9286	0.7083
(Intercept)	6.701(2.428;13.488)	2.429	0.015	813.662(11.341;720717.037)
piR_020813	-1.807(-3.818;-0.415)	-2.127	0.033	0.164(0.021;0.66)
miR-16-5p	1.803(0.678;3.517)	2.56	0.01	6.073(1.97;33.695)
piR_022437	1.248(0.259;2.713)	2.063	0.039	3.483(1.296;15.081)
model 3					0.9286	0.75
(Intercept)	6.803(2.115;13.645)	2.397	0.016	900.756(8.292;844027.632)
piR_004307	2.748(0.756;5.757)	2.23	0.025	15.613(2.13;316.434)
piR_019914	-2.168(-4.702;-0.446)	-2.046	0.04	0.114(0.009;0.639)
miR-93-5p	0.881(0.209;1.811)	2.248	0.024	2.415(1.233;6.12)
model 4					0.8571	0.875
(Intercept)	4.319(1.219;8.825)	2.29	0.021	75.133(3.384;6806.095)
piR_020813	-1.541(-3.308;-0.29)	-2.051	0.04	0.214(0.036;0.747)
piR_004307	1.129(0.129;2.466)	1.951	0.051	3.093(1.138;11.782)
miR-16-5p	1.327(0.374;2.65)	2.364	0.018	3.773(1.454;14.161)
model 5					0.8571	0.75
(Intercept)	7.585(2.542;14.208)	2.603	0.009	1969.531(12.709;1481859.036)
piR_009295	-1.37(-2.869;-0.194)	-2.061	0.039	0.254(0.056;0.823)
piR_004307	1.901(0.452;3.84)	2.276	0.022	6.698(1.572;46.542)
miR-93-5p	0.772(0.174;1.53)	2.3	0.021	2.164(1.191;4.619)
model 6					0.7857	0.875
(Intercept)	6.903(2.61;12.621)	2.769	0.005	995.67(13.604;302894.002)
piR_020813	-1.077(-2.305;-0.093)	-1.965	0.049	0.34(0.099;0.91)
miR-93-5p	1.071(0.382;1.977)	2.697	0.006	2.92(1.466;7.221)
piR_022437	0.72(0.006;1.611)	1.819	0.068	2.054(1.006;5.008)

Models 2 and 3 in Figure 7 had the highest sensitivity (93%) in identifying PGR-negative tumors.

2.8. Logistic regression models for the diagnosis of progesterone receptor-negative serous ovarian tumor by the level of tumor-associated mшRNAs and piRNAs circulating in the blood of patients

To diagnose PGR-negative serous ovarian tumors before surgical and chemotherapeutic treatment, the content of miRNA and piRNA, forming models 2 and 3 of Figure 7, was analyzed by real-time RT-PCR in the blood plasma of 38 patients (ID number in Table 1). "-ΔCt" values were obtained when using hsa_piR_004308 as a reference molecule to quantify hsa_piR_020813, hsa_piR_022437, hsa_piR_004307, hsa_piR_019914, and when using hsa-miR-30d-5p as a reference molecule to analyze hsa-miR-16-5p and hsa-miR-93-5p levels. Logistic regression models were developed based on the obtained data and presented in Figure 8. Parameters of the developed models are presented in Table 7.

Model 1 in Figure 8 has the highest sensitivity (86%) in identifying PGR-negative tumors. The formula describing this model is presented below:

$$\frac{1}{1+e^{4.3-1.29x_1-0.75x_2}}$$

where x1 – “-ΔCt” for hsa_ piR_022437, x2 – “-ΔCt” for hsa_miR-16-5p.

2.9. Functional Significance of RNA-markers associated with PGR-negative serous ovarian tumor phenotype

Potential targets of piRNAs from Figure 7 and Figure 8 (piR_020813, piR_004307, piR_022437, piR_019914, piR_009295) were predicted as described in our recent manuscript [37]. The list of RNA targets for these piRNAs in the form of RefSeq mRNA accessions is presented in Table S4 (Sheet 1), followed by conversion to gene symbols using the bioDBnet database (https://biodbnet-abcc.ncifcrf.gov/db/db2db.php, last accessed on 15 Marсh 2023). Potential target mRNAs for hsa-miR-16-5p and hsa-miR-93-5p were identified using miRtargetlink database (https://ccb-web.cs.uni-saarland.de/mirtargetlink/, last accessed on 15 Marсh 2023). The total list of miRNA and piRNA gene-targets associated with the PGR-negative serous ovarian tumor phenotype is presented in Table S4 (Sheet 2). Functional significance of the target genes was assessed by the FunRich3.1.3 (http://www.funrich.org/download) functional enrichment analysis tool (Table, Sheet 3—Sheet 12). We found that among 3166 gene-targets for miRNA and piRNA from Figure 7 and Figure 8 logistic regression models (Table S4, Sheet 2), 1421 genes are implicated in ovarian cancer (Table S4, Sheet8) and 101 genes that are involved in the pathogenesis of cancer when mutations occur in them according to the Cancer Gene Census from COSMIC database (Table S4, Sheet 11). These three gene lists were compared by Venny 2.1 (https://bioinfogp.cnb.csic.es/tools/venny/), the common one being a list of 72 genes (Figure 9A, Table S4, Sheet 12). Analysis of the expression site of 72 genes in the FunRich3.1.3 tool revealed experimentally proven involvement of these genes not only in ovarian cancer but also in breast, colorectal, liver, stomach, cervical, lung, endometrial, pancreatic, and thyroid cancers (Figure 9B).

3. Discussion

Ovarian cancer is the most lethal gynecologic cancer with a 5-year survival rate of less than 50 percent [2]. The biological properties of the tumor (increased invasiveness, profound chromosomal instability in the cancer cells, and chemoresistance) dictated by the tumor microenvironment, especially cancer-associated fibroblasts, natural killer (NK) cells, Th2 cells [14,32,35,36]. It was revealed that hormone receptor status defines tumor invasive properties [15,18,38,39,40]. Increased expression levels of estrogen and progesterone receptors provides the better outcome for patients with ovarian cancer than reduced hormone receptors levels. In particular, aggressive form of the disease is characterized by low levels of PGR expression and worse outcome. Moreover, several polymorphisms in the hormone-binding domain in the PGR gene has been associated with increased risk of ovarian cancer [41,42,43]. The protective effect of progesterone may be due to PGR mediated suppression of progesterone receptor membrane component-1 (PGRMC1), increasing the sensitivity of ovarian cancer cells to platinum-based chemotherapy [19]. Thereby, targeted depletion of PGRMC1 could be used as an additional therapy to cisplatin. Recently, chemopreventive effect of synthetic progestin Norethindrone was demonstrated in epithelial ovarian cancer cells SKOV3 by upregulation of TP53 expression and downregulation of VEGF, HIF-1α, COX-2, and PGRMC1 expression, leading to significantly reduced SKOV3 cell growth, increased apoptosis and necrosis, and inhibiting cell migration [44]. The present study was aimed to explore and evaluate the candidate factors, namely miRNA, piRNA, mRNA, associated with the PGR-negative serous ovarian tumors, to develop a non-invasive test for identification this aggressive and chemoresistant tumor phenotype before any treatment and to and be prepared to use targeted therapies in addition to conventional treatment regimens.

We found a significant inverse correlation between the PGR expression level in ovarian tumor tissues and the CA125 serum level, which in turn significantly inversely correlated with the expression level of its potential regulators miR-199a-5p, miR-214-3p, miR-424-3p, miR-424-5p, miR-125b-5p, whose target genes are MUC16 and/or MMP7 and/or MMP9 according to miRWalk database. It was found that the downregulation of miR-199a-5p expression might be implicated in ovarian cancer progression in connection with the revealed negative correlation with tumor infiltration, tumor size, lymphatic metastasis and TNM stage in ovarian cancer [45]. Moreover, downregulated expression of miR-199a-5p was found in ascites-derived spheroids of the primary tumour site origin which form new metastatic niches due to high invasive capability [46]. Mir-214-3p has a suppressive effect on CDK6 [47] and on MAPK1 [48], and mir-214-3p downregulation facilitates the cell-cycle progression, proliferation, migration, and invasion of ovarian cancer cells. miR-424-3p has a proven capability to sensitize ovarian cancer cells to cisplatin by decreasing the expression of an anti-apoptotic protein galectin-3; in the case of galectin-3 overexpression, chemoresistance occurs [49]. In ovary cancer cells, miR-424/503 cluster is silenced by DNA hypermethylation, thereby canceling suppression of the expression of kinesin family member 23 by miR-424-5p and promoting cell proliferation and migration [50]. It is important to note that kinesins play an important regulatory role in the formation of spindles, separation of chromosomes, and cytokinesis, and in case of abnormal expression/function of kinesins, daughter cells become aneuploidic, thereby resulting in tumorigenesis [51]. In particular, elevated levels of KIF23 in ovarian, breast, lung cancer are associated with adverse outcomes [52,53,54]. miR-125b-5p may be considered as a platinum chemoresistance marker in connection with its strong downregulation in the tumorspheres originated from the SKOV-3 platinum-resistant cell line [55]. This relationship between chemoresistance and expression level of miR-125b-5p may be explained by direct targeting of BCL2 mRNA by miR-125b-5p thereby increasing the sensitivity of cancer cells to cisplatin treatment as it has been demonstrated for gallbladder cancer where low miR-125b-5p expression and high expression of Bcl2 was correlated with poor prognosis [56].

In the present study we found that MMP7 mRNA level was significantly directly correlated with the level of MUC16 mRNA and MMP9 mRNA, which, in turn, was inversely correlated with the level of PGR. We have found a direct correlation between the PGR and the miR-199a-5p, miR-214-3p, miR-424-3p, miR-424-5p, miR-125b-5p expression levels, and in the case of downregulation of these miRNAs in a PGR-negative serous tumor, the levels of MUC16, MMP7, and MMP9 can be increased, enhancing the ability of ovarian cancer cells to migrate, adhere to secondary sites and form metastasis [27], and providing tumor chemoresistance as discussed above. The maintenance of MMP7 concentration at a high level may be due to the induction of MMP-7 expression via a p38 mitogen-activated protein kinase (MAPK)-dependent pathwayas the result of the interaction of the carboxy-terminal portion of the MUC16/CA125 protein and mesothelin present on mesothelial cells lining the peritoneum [28,29]. Because high mesothelin correlates with chemoresistance and poor survival in epithelial ovarian carcinoma [57], several clinical trials are underway to evaluate the safety and efficacy of mesothelin-targeted drugs in platinum-resistant ovarian cancer [58,59].

In the present study, the tissue content of miR-16-5p, miR-17-5p, miR-20a-5p, miR-93-5p, responsible for epithelial-mesenchymal transition (EMT) of the cell [60,61,62,63,64], significantly directly correlated with the blood serum CA125 concentration and significantly inversely correlated with the PGR expression level in the tumor tissue. Found here increased levels of miR-20a-5p, miR-16-5p, miR-17-5p and miR-93-5p in tumor tissues and and peripheral blood from patients with HGSOC serous ovarian cancer are in good agreement with the literature data on the same direction of changes in the biological samples from patients with serous ovarian tumors [16,65,66,67,68]. It should be noted that a sharp significant decrease in the level of PGR mRNA was found in the present study in the LGSOC and HGSOC group, with a more pronounced drop in the case of the latter. This observation, along with a significant increase in the expression level of miR-20a-5p, miR-16-5p, miR-17-5p and miR-93-5p in the HGSOC group, probably indicates the most aggressive carcinoma phenotype among other types of serous tumors.

Levels of the EMT-associated miRNAs significantly directly correlated with the content of hsa_piR_022437, hsa_piR_009295, hsa_piR_020813, hsa_piR_004307, hsa_piR_019914 in serous ovarian tumor tissues. Among them, expression level of hsa_piR_004307 significantly inversly correlated with the PGR expression level in the tumor. One of the functions of piRNA is to regulate the stability of the cell genome by interacting in the nucleus with the retrotransposon transcript as part of the RISC complex, which contacts with histone deacetylase and methyltransferase, and DNA methyltransferase, which block further transcription of the retrotransposon, thereby preventing its activity and integration into different parts of the genome [30]. In addition to the suppressor activity of piRNAs against transposons, their regulatory effects on various signaling pathways in the cell through both destabilization of the target mRNA and inhibition of translation or stabilization of the target mRNA and activation of translation are also known to exist [30,31,69,70]. When analyzing potential target genes of hsa_piR_022437, hsa_piR_009295, hsa_piR_020813, hsa_piR_004307, hsa_piR_019914 and two miRNAs, miR-16-5p and miR-93-5p, involved in the identification of PGR-negative serous ovarian tumors by logistic regression method, the participation of 72 genes in pathogenesis of ovarian, breast, colorectal, liver, stomach, cervical, lung, endometrial, pancreatic, and thyroid cancers was proved, including the occurrence of mutations in these genes according to the Cancer Gene Census. Perhaps there are common pathogenetic mechanisms of epithelial cancers by altering the functional activity of this group of genes, in particular, under the action of miRNAs and piRNAs, used as prognostic parameters in developed logistic regression models.

The genome sequences, coding certain piRNA classes, are located within protein-coding genes, and functional piRNA molecules are mainly produced from the 3'-untranslated regions (3'-UTRs) of mRNA during their translation [71]. In addition, mRNAs 3'-UTR, producing piRNAs, usually contain transposon sequences whose activity is controlled by these piRNAs at the post-transcriptional level [72]. It is important to note that some piRNAs associated with the PGR-negative tumor phenotype are located in the loci of protein-coding genes and/or transposon, namely: the hsa_piR_022437 DNA sequence is related to retrotransposon SINE and located in the SUN1 gene involved in directed cell migration [73], and hsa_piR_009295 DNA sequence is located in the centrosome linker protein rootlein (encoding for CROCC gene with alternative name – TAX1BP2), whose overexpression inhibites centrosome duplication, whereas depletion results in centrosome hyperamplification [74] that is linked to oncogenesis [75,76]. Centrosome aberrations may induce the dissemination of metastatic cells and contribute to aggressive cancer subtypes [77].

After examining the tumor response to adjuvant chemotherapy, we found that in PGR-negative LGSOC and HGSOC, there was disease stabilization/progression in 75% and 100% cases, respectively. Thus, a substantial proportion of patients with primary PGR-negative LGSOC and HGSOC were resistant to platinum-based treatment that requires novel therapeutic and preventive approaches. In the present research, we developed non-invasive method based on regression analysis of hsa_piR_022437 и miR-16-5p content in the blood plasma from patients which can be used to diagnose PGR-negative tumor phenotype with 86% sensitivity before surgery and chemotherapy to choose the treatment strategy for this most aggressive type of ovarian cancer. Recently, a trial to investigate the effect of the addition of hyperthermic intraperitoneal chemotherapy (HIPEC) to interval cytoreductive surgery was conducted among patients with stage III epithelial ovarian cancer who were receiving neoadjuvant chemotherapy, and this procedure resulted in longer recurrence-free survival and overall survival than surgery alone [78]. It may be so that besides administration of HIPEC with cisplatin it is necessary to add targeted therapy aimed at depletion of PGRMC1, decreased mesothelin expression and Th2 infiltration to increase the effectiveness of treatment and longevity of the patient with highly aggressive serous ovarian tumor as PGR-negative type.

4. Materials and Methods

4.1. Patients enrolled in the study

38 women enrolled in the study aged between 29 and 71-years-old with primary ovarian tumors were referred to the National Medical Research Center for Obstetrics, Gynecology, and Perinatology, named after Academician V.I. Kulakov of the Ministry of Healthcare of the Russian Federation for clinical and instrumental additional examination and surgical intervention in the volume depending on the stage of the disease, which was histologically verified according to FIGO. After adjuvant chemotherapy (carboplatin AUC 6 + paclitaxel 175 mg/m2), tumor response was evaluated according to the RECIST 1.1 criteria. The following groups were formed: benign serous cystadenoma, n = 10; borderline serous cystadenoma, n = 6; low-grade serous ovary cancer, n = 8; high-grade serous ovary cancer, n = 14.

4.2. RNA Isolation from Peripheral Blood Plasma

S-MONOVETTE tubes containing EDTA KE (Sarstedt AG&Co., Ltd., Nümbrecht, Germany, cat. No. 04.1915.100) were used to sample venous blood from patients entrolled in the study. 200 µL of blood plasma collected after two-step centrifugation for 20 min at 300× g (4 °C) and for 10 min at 16,000× g, were used for RNA extraction applying an miRNeasy Serum/Plasma Kit (Qiagen, Germany, cat. No. 217184).

4.3. RNA isolation from ovarian tumors

Ovarian tumor samples were collected during surgery and immediately frozen in liquid nitrogen or embedded in paraffin blocks after fixation with neutral formalin, for subsequent total RNA extraction in the former case using the miRNeasy Micro Kit (Qiagen, Hilden, Germany, catalog No. 217084), followed by the RNeasy MinElute Cleanup Kit (Qiagen, Germany, catalog No. 74204) or in the latter case using deparaffinization solution (Qiagen, Hilden, Germany, catalog No.19093) and RNeasy FFPE Kit (Qiagen, Hilden, Germany, catalog No. 73504). Qubit fluorometer 3.0 (Life Technologies, Petaling Jaya, Malaysia, cat.Q33216) was used for RNA concentration measurement. Total RNA quality was examined on the Agilent Bioanalyzer 2100 (Agilent, Waldbronn, Germany, cat. No G2939A) using the RNA 6000 Nano Kit (Agilent Technologies, Santa Clara, CA, USA, cat. No. 5067-1511).

4.4. Small RNA Deep Sequencing

500 ng of total RNA from frozen tumor tissues were used for cDNA libraries synthesis applying the NEBNext® Multiplex Small RNA Library Prep Set for Illumina® (Set11 and Set2, New England Biolab®, Frankfurt am Main, Germany, cat. No. E7300S, E7580S). After amplification for 14 PCR cycles and purification in the 6% polyacrylamide gel, cDNA libraries were sequenced on the NextSeq 500 platform (Illumina, San Diego, AC, USA, cat. No. SY-415-1001). Deep sequencing data were processed as described in our previous publication [37], using Cutadapt to remove adapters, bowtie aligner [79] to map all trimmed reads in the range of 16 bp - 50 bp to the GRCh38.p15 human genomes, miRBase v21, and piRNABase, featureCount tool from the Subread package [80] to count aligned reads, and the DESeq2 package [81] to carry out differential expression analysis of the sncRNA.

4.5. Reverse Transcription and Quantitative Real-Time PCR of small noncoding RNA

Seven microliters of total RNA obtained in 4.2. or 250 ng of total RNA from FFPE samples obtained in 4.3. were converted into cDNA in accordance with the miScript® II RT Kit protocol (Qiagen, Germany, cat. No. 218161). After completion of the reaction and dilution of the sample by 20 times, cDNA (2 µL) was amplified during real-time PCR using a forward primer specific for the studied RNA (Table S5) and the miScript SYBR Green PCR Kit (Qiagen, Germany, cat. No. 218075). The following PCR conditions for miRNA and piRNA amplification were used: (1) 15 min at 95 °C and (2) 40 cycles at 94 °C for 15 s, an optimized annealing temperature (46.2–62 °C) for 30 s and 70 °C at 30 s in a StepOnePlusTM thermocycler (Applied Biosystems, Waltham, MA, USA, cat. No. 4376600). miR-30d-5p was used as the reference RNA to quantify miRNA and hsa_piR_004308 was used as the reference RNA to quantify piRNA in the blood plasma sample by the ∆Ct method. The relative expression of miRNA and piRNA in the FFPE ovarian tumor sections was determined by the ∆Ct method using SNORD68 as the reference RNA.

4.6. Reverse Transcription and Quantitative Real-Time PCR of mRNA

125 ng of total RNA from FFPE ovarian tumor sections, were converted into cDNA in a reaction mixture (25 µL) containing 10 µM random hexameric primer (Evrogen, Russia), 1x M-MLV RT buffer (M531A, Promega), 1x dNTP mix (0.2 mM each, Evrogen, Russia), 200 U M-MLV reverse transcriptase (M1708, Promega) at 37 ℃ over 60 min, followed by incubation at 95 ℃ over 10 min; then, the sample volume was adjusted with deionized water to 100 µL. The synthesized cDNA (2 µL) was used as a template for real-time PCR in a reaction mixture (20 µL) containing 150 nM each of the forward and reverse primers specific for the studied mRNA (Table S5) in a 1x qPCRmix-HS SYBR+HighROX (Evrogen, Russia). The following PCR conditions were used: (1) 5 min at 95 °C and (2) 40 cycles at 95 °C for 20 s, an optimized annealing temperature (48.9 – 63.8 °C) for 20 s and 72 °C at 30 s in a CFX96 Real-Time System (C1000 Touch Thermal Cycler plus CFX96 Optics Module, BioRad, Singapore). The relative expression of mRNA was determined by the ∆Ct method using geometric mean of ACTB, TUBA, GAPDH as the reference RNAs.

4.7. Immunohistochemistry

PGR immunohistochemical staining of formalin-fixed paraffin-embedded specimens of serous ovarian tumors was performed as described in our previous publication [16] according to the Allred scale [82].

4.8. Statistical Analysis of the Obtained Data

Scripts written in R language [80] and RStudio [83] were used for statistical processing as described in our previous publication [16] applying Shapiro–Wilk test, the Mann–Whitney test for paired comparison, Spearman’s nonparametric correlation test, logistic regression analysis. Study results were considered reliable if the value of statistical significance (p) was less than 0.05.