Plasma miRNAs as biomarkers for endometriosis.

STUDY QUESTION: Can plasma miRNAs be used for the non-invasive diagnosis of endometriosis in infertile women? SUMMARY ANSWER: miRNA-based diagnostic models for endometriosis failed the test of independent validation. WHAT IS KNOWN ALREADY: Circulating miRNAs have been described to be differentially expressed in patients with endometriosis compared with women without endometriosis, suggesting that they could be used for the non-invasive diagnosis of endometriosis. However, these studies have shown limited consistency or conflicting results, and no miRNA-based diagnostic test has been validated in an independent patient cohort. STUDY DESIGN, SIZE, DURATION: We performed genome-wide miRNA expression profiling by small RNA sequencing to identify a set of plasma miRNAs with discriminative potential between patients with and without endometriosis. Expression of this set of miRNAs was confirmed by RT-qPCR. Diagnostic models were built using multivariate logistic regression with stepwise feature selection. In a final step, the models were tested for validation in an independent patient cohort. PARTICIPANTS/MATERIALS, SETTINGS,
METHODS: Plasma of all patients was available in the biobank of the Leuven Endometriosis Centre of Excellence. Biomarker discovery and model development were performed in a discovery cohort of 120 patients (controls = 38, endometriosis = 82), and models were tested for validation in an independent cohort of 90 patients (controls = 30, endometriosis = 60). RNA was extracted with the miRNeasy Plasma Kit. Genome-wide miRNA expression analysis was done by small RNA sequencing using the NEBNext small RNA library prep kit and the NextSeq 500 System. cDNA synthesis and qPCR were performed using the Qiagen miScript technology. MAIN RESULTS AND THE ROLE OF CHANCE: We identified a set of 42 miRNAs with discriminative power between patients with and without endometriosis based on genome-wide miRNA expression profiling. Expression of 41 miRNAs was confirmed by RT-qPCR, and 3 diagnostic models were built. Only the model for minimal-mild endometriosis (Model 2: hsa-miR-125b-5p, hsa-miR-28-5p and hsa-miR-29a-3p) had diagnostic power above chance performance in the independent validation (AUC = 60%) with an acceptable sensitivity (78%) but poor specificity (37%). LIMITATIONS, REASONS FOR CAUTION: The diagnostic models were built and tested for validation in two patient cohorts from a single tertiary endometriosis centre. Further validation tests in large cohorts with patients from multiple endometriosis centres are needed. WIDER IMPLICATION OF THE
FINDINGS: Our study supports a possible biological link between certain miRNAs and endometriosis, but the potential of these miRNAs as clinically useful biomarkers is questionable in women with infertility. Large studies in well-described patient cohorts, with rigorous methodology for miRNA expression analysis, sufficient statistical power and an independent validation step, are necessary to answer the question of whether miRNAs can be used as diagnostics markers for endometriosis. STUDY FUNDING/COMPETING INTEREST(S): The project was funded by a grant from the Research Foundation - Flanders (FWO). A.V., D.F.O. and D.P. are PhD fellows from the FWO. T.D. is vice president and Head of Global Medical Affairs Fertility, Research and Development, Merck KGaA, Darmstadt, Germany. He is also a professor in Reproductive Medicine and Biology at the Department of Development and Regeneration, Group Biomedical Sciences, KU Leuven (University of Leuven), Belgium and an adjunct professor at the Department of Obstetrics and Gynecology in the University of Yale, New Haven, USA. Neither his corporate role nor his academic roles represent a conflict of interest with respect to the work done by him for this study. The other co-authors have no conflict of interest. TRIAL REGISTRATION NUMBER: Not applicable.

Introduction

Endometriosis, defined as the presence of endometrial-like tissue outside the uterus, is a chronic gynaecological disease with an overall prevalence of 2–10% in women of reproductive age and up to 50% in infertile women (Giudice, 2010). Patients with endometriosis can be asymptomatic but typically present with chronic pelvic pain, infertility or both. Due to its chronic character, lack of medical treatment and debilitating symptoms, endometriosis results in major public health care costs, comparable to that of Crohn’s disease (Simoens ).

Endometriosis appears as superficial peritoneal lesions, adhesions, ovarian endometriotic cysts and deeply infiltrative disease and is classified according to the staging system of the American Society for Reproductive Investigation (ASRM) into four stages (minimal, mild, moderate and severe disease) (1997). Ultrasound can detect ovarian endometriotic cysts and deep endometriotic nodules but does not rule out peritoneal endometriosis or endometriosis-associated adhesions (Nisenblat ). At present, there is no non-invasive test for endometriosis. Hence, the gold standard for diagnosis is still laparoscopic visualization of lesions confirmed by histology (Dunselman ). Several studies have reported diagnostic delays in endometriosis averaging between 4–10 years (Dunselman ). The lack of a non-invasive test is an important contributing factor to this diagnostic delay. Non-invasive diagnosis of endometriosis would allow early diagnosis and treatment, with the potential to improve quality of life and reduce the costs related to endometriosis. It has therefore been selected as a research priority by the World Endometriosis Society and the World Endometriosis Research Foundation (Rogers ).

MicroRNAs (miRNAs) are a class of small non-coding regulatory RNAs estimated to influence the translation of mRNAs in 30% of all genes in animals (Santamaria and Taylor, 2014). At present, more than 2400 human miRNAs have been identified and registered in the miRBase database (Santamaria and Taylor, 2014). miRNAs control gene expression post-transcriptionally by inhibiting translation or promoting mRNA degradation in the cytoplasm (Moreno-Moya ). miRNAs are expressed in all tissues and regulate a wide spectrum of processes, including cellular differentiation, proliferation and apoptosis (Moreno-Moya ). Although most miRNAs are localized inside the cell, a significant number of miRNAs have also been detected in extracellular body fluids such as serum, plasma, spinal fluid, follicular fluid, saliva and urine (Traver ).

Several studies have indicated that a large number of miRNAs are involved in endometriosis, showing differential expression between eutopic and ectopic endometrial tissues (Teague , Braza-Boils , Nothnick 2017, Haikalis , Panir , Rekker , Zhao ). Several of these differentially expressed miRNAs are known to regulate genes involved in processes that are crucial for the establishment and progression of endometriosis, such as angiogenesis, inflammation and immune regulation (Santamaria and Taylor, 2014). Further investigation of the expression profiles of miRNAs could lead to new insights in the pathophysiology of endometriosis and major improvements in the management of endometriosis (Santamaria and Taylor, 2014).

MiRNAs have also become a major focus of research aimed at identifying new nucleic acid–based biomarkers for human disease. MiRNAs present in serum or plasma are bound to (lipo) proteins or localized inside exosomes, which protects them from endogenous RNase activity (Moreno-Moya ). Several studies have evaluated the potential of circulating miRNAs as diagnostic marker for endometriosis (Jia ; Suryawanshi ; Wang ; Hsu ; Cho ; Rekker ; Cosar ; Nothnick ; Wang ; Bashti ; Maged ; Pateisky ). However, these studies have shown limited consistency and conflicting results.

The goal of this study was to discover and validate plasma miRNAs that can be used for the non-invasive diagnosis of endometriosis. We performed a genome-wide miRNA expression analysis in a large discovery cohort to identify a set of plasma miRNAs with discriminative potential between patients with and without endometriosis. This set of miRNAs was used to build diagnostic models for endometriosis then tested for validation in an independent patient cohort.

Materials and Methods

The project consisted of three parts: (i) biomarker discovery based on a genome-wide miRNA expression profiling by small RNA sequencing, (ii) development of diagnostic models based on targeted expression profiling using RT-qPCR for selected miRNAs and (iii) validation tests of the diagnostic models in an independent cohort (Fig. 1). All patients had signed a written informed consent prior to recruitment, and the study protocol was approved by the Institutional Ethical and Review Board of University Hospital Leuven (Belgium).

Figure 1

Overview of the project.

MiRNA biomarker discovery

Discovery cohort

In the first part of our study, we analysed the plasma miRNA expression profile of 120 patients from the biobank of the Leuven Endometriosis Center of Excellence. All patients had laparoscopically and histologically proven presence (n = 82) or absence (n = 38) of endometriosis and had not used any hormonal medication in the 3 months prior to laparoscopy. Patients who had a visual diagnosis of endometriosis during laparoscopy that was not confirmed by histology were not selected. In the control group, we only selected patients without endometriosis lesions during laparoscopic inspection. Blood had been collected in EDTA tubes on the day of laparoscopy before induction of anaesthesia. Within 1 h, samples were centrifuged at 1400 ×g for 10 min at 4°C, then aliquoted, labelled and stored at −80°C until analysis. The plasma samples were obtained during the different phases of the menstrual cycle (menstrual n = 25, follicular n = 43, luteal n = 52). The cycle phase identification was based on histological examination of an endometrial biopsy taken during laparoscopy. Women with endometriosis were classified as minimal/mild or moderate/severe disease (ASRM endometriosis stages: I–II n = 41 and III–IV n=41) (1997). Patient characteristics are summarized in Table I.

Table I

Patient characteristics of the discovery and validation cohort.

	Discovery cohort (n = 120)				Validation cohort (n = 90)
	Controls (n = 38)	Endometriosis (n = 82)	All (n = 120)	1 P-value	Controls (n = 30)	Endometriosis (n = 60)	All (n = 90)	1 P-value	2 P-value
Age (years, mean ± SD)	30 ± 4	32 ± 4	31 ± 4	0.007	31 ± 5	31 ± 4	31 ± 4	0.855	0.324
BMI (kg/m², mean ± SD)	23.6 ± 3.6	22.5 ± 5.1	22.8 ± 4.7	0.224	24.0 ± 3.5	21.6 ± 2.7	22.4 ± 3.2	0.001	0.495
Smoking, n (%)	11 (29%)	11 (13%)	22 (18%)	0.097	11 (37%)	13 (22%)	24 (27%)	0.193	0.124
Cycle phase, n (%)
Menstrual	8 (21%)	17 (21%)	25 (21%)	0.597	4 (13%)	7 (12%)	11 (12%)	0.820	0.459
Follicular	13 (34%)	30 (37%)	43 (36%)	0.978	13 (43%)	26 (43%)	39 (43%)	1.000	0.297
Luteal	17 (45%)	35 (43%)	52 (43%)	0.586	13 (43%)	27 (45%)	40 (44%)	0.881	0.660
Infertility, n (%)
All infertility	35 (92%)	78 (95%)	113 (94%)	0.928	30 (100%)	59 (98%)	89 (99%)	0.477	0.120
Primary	21 (55%)	61 (74%)	82 (68%)	0.036	17 (57%)	45 (75%)	62 (69%)	0.126	0.932
Secondary	14 (37%)	17 (21%)	31 (26%)	0.061	13 (43%)	14 (23%)	27 (30%)	0.051	0.504
Pain symptoms, n (%)
Any pain symptoms	30 (79%)	73 (89%)	103 (86%)	0.161	21 (70%)	47 (78%)	68 (76%)	0.386	0.062
Dyspareunia	17 (45%)	16 (20%)	33 (28%)	0.011	4 (13%)	14 (23%)	18 (20%)	0.264	0.130
Dysmenorrhea	25 (66%)	68 (83%)	93 (78%)	0.037	19 (63%)	45 (75%)	64 (71%)	0.129	0.486
Dyschezia	2 (5%)	5 (6%)	7 (6%)	0.407	4 (13%)	1 (2%)	5 (6%)	0.015	0.652
Chronic pelvic pain	5 (13%)	12 (15%)	17 (14%)	0.829	2 (7%)	10 (17%)	12 (13%)	0.681	0.779
ASRM stage, n (%)
Stage I-II	NA	41 (50%)	NA	NA	NA	41 (68%)	NA	NA	0.029
Stage III-IV	NA	41 (50%)	NA	NA	NA	19 (32%)	NA	NA	0.029
Ultrasound-negative endometriosis n (%)	NA	61 (74%)	NA	NA	NA	50 (83%)	NA	NA	0.203

NA = not applicable.

A t-test was used for comparison of continuous variables (age & BMI) and chi-square test for categorical variables

1Controls versus endometriosis.

2Discovery cohort versus validation cohort.

Small RNA sequencing

RNA was extracted from 200 μl of plasma using the miRNeasy Serum/Plasma Kit (Qiagen, Germany) according to the manufacturer’s protocol. Libraries for small RNA sequencing were prepared using the NEBNext Small RNA Library Prep Set for Illumina (New England Biolabs, USA). Briefly, 6 μl of total RNA was used as the input for RNA adapter ligation (using 3′ and 5′ RNA adapters) prior to reverse transcription and PCR amplification with bar-coded primers. The PCR products were pooled based on equal volumes then used for size selection on a Pippin Prep (Sage Science, USA) to recover the fractions containing mature miRNAs. The resulting small RNA libraries were concentrated by ethanol precipitation and quantified using a Qubit 2.0 Fluorometer (Thermo Fisher Scientific, USA) prior to sequencing on a NextSeq 500 sequencer (Illumina, USA) with read lengths of 75 base pairs and 15 million single-end reads per sample, on average.

Sequencing reads were processed using the data processing pipeline Cobra (Biogazelle, Belgium). First, reads were filtered based on stringent read quality control. After adapter trimming, reads were collapsed and mapped to the reference genome (Homo sapiens, GRC38) using Bowtie (Langmead ). Mapped reads were subsequently annotated to different contaminants (tRNA fragments, rRNA, sn(o) RNA, piRNA, etc.) and mature miRNAs using genome annotation data from Ensembl 76, UCSC and miRBase v21. We used the sum of all isomiR reads from the canonical mature miRNA locus for miRNA expression analysis. Prior to normalization, data were filtered using a cut-off of 4fourreads (i.e. only those miRNAs with four or more reads were considered expressed). miRNA expression data were normalized based on the total read count per sample, i.e. the read count for each miRNA was divided by the total read count in that sample and multiplied by the median total read count across all samples. After this normalization step, the data were log2 transformed. To assess possible hemolysis, we calculated the miR-Ratio, i.e. the fold change between hsa-miR-451a and hsa-miR-23a-3p (Blondal , Shah ).

Signature selection

Differentially expressed miRNAs were identified using a Mann–Whitney U test on the normalized miRNA expression data, followed by multiple testing correction according to the Benjamini–Hochberg method. Only miRNAs that were expressed in 80% of samples of either one of both analysed patient subgroups were included.

Multivariate logistic regression with stepwise feature selection was used to select miRNAs with discriminative potential between disease and controls. We performed three comparisons (all disease versus all controls, all stage I–II versus all controls and all stage III–IV versus all controls) in all cycle phases combined and in the separate cycle phases (menstrual, follicular and luteal), resulting in a total of 12 analyses. Only miRNAs that were expressed in 80% of samples of either one of both analysed patient subgroups were included. We selected the top 10 discriminative miRNAs from the 12 analyses for targeted miRNA profiling with RT-qPCR. Due to significant overlap between the different analyses, this resulted in a panel of 42 miRNAs.

Development of the diagnostic model

In the second part of our study, we performed a targeted plasma miRNA expression profiling with RT-qPCR for the set of 42 miRNAs selected in the first part of the project. We used the same discovery cohort (n = 120) from the biobank of the Leuven Endometriosis Center of Excellence. RNA had already been extracted prior to small RNA sequencing profiling.

To identify candidate reference miRNAs, stably expressed in the sample cohort, we selected miRNAs expressed in all samples based on the available small RNA sequencing data. Subsequently, we calculated the coefficient of variation (CV) for each miRNA using the normalized miRNA expression data. The top 10 miRNAs with lowest CV, not belonging to the same miRNA cluster or located at the same genomic locus, were selected for geNorm analysis (Hellemans et al., 2007). The geNorm analysis was performed using Biogazelle’s qbase+ software (www.qbaseplus.com) using log2-transformed miRNA count data. This geNorm study was characterized by medium reference target stability (average geNorm M ≤ 1.0). The optimal number of reference targets in this experimental situation is five (geNorm V < 0.15 when comparing a normalization factor based on the four or five most stable targets). As such, the optimal normalization factor can be calculated as the geometric mean of reference targets hsa-miR-361-3p, hsa-miR-423-3p, hsa-miR-28-3p, hsa-miR-191-5p, hsa-miR-425-5p. These miRNAs were considered as candidate reference miRNAs.

For cDNA preparation, 1.5 μl of total RNA was reverse transcribed using the miScript II RT Kit (Qiagen, Germany). cDNA was pre-amplified in a 12-cycle PCR reaction using the miScript PreAMP PCR Kit (Qiagen, Germany). Pre-amplified cDNA was diluted and used as input for a 40-cycle qPCR reaction, quantifying 42 candidate miRNAs of interest, the five candidate reference miRNAs (hsa-miR-361-3p, hsa-miR-423-3p, hsa-miR-28-3p, hsa-miR-191-5p and hsa-miR-425-5p) and two controls (positive PCR control and miRNA reverse transcription control) using miScript Primer Assays (Qiagen, Germany) with the miScript SYBR Green PCR Kit (Qiagen, Germany). All reactions were performed in 384-well plates on the CFX instrument (Bio-Rad, USA) using the gene maximization strategy. Cq-values were determined based on the single threshold algorithm and filtered using a detection cut-off of 29 cycles.

To confirm candidate reference miRNAs, a geNorm pilot experiment was performed using the RT-qPCR data of five candidate reference miRNAs on all samples. Data analysis was performed using qbase+ software (www.qbaseplus.com) (Biogazelle, Belgium). This geNorm study was characterized by high reference target stability (average geNorm M ≤ 0.5). The optimal number of reference targets in this experimental situation is two (geNorm V < 0.15 when comparing a normalization factor based on the two or three most stable targets). As such, the optimal normalization factor can be calculated as the geometric mean of reference targets hsa-miR-28-3p and hsa-miR-423-3p.

Model building

Diagnostic models were built using multivariate logistic regression with stepwise feature selection and a 5% significance level for variables to enter or stay in the model. We built three diagnostic models: the first model to discriminate between controls and all patients with endometriosis and the second and third models to discriminate between controls and patients with minimal–mild (stage I–II) or moderate–severe (stage III–IV) endometriosis, respectively. To determine the optimal cut-off for the diagnostic models, the Youden index was used.

Validation in independent cohort

In the final part of our study, we validated the diagnostic models in an independent validation cohort of 90 patients from the biobank of the Leuven Endometriosis Center of Excellence. All patients had laparoscopically and histologically proven presence (n = 60) or absence (n = 30) of endometriosis and had not used any hormonal medication in the 3 months before laparoscopy. Plasma was collected on the day of laparoscopy before induction of anaesthesia. The plasma samples were obtained during the different phases of the menstrual cycle (menstrual n = 11, follicular n = 39, luteal n = 40). The cycle phase was based on histological examination of an endometrial biopsy taken during laparoscopy. Women with endometriosis were classified as minimal/mild or moderate/severe disease (ASRM endometriosis stages: I–II n = 41, III–IV n = 19) (1997). Patient demographics are summarized in Table I. RNA extraction, cDNA synthesis and qPCR were performed with the same technology and protocols as described for the discovery cohort. The model algorithm trained in the discovery cohort was tested in this independent test set.

Results

Description of the discovery and validation cohort

The clinical characteristics of the patients from the discovery and validation cohort are presented in Table I. There were no significant differences between the discovery and validation cohort, except for the significantly larger proportion of moderate/severe endometriosis in the discovery cohort (P = 0.029).

In the discovery cohort, the endometriosis group was older (32 ± 4) than the control group (30 ± 4). Although statistically significant (P = 0.007), it is unlikely that this small age difference is clinically relevant. In the validation cohort, the endometriosis group had a significantly lower BMI (21.6 ± 2.7) than the control group (24.0 ± 3.5) (P = 0.001). In both patient cohorts, there were more patients with primary infertility in the endometriosis group than in the control group and thus an opposite distribution in secondary infertility. Only in the discovery cohort this reached significance (P = 0.036). In both patient cohorts, the vast majority of patients had some pain symptoms suggestive of endometriosis. In the discovery cohort, the endometriosis group had significantly less dyspareunia (P = 0.011) and more dysmenorrhea (P = 0.037). In the validation cohort, the control group had significantly more dyschezia (P = 0.015).

Quality control of plasma samples

Visual inspection of the plasma samples to assess possible hemolysis identified 16 out of 120 samples (13%) in the discovery cohort that showed a discolouration: 10 samples (8%) were light orange and 6 samples (5%) were dark orange or red. In the validation cohort, 20 out of 90 samples (22%) were discoloured: 16 samples light orange (18%) and 4 samples (4%) dark orange or red (4%). Using the miR-Ratio thresholds of 5 and 7 for assessment of hemolysis in our discovery cohort, we found 39 samples had a low risk of hemolysis (miR-Ratio <5), 31 samples had an intermediate risk (5< miR-Ratio <7) and 50 samples had a high risk of hemolysis (miR-ratio >7) (Shah ). There was no significant difference in the fraction of samples at intermediate and/or high risk of hemolysis in controls versus endometriosis patients (P = 0.44, Chi square test) (Supplementary Fig. S1). Since these thresholds set for the detection of sample hemolysis are based on RT-qPCR data, it is unclear if extrapolation to next generation sequencing (NGS) data is possible. Therefore, we analysed expression of two hemolysis sensitive miRNAs: miR-16-5p and miR-15b-5p. We did not see an increased expression of hsa-miR-16-5p and hsa-miR-15b-3p in the samples with a high risk for hemolysis based on the miR-Ratio (P = 0.56 and P = 0.49, Kruskal–Wallis test; see Supplementary Fig. S1).

We did not see any effect of biobanking time on RNA concentration (Supplementary Fig. S2), and there were no significant differences in RNA concentrations between cases and controls in either the discovery or validation cohorts (P = 0.72 and P = 0.26, Mann–Whitney U test; Supplementary Fig. S3).

MicroRNA biomarker discovery

Two samples of control patients from the discovery cohort were excluded from the analysis because of poor small RNA library quality. Univariate analysis was performed for six comparisons using a Mann–Whitney U test: all stages versus all controls, all stage I–II versus all controls, all stage III–IV versus all controls, all disease versus all controls in the luteal phase, all disease versus all controls in the follicular phase and all disease versus all controls in the menstrual phase. This resulted in 39 (28 up, 11 down), 21 (15 up, 6 down), 49 (34 up, 15 down), 24 (19 up, 5 down), 2 (2 up) and 31 (18 up, 13 down) differentially expressed (P < 0.05) miRNAs, respectively (Supplementary Table SI). However, after correction for multiple testing according to the Benjamini–Hochberg method, no miRNAs remained differentially expressed (adjusted P < 0.05).

To identify miRNAs with discriminative potential between patients with and without endometriosis, irrespective of their differential expression in the univariate analysis, we performed multivariate logistic regression. Firstly, we compared cases and controls in all cycle phases combined and in the separate cycle phases (menstrual, follicular and luteal). Secondly, we performed a sub-analysis for minimal–mild and moderate–severe endometriosis, resulting in a total of 12 analyses. miRNAs were ranked according to their discriminative power, and we selected the top 10 miRNAs from the 12 analyses for confirmation of miRNA expression with RT-qPCR. Because of important overlap between the top-ranked miRNAs in the different analyses, this resulted in a panel of 42 miRNAs, as listed in Table II. Of the 42 miRNAs, 5 were differentially expressed in the univariate analysis (i.e. P < 0.05 in minimally 1 of the 6 comparisons using Mann–Whitney U test, before multiple testing correction) (Table II).

Table II

The 42 miRNAs selected by multivariate logistic regression on the RNA-seq data.

The panel of 42 miRNAs selected for validation
hsa-let-7a-5p	hsa-miR-103a-3p	hsa-miR-15b-5p	hsa-miR-199b-5p	hsa-miR-23a-3p	hsa-miR-29b-3p
hsa-let-7b-5p	hsa-miR-106a-5p	hsa-miR-16-5p	hsa-miR-19a-3p	hsa-miR-24-3p	hsa-miR-30a-3p
hsa-let-7c-5p	hsa-miR-107	hsa-miR-17-3p*	hsa-miR-19b-3p	hsa-miR-25-3p	hsa-miR-30a-5p
hsa-let-7d-5p	hsa-miR-10a-5p	hsa-miR-17-5p	hsa-miR-20a-5p	hsa-miR-26a-5p	hsa-miR-33a-5p
hsa-let-7e-5p	hsa-miR-125b-5p	hsa-miR-182-5p*	hsa-miR-21-5p*	hsa-miR-26b-5p	hsa-miR-92a-3p
hsa-let-7f-5p	hsa-miR-148a-3p	hsa-miR-18a-5p	hsa-miR-210-3p*	hsa-miR-28-5p	hsa-miR-95-3p
hsa-miR-101-3p	hsa-miR-15a-5p	hsa-miR-199a-5p	hsa-miR-22-3p	hsa-miR-29a-3p	hsa-miR-98-5p

*Significant in the univariate analysis (before adjustment for multiple testing).

Development of diagnostic models

In the second part of the project, we performed a targeted miRNA profiling for the panel of 42 miRNAs using RT-qPCR in the same discovery cohort. However, for hsa-miR-17-5p and hsa-miR-106a-5p, no separate qPCR primer is available in the miScript assays of Qiagen because only one base is different between the two miRNA sequences. Therefore, miRNA expression was confirmed for only 41 miRNAs. The miScript Assay for hsa-miR-17-5p was used and quantified both hsa-miR-17-5p and hsa-miR-106a-5p as a single miRNA. As a result, we had 41 features (miRNAs) in the downstream data analysis.

We used multivariate logistic regression with stepwise feature selection to build three diagnostic models for endometriosis: model 1 to discriminate between controls and patients with all stages of endometriosis and models 2 and 3 to discriminate between controls and patients with minimal–mild (stage I–II) or moderate–severe (stage III–IV) endometriosis, respectively (Table III). All models consisted of three miRNAs with hsa-miR-28-5p present in all models and hsa-miR-21-5p present in model 1 and model 3. Model 1 (hsa-let-7d-5p, hsa-miR-21-5p and hsa-miR-28-5p) had an AUC of 73% with a high sensitivity (85%) but a low specificity (49%) at optimal cut-off (Table III). Model 2 (hsa-miR-125b-5p, hsa-miR-28-5p and hsa-miR-29a-3p) performed slightly better with an AUC of 77% and could detect minimal–mild endometriosis with a sensitivity of 80% and specificity of 62% at optimal cut-off (Table III). Model 3 (hsa-miR-21-5p, hsa-miR-28-5p and hsa-miR-30a-5p) had the best AUC (81%) with a high specificity (89%) for stage III–IV endometriosis and a sensitivity of 63% at optimal cut-off (Table III).

Table III

Overview of the logistic regression models.

Model 1: All endometriosis stages	Coefficient	SE	P-Value	Optimal cut-off	Diagnostic performance
Intercept	1.000	0.236	<0.0001	0.5556	AUC: 73% (64–83%) Sensitivity: 85% Specificity: 49%
hsa-let-7d-5p	3.219	1.363	0.0182
hsa-miR-21-5p	−3.758	1.093	0.0006
hsa-miR-28-5p	5.967	2.473	0.0158
Model 2: Stage I-II endometriosis	Coefficient	SE	P-Value	Optimal cut-off	Diagnostic performance
Intercept	0.321	0.271	0.2370	0.4759	AUC: 77% (66–87%) Sensitivity: 80% Specificity: 62%
hsa-miR-125b-5p	3.551	1.221	0.0036
hsa-miR-28-5p	11.356	3.431	0.0009
hsa-miR-29a-3p	−6.906	2.030	0.0007
Model 3: Stage III-IV endometriosis	Coefficient	SE	P-Value	Optimal cut-off	Diagnostic performance
Intercept	0.198	0.271	0.4664	0.6247	AUC: 81% (71–90%) Sensitivity: 63% Specificity: 89%
hsa-miR-21-5p	−4.794	1.357	0.0004
hsa-miR-28-5p	11.343	3.270	0.0005
hsa-miR-30a-5p	6.874	2.610	0.0084

Validation of diagnostic models

In the final step of our project, we validated the diagnostic models in an independent validation cohort (Fig. 2). Models 1 and 3 experienced a drop in diagnostic accuracy to levels of chance performance with AUCs of 50% (Fig. 2). Model 2 had an AUC of 60% in the validation cohort and could, at optimal cut-off, detect minimal–mild endometriosis with an acceptable sensitivity of 78% but with a very low specificity of 37% (Fig. 2).

Figure 2

Validation of the diagnostic models in an independent patient cohort. Black curve = Discovery cohort; Red curve = validation cohort. Data in the contingency tables are from the validation cohort.

Effect of menstrual cycle phase on miRNA levels

To assess a potential effect of the phase of the menstrual cycle on circulating miRNA levels, we compared the levels of the 41 miRNAs in the luteal, follicular and menstrual phase (Supplementary Tables SII and SIII). In the discovery cohort, we found 16 miRNAs with a significant difference between the levels in either the luteal, menstrual or follicular phase: hsa-let-7c-5p (P = 0.027), hsa-miR-101-3p (P = 0.004), hsa-miR-103a-3p (P = 0.009), hsa-miR-10a-5p (P = 0.020), hsa-miR-125b-5p (P = 0.008), hsa-miR-17-3p (P = 0.006), hsa-miR-182-5p (P = 0.011), hsa-miR-18a-5p (P = 0.016), hsa-miR-199b-5p (P = 0.009), hsa-miR-210-3p (P = 0.008), hsa-miR-24-3p (P = 0.013), hsa-miR-25-3p (P = 0.042), hsa-miR-26a-5p (P = 0.024), hsa-miR-26b-5p (P = 0.030), hsa-92a-3p (P = 0.002) and hsa-miR-95-3p (P = 0.008). In the validation cohort only four miRNAs had a possible cycle effect: hsa-miR-let-7d-5p (P = 0.009), miR-15b-5p (P = 0.05), hsa-miR-199a-5p (P = 0.031), hsa-miR-30a-5p (P = 0.022). No miRNA had significant cycle changes in both the discovery and validation cohort. Therefore, it is more likely that the observed changes are due to biological variability rather than a true menstrual cycle effect.

Discussion

In this study, we investigated whether plasma miRNAs can be used for the non-invasive diagnosis of endometriosis. First, we identified a set of 42 plasma miRNAs with discriminative potential between patients with and without endometriosis based on a genome-wide miRNA expression profiling. Secondly, we validated miRNA expression with RT-qPCR and developed three diagnostic models in the discovery cohort of 120 patients. Finally, the models were validated in an independent cohort of 90 patients. Two out of three models did not have a diagnostic potential above chance performance in the validation cohort. The diagnostic model for minimal–mild endometriosis had an acceptable sensitivity of 78% but with a very low specificity of 37% and an AUC of only 60%.

The strength of our study lies in its rigorous methodological design. To minimize potential inaccuracies due to sample processing, our biobank has adopted strict standard operating procedures (Fassbender ). All patients were operated in the Leuven Endometriosis Center of Excellence, a tertiary care unit for endometriosis patients in Belgium, receiving pre- and post-operative care according to a strict clinical protocol (Meuleman ). Staging of endometriosis was done by two experienced endometriosis surgeons, again minimizing potential variability. The use of standard protocols for both clinical management and sample collection allowed us to minimize the clinical and technical variability in the two patient cohorts.

Secondly, we used a workflow with highly validated standard operating procedures for miRNA expression analysis and data processing to maximize the accuracy of our results. In the discovery phase, we used next generation sequencing for miRNA expression analysis. This state-of-the art technology is ideal for unbiased genome-wide miRNA expression profiling. For the development and validation of the diagnostic models, we used RT-qPCR that has high specificity and linear dynamic range of quantification, making it the gold standard for expression profiling of even low abundant miRNAs. Normalization of miRNA expression data is challenging since no standard reference or housekeeping miRNAs have been identified. Therefore, we used qbase+ software for data processing and normalization. This software is built upon an established quantification model including PCR efficiency correction, multiple reference gene normalization and inter-run calibration (Hellemans ).

Finally, we used an approach based on multivariate logistic regression with stepwise feature selection for development of the diagnostic models. The multivariate statistical approach allowed us to model the relationship between the diagnostic categories and all of the miRNAs simultaneously, while taking into account the correlation that may exist between the miRNAs. We specifically chose to build diagnostic models without taking the required number of miRNAs into account since the combination of multiple biomarkers might be necessary to capture the complex dynamics of the endometriosis disease process. We did so because in practice, a substantial part of the workload for miRNA expression analysis relates to RNA extraction, which is independent of the number of miRNAs to be analysed. The RT-qPCR can be fully automated, and the number of analysed miRNAs has relatively little impact on cost and time.

While our study focused on minimizing potential variability and bias, several limitations must be addressed. Firstly, this is not a prospective study but a retrospective, biobank-based cohort study. All patients were selected from the biobank of the Leuven University Fertility Center, and, as stated in the methods’ section, we excluded samples from patients receiving hormonal suppression therapy at the time of surgery. We included samples from patients with laparoscopically proven presence or absence of endometriosis, with three possible indications for surgery: surgical treatment of endometriosis diagnosed on imaging or diagnostic laparoscopy, investigation of pain symptoms suggestive of endometriosis or a diagnostic laparoscopy as part of the diagnostic work-up for infertility. In routine practice in our centre, laparoscopies in women with (unexplained) infertility are predominantly performed when they also have pain symptoms. This explains the high pain prevalence in both cases and controls (Table I). Furthermore, the proportion of patients with only pain symptoms and no infertility was rather low (Table I), because patients with pain suggestive of endometriosis and without intention to conceive in the near future, typically receive hormonal therapy as first line treatment. Since we excluded samples from patients using hormonal medication at the time of surgery, we predominantly selected samples from infertile patients with pain symptoms suggestive of endometriosis. Finally, excluding patients under hormonal therapy resulted in selection of a high number of patients with ultrasound negative endometriosis because, in our hospital, it is common practice that patients with evidence of endometriosis on imaging (ultrasound or MRI) are given hormonal suppression (oral contraceptives or GnRH analogues) in the pre-operative period. Overall, our study was done in a highly selected population of predominantly infertile women with pain symptoms suggestive of endometriosis and a normal pelvic ultrasound. More research is needed to confirm our data in patient populations with pain and without infertility, but this research effort will be confounded by the fact that many patients with endometriosis-associated pain receive hormonal suppression therapy.

In order to address the challenge in the selection of a control group in endometriosis biomarker studies, both our study and control samples were obtained from a similar and well-characterized patient population (Table I), i.e. women with subfertility and pain symptoms suggestive of endometriosis. Since a diagnostic laparoscopy is often considered in these patients, they represent a relevant patient group for which a non-invasive diagnostic test for endometriosis would be useful (Vodolazkaia , Fassbender ). To show true promise as an aid to diagnosis, a biomarker for endometriosis needs to distinguish women with endometriosis from unaffected women with a similar presentation (May ), In our study, both cases and controls indeed had a similar clinical presentation and indication for surgery (nearly all had infertility and a majority had pain symptoms suggestive of endometriosis) and only differed with respect to the presence or absence of macroscopic endometriosis confirmed by histology.

Another important limitation of our study is the relatively small number of patients in the two study cohorts. Although 210 patients represent the largest cohort in which the potential of circulating miRNAs as biomarkers for endometriosis has been studied so far, this number still represents a relatively small population for model building. However, the large variability and conflicting results in the existing literature on differential expression of circulating miRNAs between patients with and without endometriosis did not allow us to make a reliable estimation of an expected effect size, and therefore, a formal power calculation was not possible at the start of our study. We attempted to address the issue of statistical power in the methodological design of the project. We anticipated a certain degree of overfitting in the multivariate analysis in the discovery phase. However, multiple comparison or overfitting of a multivariate analysis does not result in loss of significance; it results in more type 1 error. We accepted this possibility of false positive selection of features, since we considered the discovery phase as exploratory.

In the model building phase, we used multivariate logistic regression with stepwise feature selection and a 5% significance level for variables to enter or stay in the model. This technique selects the largest effects first (i.e. the miRNAs with the best discriminative potential) and then adds smaller effects. We recognize that the relatively small sample size in the model building phase might have failed to detect small effects. However, since a model combining the three largest effects did not withstand independent validation, it seems unlikely that addition of smaller effects (i.e. miRNAs with very limited additional discriminative power) would have changed this observation. Ideally, the validation cohort would have been larger. However, a small validation set does not result in an underestimation of the AUC; it merely results in a wider confidence interval. In our data set, a narrower confidence would not have changed the overall conclusion for any of the three diagnostic models.

Several studies investigating the potential of circulating miRNAs as biomarkers for endometriosis have been published (Jia , Suryawanshi , Wang , Hsu , Cho , Rekker , Cosar , Nothnick , Wang , Bashti , Maged , Pateisky ). These studies have shown differential expression of circulating miRNAs between patients with and without endometriosis. However, these studies have shown limited consistency or conflicting results, and no miRNA-based diagnostic test has been validated in an independent patient cohort so far. The disparities between the different studies and the inability to validate initial findings can be explained by multiple sources of variability in miRNA research related to endometriosis. There are pre-analytical factors such as the source of the miRNAs (serum or plasma), type of blood collection tubes (EDTA versus heparin), hemolysis, sample processing protocol, etc. Secondly, there are technical factors related to the method used for RNA extraction, miRNAs expression analysis (micro-array, RT-qPCR and NGS) and the strategy for normalization of miRNA expression data. Finally, there are biological factors: genetic background of the study cohort, the control population (self-reported healthy versus laparoscopically proven absence of endometriosis), extent of endometriosis (stage I versus stage IV), etc.

In the discovery phase of the project, we applied a genome-wide approach, and when we selected the panel of 42 candidate miRNAs we noted some interesting overlap with previously published studies: 8 out 42 miRNAs in our panel had been reported as potential circulating biomarkers for endometriosis in other studies (Jia , Wang , Cho , Cosar ). From the miRNAs included in the three diagnostic models (hsa-let-7d-5p, hsa-miR-125b-5p, hsa-miR-21-5p, hsa-miR-28-5p, hsa-miR-29a-3p and hsa-miR-30a-5p), two miRNAs have been suggested as circulating biomarkers for endometriosis (hsa-let-7d-5p and hsa-miR-125b-5p) (Cho , Cosar ), and three miRNAs have been implied to play a role in the pathophysiology of endometriosis (hsa-miR-30a-5p, hsa-miR-125b-5p and hsa-miR-21-5p) (Chang , Aoyagi , Haikalis , Liu , Park ). Our study supports a possible link between these miRNAs and endometriosis. However, as shown in our study, the potential of these miRNAs as clinically useful biomarkers is questionable.

In our project, we tried to minimize potential biological and technical variability by selecting the samples from the same biobank and using identical methodology for miRNA expression analysis. Nevertheless, we failed to validate our diagnostic models. Our data suggest an important natural variation in circulating miRNA levels between individuals, with only very modest or minimal changes related to endometriosis. Interestingly, we found that, in the independent validation cohort, the diagnostic miRNA model for stage I–II endometriosis performed slightly better than the model for stage III–IV endometriosis. Intuitively, one would expect a higher difference between stage III–IV endometriosis and controls. The better diagnostic performance from the models for stage I–II endometriosis might reflect an influence of endometriosis phenotype on miRNA expression. The stage III–IV group is more heterogenous since it consists of patients with extensive adhesions, endometriomas, deep endometriotic nodules, etc., whereas patients in the stage I–II group almost exclusively have peritoneal endometriosis and sometimes limited adhesions.

Conclusion

In this study, we investigated the potential of plasma miRNAs as biomarkers for endometriosis in women with infertility and pain. Our study supports a possible biological link between certain miRNAs and endometriosis. The potential of these miRNAs as clinically useful biomarkers is questionable given the important individual natural variation in miRNA levels. Large studies in well-described patient cohorts, with rigorous methodology for miRNA expression analysis, sufficient statistical power and an independent validation step, are necessary to answer the question of whether miRNAs can be used for the non-invasive diagnosis of endometriosis.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.

Click here for additional data file.