Steroid hormones regulate genome-wide epigenetic programming and gene transcription in human endometrial cells with marked aberrancies in endometriosis.

Programmed cellular responses to cycling ovarian-derived steroid hormones are central to normal endometrial function. Abnormalities therein, as in the estrogen-dependent, progesterone-"resistant" disorder, endometriosis, predispose to infertility and poor pregnancy outcomes. The endometrial stromal fibroblast (eSF) is a master regulator of pregnancy success. However, the complex hormone-epigenome-transcriptome interplay in eSF by each individual steroid hormone, estradiol (E2) and/or progesterone (P4), under physiologic and pathophysiologic conditions, is poorly understood and was investigated herein. Genome-wide analysis in normal, early and late stage eutopic eSF revealed: i) In contrast to P4, E2 extensively affected the eSF DNA methylome and transcriptome. Importantly, E2 resulted in a more open versus closed chromatin, confirmed by histone modification analysis. Combined E2 with P4 affected a totally different landscape than E2 or P4 alone. ii) P4 responses were aberrant in early and late stage endometriosis, and mapping differentially methylated CpG sites with progesterone receptor targets from the literature revealed different but not decreased P4-targets, leading to question the P4-"resistant" phenotype in endometriosis. Interestingly, an aberrant E2-response was noted in eSF from endometriosis women; iii) Steroid hormones affected specific genomic contexts and locations, significantly enriching enhancers and intergenic regions and minimally involving proximal promoters and CpG islands, regardless of hormone type and eSF disease state. iv) In eSF from women with endometriosis, aberrant hormone-induced methylation signatures were mainly due to existing DNA methylation marks prior to hormone treatments and involved known endometriosis genes and pathways. v) Distinct DNA methylation and transcriptomic signatures revealed early and late stage endometriosis comprise unique disease subtypes. Taken together, the data herein, for the first time, provide significant insight into the hormone-epigenome-transcriptome interplay of each steroid hormone in normal eSF, and aberrant E2 response, distinct disease subtypes, and pre-existing epigenetic aberrancies in the setting of endometriosis, provide mechanistic insights into how endometriosis affects endometrial function/dysfunction.

Introduction

Endometrium is a dynamic tissue whose cellular components undergo cyclic proliferation and differentiation, preparing for embryo implantation by highly coordinated spatiotemporal actions of ovarian-derived estradiol (E2) and progesterone (P4) [1,2]. These hormones bind cognate receptors [estrogen receptor (ER) and progesterone receptor (PR)], whose activities are tightly regulated by post-translational modifications and interactions with cell- and tissue-specific co-regulators [3-5]. Binding ER and PR leads to their nuclear translocation, complexing with nuclear response elements, remodeling chromatin by co-modulator recruitment [3], and alteration of the transcriptional machinery. In endometrium, dynamic circulating E2 and P4 levels drive the normal functionality of the tissue. Moreover, environmental and inflammatory signals can alter steroid hormone-driven endometrial gene transcription and cellular function resulting in tissue dyshomeostasis [6], including endometrial hyperplasia and cancer, endometrial-based infertility, endometriosis, and poor pregnancy outcomes [7,8]. While changes in chromatin accessibility, PR targets and changes in histones and gene expression in eSF decidualization by E2, cAMP and MPA have been shown [9-11], how E2 and P4 individually interact with the endometrial epigenome normally or in inflammatory disorders that compromise endometrial function, e.g., as in the disorder endometriosis, are incompletely understood. We hypothesized that these steroid hormones induce unique genome-wide signatures in normal human endometrial stromal fibroblasts with aberrant signatures in endometrial cells from women with endometriosis, and their effect on the epigenome is directed by specific genomic sequences and locations.

Endometriosis is a common, chronic disorder wherein endometrial tissue, shed into the pelvis at menses, elicits an inflammatory response, neuroangiogenesis and fibrosis, resulting in infertility and chronic pelvic pain [12]. Hallmarks of the disorder are its dependence on E2 for growth, disrupted P4 signaling caused by chronic inflammation in endometriosis lesions and in the eutopic endometrium (uterine lining) [13], and epigenetic chromatin changes that determine endometrial cellular responses to mitogenic and differentiative signals [6,12,14]. Normally, the eutopic endometrial DNA methylome varies according to the hormonal milieu, with greatest differences in the E2-dominant (proliferative) versus P4-dominant (secretory) phase of the cycle and associated with gene expression changes [10,13-16]. Chronic inflammation affects the chromatin landscape in endometrium of women and animal models of endometriosis [6]. These observations on bulk tissue raise fundamental questions about steroid hormone-epigenome interactions in cellular components of the endometrium normally and in women with disease, how steroid hormones affect the epigenome, how the epigenome affects steroid hormone response and action, and if there are epigenetic differences in endometrium of women with endometriosis and different stages of disease, and if so, what role they play in these processes.

Herein, we studied responses of endometrial stromal fibroblasts (eSF) isolated from normal women (controls) and those with endometriosis. eSF comprise a major endometrial cell type whose programmed responses to E2 and P4 are essential for pregnancy success and whose responses are compromised in inflammatory disorders [17], including endometriosis [18,19]. Given the centrality of a normal eSF P4 response for pregnancy and that women with endometriosis have infertility and poor pregnancy outcomes reported by some to be due, in part, to altered eSF basal gene expression and abnormal response to P4 [18-21], understanding steroid hormone-signaling and regulation of transcription in this cell type is paramount. Moreover, controls for this study were women with no known gynecologic or systemic disorders, thereby enabling establishment of a normative platform for steroid hormone effects on the epigenome and gene transcription in this cell type in endometrium, the tissue that is the anatomic pre-requisite for continuation of the species.

Results

Distinct Effects of Ovarian Steroid Hormones in Normal Human Endometrium

We assessed the genome-wide effect of individual steroid hormones, E2 and P4, and their combination (E2+ P4) on endometrial stromal fibroblasts (eSF) after 14 days of exposure mimicking the timeframe in the menstrual cycle for maximal hormone responsiveness. Since the effects of E2 and P4 individually and together on the hormone-epigenome interplay in normal endometrial cells was unknown and was a main aim of this study, we applied stringent criteria and utilized only endometrial samples from extensively screened volunteers without any gynecologic disorders and no uterine pathology (NUP), with confirmed in vitro eSF progesterone responsiveness (see methods and ).

Interrogation of 485,577 methylation targets across the genome revealed that E2 and P4 and their combination affected the DNA methylome to different extents and with distinct patterns in eSFnormal. (Note: throughout the Results section steroid responses of eSF DNA methylomes are compared to untreated (vehicle) cells for each group.) E2 induced the most extensive changes in the eSFnormal DNA methylome (2047 CpG sites), followed by combined E2+P4 (569 CpG sites), and P4 alone had the least effect (505 CpG sites) (). Importantly, combined E2+P4 resulted in dramatically reduced numbers of differentially methylated loci compared to E2 alone (). While individual hormone treatments elicited hormone-specific DNA methylome changes, the simultaneous presence of both hormones altered their individual effects. The pattern of loss and gain of methylation is also distinct for each hormone. E2, and E2+P4 induced more loss than gain of methylation ( yellow: gain of methylation, blue: loss of methylation vs vehicle), while P4 induced similar numbers of loss and gain of methylation. Concordant with differential patterns and extents of methylation changes, we found minimal overlap in the differentially methylated CpG sites affected by each hormone and the majority was unique (, ). In particular, loci differentially methylated in response to E2+P4 were mostly unique compared to those in response to P4 or E2 alone and were not a combination of the response to each hormone individually ().

Pathways and biofunctions () as well as functional enrichment clustering ( were also unique to each hormone with E2 pathways enriching for gap junctions, melanogenesis, and glutamatergic and dopaminergic synapses pathways, and zinc and ion binding, cell membrane, glycoprotein and signal peptide functional clusters, with fewer and unique statistically significant pathways and functional clusters affected by P4 and E2+P4. Together these data indicate that each hormone affects different regions and E2+P4 targets are not a combination of E2 and P4. Differentially methylated loci in all hormonally treated eSFnormal involved several pathways, many important in normal endometrial function and dysfunction. Important pathways affected by each hormone and the differentially methylated genes in each pathway are shown in The data were further mined for differences in DNA methylation patterns, profiles, and genomic locations, regulatory elements, transcribed genes and biofunctions induced by each hormonal treatment in cells from normal and endometriosis women (see below).

Genomic locations, regulatory elements, CpG islands and neighborhood context

Interestingly, while the patterns, profiles, differentially methylated CpG sites, pathways and biofunctions were unique to each hormone, the genome-wide distribution of their affected CpG sites shows specific enrichments and depletions. All hormones (E2, P4 and E2+P4) were statistically significantly enriched in intergenic regions ( and in enhancers (, albeit with different extents and in gain vs loss of methylation and some variations based on hormones. These may reflect hormone binding sites in these regions, as had been reported in breast cancer cell lines [22,23]. There was a marked depletion of differential methylation in close proximity to transcription start sites (TSS) up to -200 nt upstream (TSS200) for all hormonal treatments in both of gain or loss of methylation (; for gain/loss for each hormone; ). But, CpG sites with gain of methylation in all hormonal treatments exhibited low enrichment at 5’UTRs, 1st exons and gene bodies, while in loss of methylation TSS1500 and 1st exons were less involved but gene bodies and 3’UTRs were more enriched. Greatest differences in gain versus loss of methylation () involved gene bodies, and 3’UTRs, and less at TSS1500, TSS200, 5’UTRs, and 1st exons.

CpG islands (CGI), CGI shores and shelves. For all hormones, there was low involvement of CGIs and CGI shelves and shores, compared to interrogated loci on the HM450 platform ( and ). Most DNA methylation changes in eSFnormal did not involve CGIs. Indeed, while 31% of all interrogated loci were at CGIs and 33% at CGI shores and shelves (total CGI-related 64%), only 5–7% of the differentially methylated loci for any hormonal treatment were located at CGIs (33–40% overall in CGIs, shores and shelves, compared to 64% arrayed on the platform). The majority of differentially methylated CpG sites involved non-CGIs (59–66%) versus 36% non-CGI CpG sites on the platform. However, there were more differentially methylated CpG sites at CGI shores and shelves than in CGIs, in both gain and loss of methylation. Notably, more CGIs, less CGI shelves (north and south shelves), more CGI shores (north and south shores), and less non-CGI CpG sites were involved in gain versus loss of methylation.

Changes in gene expression in response to hormones

Steroid hormones affect their target genes through various mechanisms and as such, changes in DNA methylation may not fully reflect their effect on changes in gene expression, particularly in the case of those loci whose transcriptional regulation does not involve chromatin modifiers. To elucidate a more complete effect of hormones, transcriptomic profiles were determined in the same steroid hormone-treated eSF used for DNA methylation analysis and was compared to its corresponding gene expression profiles in untreated eSF. E2 induced more up- than down-regulated genes ( top up- and down-regulated loci; full gene list), and P4 elicited similar numbers of up- and down-regulated genes, (). However, more genes were differentially expressed when E2 and P4 were combined, with more genes up- than down-regulated ().

E2+P4 induced the largest and P4 the smallest changes in gene expression. Almost all the genes differentially induced by P4 were shared with E2+P4 and some were shared with the E2 treatment. Notably, half of the E2-induced and the majority of E2+P4 induced differentially expressed genes were unique. However, in commonly up-regulated genes between E2 and E2+P4, the variable fold changes indicate inhibitory or stimulatory effects when E2 is combined with P4. For example, PGR is upregulated by both E2 and E2+P4 (FC = 4.5 vs 1.7, respectively), indicating that the addition of P4 limited up-regulation of PGR compared to E2 alone. Among other up-regulated gene in common in E2, P4, E2+P4, are IGF1 and SPARCL1 with known roles in endometrial biology. But, the FCs were different (IGF1: E2 = 21.5, P4 = 4.6, E2+P4 = 13.6; SPARCL1: E2 = 4.6, P4 = 13.2; E2+P4 = 53.1) indicating potentially different mechanisms for up- or down-regulation and for different genes, potentially affected by genomic location and other regulatory factors, modifiers and gene/region-specific mechanisms involved in hormone-induced gene expression regulation.

Pathways and biofunctions

E2 increased tissue and cellular development, growth and maintenance, and downregulated cell-to-cell signaling, immune cell trafficking, inflammatory response, apoptosis and cellular migration (). P4 elicited down-regulation of cellular regeneration and proliferation and cell-cell signaling and adhesion. E2+P4 upregulated cell death and molecular transport and downregulated cell growth and proliferation, carbohydrate metabolism and molecular transport. The genes commonly upregulated by E2, P4, E2+P4 involved catalytic activity, receptor and signal transduction, binding, transporter and structural molecule activity. The main biofunctions of differentially expressed genes that were shared with differentially methylated loci involved cell membrane and signaling in response to E2.

Association of gene expression with DNA methylation

Hormonally-induced differentially methylated CpG sites were assessed for association with differential gene expression for each corresponding locus, noting that not all transcribed loci are included in both platforms and many intergenic regions in DNA methylation platform were not represented on the gene expression array used in this study. Only loci with a strong positive or negative association (by Spearman rho, and corrected p<0.05, see Methods) were considered. There was a large number of functional gene clusters with strong association of DNA methylation and gene expression for E2 in eSFnormal (), which was not observed for P4 or E2+P4 treatments.

Effects of Ovarian Steroid Hormones in Endometrium of Women with Endometriosis

We next aimed to determine the effect of hormones on the endometrium of endometriosis patient, known to have abnormal P4 response. We applied strict criteria using eSF from patients with only endometriosis and no other uterine, pelvic or gynecologic disorders and those that show P4-resistance confirmed by microscopy and IGFBP1 assay (). Furthermore, to understand the effect of disease stage on the hormone-epigenome interplay, we used early (stage I) and late stage (stage IV) disease. Similar to normal, E2 induced the most and P4 the least DNA methylation changes in eSFstage-I and eSFstage-IV (, respectively). But, in eSFendo, the extent of E2-induced changes was less than in eSFnormal specifically in eSFstage-IV (418 CpG sites) exhibiting significantly less methylation alterations compared to eSFnormal and eSFstage-I (2047 and 1633 CpG sites, respectively) (. Opposite to that of eSFnormal the majority of changes in both stages were gain of methylation ( heatmap). The extent of E2-induced differentially methylated loci differed considerably between the two stages (): stage I showed extensive changes induced by E2, much reduced in stage IV. The considerable difference in the extent of E2- induced methylation in eSFstage-IV and in the gain/loss pattern (), indicate an aberrant response to E2 in both stages of disease and more extensively in stage IV, not previously reported. Progesterone, similar to eSFnormal induced the least DNA methylome alterations in both eSFstage-I and eSFstage-IV despite the difference in the robust decidualization response to P4 in eSFnormal and the refractory decidualization response to P4 in eSFendo (. While P4 induced similar numbers of loss and gain of methylation in eSFnormal in eSFstage-I there was more loss than gain of methylation and in eSFstage-IV more gain than loss of methylation. Interestingly, despite complete lack of decidualization, eSFstage-I exhibited more E2+P4-induced differentially methylated loci versus eSFnormal, and in both disease stages there was more loss than gain of methylation (). eSFstage-IV showed the fewest methylome changes in all three hormone treatments, suggesting extensive aberrancies in hormone-methylome interactions in late stage disease. Particularly important is the novel observation of an aberrant response to E2 and not just to P4 and E2+P4, as previously believed [24].

Hormone induced differentially methylated CpG sites in stage I eSF (Endo I) and stage IV eSF (Endo IV).

2A. Differentially methylated CpG sites induced by E2, P4 and E2+P4 versus vehicle. Heatmaps reflect the differential methylation of each sample in each hormone treatment versus its corresponding non-treated vehicle control (Δβ: Hormone treated minus vehicle control) (see Fig 1A legend for details). 2B. Number of differentially methylated CpG sites and in gain/loss of methylation for each hormone treatment in Endo I and Endo IV. C. Unique and common differentially methylated CpG sites for each hormone in Endo I (left) and Endo IV (right) indicating mostly unique loci for each hormone. 2D. Unique and common differentially methylated CpG sites across normal (NUP), Endo I and Endo IV, for each hormone: E2: left, P4: middle and E2+P4: right. 2E. Enrichment of intergenic regions; Endo I, top and Endo IV bottom charts (see Fig 1D legend for details). 2F. Enrichment of enhancers for each hormone and based on loss or gain of methylation in Endo I (left panel) and Endo IV (right panel) (see Fig 1E legend for details). 2G. Genomic distribution of all differentially methylated CpG sites in each group (Endo I, left, Endo IV right panel) and by gain or loss of methylation (see Fig 1F legend for details). 2H. Distribution of differentially methylated CpGs by CpG islands (CGI), CGI north/south shores and shelves for Endo I (left) and Endo IV (right) based on gain or loss of methylation for all hormones. For details see Fig 1G legend. N Shelf: North Shelf; S Shelf: South Shelf; N Shore: North Shore; S Shore: South Shore. Endo I: stage I; Endo IV: stage IV.

Fig 1

Hormone induced differentially methylated CpG sites in normal eSF (NUP).

1A. Differentially methylated CpG sites induced by E2, P4 and E2+P4 versus vehicle. Each heatmap reflects differential methylation of each sample in each hormone treatment versus its corresponding non-treated vehicle control (Δβ: Hormone treated minus vehicle control). Yellow heatmaps above the X-axis reflect gain of methylation vs. vehicle; blue heatmaps below the X-axis reflect Δβ loss of methylation. In each heatmap, rows show Δβ of differentially methylated loci, columns indicate samples. Y-axis shows the number of differentially methylated loci for either gain or loss of methylation for each hormone treatment. 1B. Number of differentially methylated CpG sites and in gain/loss of methylation for each hormone treatment. 1C. Venn diagram of unique and common differentially methylated CpG sites for each hormone shows little overlap between differentially methylated loci in each hormone treatment 1D. Enrichment of intergenic regions in % in each hormone treatment for all differentially methylated loci (All Loci), those with gain or loss of methylation (Gain, Loss) and by individual hormones (E2, P4, E2+P4). Enrichment is assessed by Z-test and p<0.05 are shown in parentheses. Black bar represents percentage of intergenic loci of total interrogated CpG sites. 1E. Statistically significant involvement of enhancers by hormones and gain or loss of methylation. Enrichment is assessed by Z-test and p<0.05 are shown in parentheses. Black bar represents percentage of enhancers of total interrogated CpG sites. 1F. Genomic distribution of all differentially methylated CpG sites in each hormone and by gain or loss of methylation, assessed at TSS1500, TSS200, 5’UTR, 1st exon, gene body, 3’UTR, and intergenic regions. Black line represents the percentage of interrogated CpG site at each location, green line (top panel) shows all differentially methylated loci in NUP for all hormones, yellow line (middle panel) shows all loci with gain of methylation in NUP, and blue line (bottom panel) shows all loci with loss of methylation for all hormones. Enrichment is assessed by Z-test and p<0.05 are shown in parentheses for each genomic location. 1G. Distribution of differentially methylated CpGs by CpG islands (CGI), CGI north/south shores and shelves for all loci with gain of methylation in all hormone treatments (orange lines) or with loss of methylation (blue line) in comparison to the distribution of the interrogated CpG sites in each of these genomic locations (black line). N Shelf: North Shelf; S Shelf: South Shelf; N Shore: North Shore; S Shore: South Shore. NUP: normal (no uterine pathology).

Differentially methylated loci were unique in response to different hormones in each disease stage and between the two stages. (). As in normal, E2+P4 induced methylation were mostly unique and not a combination of the response to E2 or P4 individually (), reaffirming that E2 and P4 interact differently with the epigenome when combined than when individually administered (). Moreover, the majority of loci differentially methylated in response to each specific hormone was also unique in eSFnormal vs eSFstage-I vs eSFstage-IV (). These data suggest that the hormone-DNA methylome dynamics differ under normal and disease conditions, and furthermore that the stage of disease affect the hormone-methylome response.

Despite distinct profile differences with normal, hormone-induced differentially methylated CpGs for E2, P4 and E2+P4in both stages of eSFendo were also statistically significantly enriched in intergenic regions () although the extent of this enrichment differed among hormone treatments, and between disease stages (, for gain/loss for each hormone; ). Similar to normal and in both gain and loss of methylation, there was marked enrichment of enhancers, although the extent differed with specific hormone treatments and disease stage (). In both E2 and P4 treatment, eSFnormal involved more enhancers than eSFstage-I and eSFstage-IV, but E2+P4 treatment induced involvement of more enhancers in eSFstage-IV, particularly in loss of methylation with nearly 50% of CpGs associated with enhancers ().

In eSFstage-I, the genomic distribution of loci with gain or loss of methylation differed from eSFnormal at 1st exons, gene bodies, and 3’UTRs. In eSFstage-IV the genomic distribution was mostly similar in gain and loss of methylation. These differed at 5’UTRs, gene bodies and 3’UTRs compared to eSFnormal and at 5’UTRs and 1st exons compared to eSFstage-I (). Overall, these data demonstrate that: hormone treatments regardless of disease and its stage affected CpG sites more at the 3’UTR and intergenic regions and much less at proximal promoters/TSS; genomic locations of CpG sites differentially methylated in response to hormones differed based on loss/gain of methylation; while decreased proximal promoter (TSS200) involvement and increased intergenic region involvement were common in eSFnormal, eSFstage-I and eSFstage-IV. Low involvement of promoters/TSS and increased involvement of 3’UTR and intergenic regions were remarkable, considering vast differences in patterns, profiles and loci differentially methylated in eSF under the different hormonal treatments and in normal versus disease. These observations underscore key roles for genomic locations and potentially chromatin configurations further directing hormonal effects.

CpG islands (CGI), CGI shores and shelves. There was low involvement of CGIs and CGI shelves and shores in both disease stages (and ), similar to normal. But, in both disease stages loss of methylation involved more CGIs than gain of methylation.

Pathways and biofunctions associated with differentially methylated loci

Interestingly, hormone treatments significantly enriched more pathways in eSFstage-I and eSFstage-IV versus eSFnormal () ( shows important pathways affected by each hormone in each eSF group and marking differentially methylated genes in those pathways). Thus, while there were fewer loci in eSFstage-I and eSFstage-IV compared to eSFnormal, more specific pathways were affected; whereas, with eSFnormal hormone effects did not particularly affect specific canonical pathways and likely involved broader targets across the genome. E2 affected pathways in eSFstage-I involved endometrial function/dysfunction and endometriosis (e.g. MAPK, PI3K-Akt, ErbB signaling, focal adhesion, gap junctions, among others ()). E2 elicited pathways in eSFstage-IV associated with proteoglycans and estrogen, ErbB, Ras, GnRH, and FoxO signaling ()–all relevant to endometriosis pathophysiology [8,14]. While P4 did not significantly enrich specific pathways in eSFnormal, indicating a more genome-wide effect instead of limited effect at specific canonical pathways, several statistically significant pathways were enriched in eSFstage-I and even more in eSFstage-IV. These data suggest an aberrant response to P4 in eSF from women with disease, which is enhanced in stage IV () involving specific pathways including estrogen, MAPK and ErbB signaling, confirming pathways associated with transcriptomic data [19,20]. Similarly, in response to E2+P4, there were more enriched pathways in eSFendo than in eSFnormal relevant to endometrial function (adhesion, and disease/cancer ()).

Gene functional enrichments in eSFnormal, eSFstage-I and eSFstage-IV for each hormone treatment are shown in While there was little overlap in genes or functional clusters in eSFendo compared to eSFnormal, the greatest number of genes in the same functional cluster in all eSF groups induced by E2 involved those with signal peptide, membrane, and glycoprotein functions. eSFstage-I had more gene functional clusters with specific functions in all hormone treatments compared to eSFstage-IV or eSFnormal. eSFstage-IV had the fewest functional clusters in all treatments compared to eSFstage-I and eSFnormal, and the genes affected specific pathways involved in endometriosis and cancer as observed in the pathway analysis (above). Importantly, P4 treatment affected specific gene functions in disease (adhesion, synapse, cell junction, cadherins), different from eSFnormal. E2+P4 elicited, only in eSFstage-I, several distinct functional clusters with specific functions in endometrial biology and endometriosis, including EGF/EGF-like genes, ECM receptor interaction, focal adhesion, PI3K-Akt pathway, synapse, cell junctions, spectrins and others. The most enriched cluster elicited by E2+P4 in eSFstage-IV included fibronectins (large glycoproteins in ECM that bind integrins and other matrix components with major roles in cell adhesion, growth, migration differentiation, fibrosis and cancer).

Functional enrichment differences did not show a gradual change from eSFnormal to eSFstage-I and then to eSFstage-IV, rather showed distinct enrichments, suggesting inherent differences between disease stages. These data support that stage I and stage IV belong to distinct disease subtypes.

Aberrant hormone-induced methylation in eSFendo are due to pre-existing methylation abnormalities

As patterns, profiles, pathways and gene functions differed in responses of eSFstage-I, and eSFstage-IV to E2, P4 and E2+P4 in comparison to those of eSFnormal, the question arose whether these could be due, in part, to aberrant DNA methylation signatures present prior to hormone treatment (referred to “pre-existing differences” herein). The DNA methylation status of untreated (vehicle) eSFendo from women with endometriosis were assessed for loci with aberrant methylation changes in response to each hormone, compared to eSFnormal and whether they differed from untreated (vehicle) eSFnormal (see Methods). Numerous aberrantly differentially methylated loci in disease were found to be due to pre-existing DNA methylation differences across the genome (), including up to 53% of aberrant E2+P4 loss of methylation in eSFstage-IV, showing an aberrant methylation pattern from that of the untreated normal eSF, further resulting in aberrant response to hormone treatments. These data are supported by previous gene expression analysis of eSFnormal and eSFendo at t = 0 in culture [19], demonstrating intrinsic and pre-existing abnormalities in the eSFendo cells, although unclear whether these aberrancies are due to disease or are contributing to its progression/pathogenesis.

As in eSFnormal, transcriptomic profiles were determined for both stages of disease for each hormone (versus control). Whether these responses were abnormal was investigated compared to normal. In eSFendo, all treatments resulted in fewer differentially expressed genes versus eSFnormal (. P4 alone induced the fewest differentially expressed genes in eSFstage-I, and eSFstage-IV, consistent with the DNA methylation changes, an aberrant response to P4 and abnormal decidualization in eSF from women with endometriosis () [18,19]. In general, genes differentially expressed in response to each hormonal treatment were different in each disease stage versus normal, although some were in common (). Of interest to endometrial function and in disease, PGR was also upregulated in response to E2 in eSFstage-I, and eSFstage-IV similar to eSFnormal () and was also among the top up-regulated genes in all groups despite the aberrant and very limited P4 response observed in the DNA methylation, gene expression and IGFBP1 production in eSF from women with disease (). E2, P4 and E2+P4 up-regulated IGF1 and IL1R1 in all eSF and SPARCL1 was up-regulated in all E2+P4 treated eSF ( (full list and Venn diagrams of unique and common genes within each group and each hormone treatment across groups)). Gene expression profiles in response to hormones, similar to the DNA methylome, demonstrated distinct and aberrant molecular signatures in eSFstage-I versus eSFstage-IV and compared to eSFnormal. Moreover, these were not limited to P4 and E2+P4 treatments and, importantly, included an abnormal response of eSFstage-I and eSFstage-IV to E2.

Changes in gene expression in stage I disease

In eSFstage I, similar to eSFnormal, all hormone treatments resulted in more gene up-regulation than down-regulation, but unlike eSFnormal, all hormones, including E2 affected fewer differentially expressed genes, with a marked minimal effect with P4 (, full gene list; , common genes and Venn diagrams). Also, in E2 and in E2+P4 half and the majority of the genes, respectively, were unique, but in P4, the majority of differentially expressed genes were in common with E2+P4. Note that the number of P4-induced differentially expressed genes were very limited in stage I, while the combination of E2+P4 in stage I disease resulted in more differentially expressed genes than with E2 or P4 treatments alone (). Similar to eSFnormal, where addition of E2 minimally affected P4 target genes, E2 combined with P4 affected the target genes of E2 alone. While there were commonly upregulated P4 target genes eSFstage I and eSFnormal including IGF1, GREB1, and PGR, many key genes were missing in eSFstage I, such as, SPARCL1 which was upregulated in normal but not in stage I disease, further indicating an aberrant E2 response in disease.

Changes in gene expression in stage IV

Similar to eSFstage I the number of differentially expressed genes in response to E2 as well as to P4 treatments were far fewer than what was observed in eSFnormal (, full gene list). Similar to eSFnormal and eSFstage I there were more differentially expressed genes by E2+P4. Among commonly up-regulated E2+P4 induced eSFstageIV and eSFnormal were SPARCL1, IGF1, and LAMA3. Among the 103 genes commonly down-regulated were CCL2, RGS4, RGS5, IL-6, MEST, KRT19, KRT18 and H19. The overlap in differential expression of specific up- and down-regulated genes with E2+P4 treatment of eSF from women with and without endometriosis is remarkable, since stage IV disease eSF cells did not decidualize and are not considered to be P4-responsive [24,25].

Pathways and biofunctions

Pathways and biofunctions, derived from the gene expression data, underscored distinct differences between eSFendo and eSFnormal and between stages of disease, similar to the DNA methylation data. In eSFstage-I, with more limited E2 effects, pathways included activation of cellular proliferation and viability (), and eSFstage-IV involved increased tissue and cellular development (as with eSFnormal), proliferation, cell-cell signaling and adhesion (unlike eSFnormal). Note that E2-induced biofunctions and pathways were different in eSFstage-I and eSFstage-IV and both differed from eSFnormal (), consistent with the DNA methylation data. There were no enriched pathways in eSFstage-I and moderately enriched (Z score = 1.9) up-regulation of cell growth and proliferation in eSFstage-IV in reponse to P4 (far fewer loci). Similar to eSFnormal, carbohydrate metabolism and molecular transport was also seen in eSFstage-I in response to E2+P4, which also showed up-regulation of cell invasion and viability. Importantly, E2+P4 increased cell survival, cell movement and invasion, cell-to-cell signaling and adhesion, and downregulated cellular proliferation and growth in eSFstage-IV. Genes involved in these pathways and their upstream regulators are shown in . Similar to normal, the main biofunctions of differentially expressed genes that were shared with differentially methylated loci involved cell membrane and signaling in response to E2 in eSFstage-IV.

In loci with a strong positive or negative association of DNA methylation and gene expression (by Spearman rho, and corrected p<0.05) distinct differences were found in eSFstage-I and eSFstage-IV and versus eSFnormal (full lists, for unique and common loci between each group). Functional enrichment analyses revealed distinct differences in numbers and types of gene functional clusters in each stage of disease versus normal. While there was a large number of functional gene clusters with strong association for E2 in eSFnormal (), eSFstage-I and eSFstage-IV showed different and more limited functional clusters. This result further indicates that the E2 response is aberrant in eSFendo compared to normal eSF.

There were also multiple differences in response to P4, and E2+P4 among the eSF groups, further highlighting distinct molecular signatures in each disease stage. Importantly, eSFstage-I showed distinct clusters in response to P4 and to E2+P4, including functions characteristic of endometrial biology and endometriosis pathophysiology ().

While there were no statistically significantly enriched gene functional clusters in response to P4 in stage IV disease, a moderate enrichment of fibronectins, cell adhesion and secreted proteins were noted. These are consistent with the important role for cell adhesion in stage I and stage IV disease. In response to E2+P4, eSFstage-IV showed enrichment of calcium channels and integrins. Together these data suggest that the responses to all hormone treatments are aberrant in eSF derived from women with stage I and stage IV disease and are specific to each stage, supporting distinct disease subtypes.

Comparison of in vitro versus in vivo data

Herein, eSF hormone treatments in vitro were chosen to approximate the hormonal milieu in vivo (E2-dominant proliferative phase endometrium (PE) and E2+P4-dominant mid-secretory phase endometrium (MSE)). Comparing in vitro hormone eSF transcriptome data to corresponding phases in bulk endometrial tissue in normal versus disease [26] and FACS-isolated eSFendo and eSFnormal [20] revealed great overlap of differentially expressed genes (). GO functional analysis of genes differentially expressed in E2 treated eSFstage-IV and PE tissue revealed many genes in common involved in regulation of cell migration and motility, proteolysis, negative regulation of cell death, regulation of fibroblast proliferation and others. Regulation of inflammatory response, cell migration/motility, transport, protein import to nucleus, signal transduction, wound healing, and epithelial development and others were noted in stage I ().

Comparing MSE tissue and E2+P4-treated eSFstage-I and eSFstage-IV vs normal, pathway analysis revealed common signaling pathways involving PI3K-Akt, Rap1 and Ras, and cancer ().

Comparing transcriptomes of cultured eSFendo and eSFnormal and freshly isolated (uncultured) FACS-sorted eSFendo and eSFnormal from human eutopic endometrium [20] revealed many shared genes (; ). Note that FACS-sorted eSFnormal and eSFendo included samples from various cycle phases, different disease stages [20], and a more limited sample number compared to the bulk tissue study. Thus, the number of overlapping loci in cultured and freshly isolated eSF is expectedly smaller than those shared with whole tissue.

Together, the extent of overlap with whole tissue samples and FACS isolated eSF indicates the in vitro hormonal treatment of eSF, the predominant cell type in endometrium, is a good model and reflects a persistent eSF signature in the whole tissue.

Histone H3K27me3 and H3K27ac modifications in response to E2

Since E2 induced the largest changes in the DNA methylome of eSF, we sought to assess its effect on the histone marks to better understand how E2 affected the regulatory function of the epigenetic machinery. We assessed silencing and activating histone modifications, H3K27me3, and H3K27ac, using chromatin immunoprecipitation followed by deep sequencing (ChIP-Seq). Modifications of H3K27me3, and H3K27ac have been found in loci involved in eSF decidualization [10,11,27]. In response to E2 we observed more differential peaks in H3K27ac than H3K27me3 in line with our observation of more loss of DNA methylation corresponding to a more open chromatin state induced by E2 ( for peaks associated with each histone mark). GO gene functional analysis for each histone modification renriched pathways related to regulation of signaling, cell morphogenesis and differentiation, G-protein coupled receptor signaling, regulation of mitotic cell cycle and intracellular protein transport among others, many of which are shared with DNA methylation data ( for pathways for each histone mark).

Association of PGR target loci identified in E2+cAMP+MPA decidualization and E2+P4 induced DNA methylation

Increased binding of PGR to open chromatin was shown previously in decidualizing cells by ChIP-Seq experiments [10,27] and that the presence of PGR binding site and its putative co-regulator FOSL2 in a genomic location is associated with open chromatin during decidualization [10,11]. Moreover, direct PGR targets in eSF treated for 72hrs with E2+MPA+cAMP were identified by Mazur et. al., using ChIP-Seq and RNA-Seq [10]. We assessed the overlap of the E2+P4 induced differentially methylated CpG sites associated with genes in normal, stage I and stage IV disease to genes with PGR binding sites present in the Mazur et. al. study within the extended promoter region (as defined to be -7500bp and +2500bp from TSS) and intervals within ±10KB, as well as ±25KB from transcriptional start/stop site in normal eSF. We found a small subset of genes overlapping in normal eSF (); however, these common genes were enriched for biofunctions that are involved in cell morphogenesis, differentiation and cell projections, endosome organization and cytoskeletal organization (). These are important during decidualization as eSF decidualization is characterized by morphological changes, expansion/restructuring of extracellular matrix, surface projections and expansion of endoplasmic reticulum. Interestingly, a larger number of genes in stage I and IV overlapped with PGR binding sites than normal eSF (). Stage IV and normal shared more common genes with PGR binding sites than they did with stage I (). Biofunction analysis showed more biofunctions involved in stage I than normal, such as tissue morphogenesis, response to TGF-beta signaling, response to growth factor and extracellular matrix among others (). In stage IV, the biofunctions involved negative regulation of Wnt signaling, and intracellular signal transduction () known to be affected in endometriosis. These data further support the notion of aberrant P4 response, rather than P4-resistance in endometriosis.

Discussion

Unique steroid hormone effects on normal endometrial stromal fibroblasts

The eSF is the most abundant steroid hormone-responsive cell in endometrium and is a master regulator of tissue function and pregnancy success, and thus how the steroid hormones E2 and P4 regulate the epigenome and transcriptional machinery in this cell type in a timeline similar to in vivo exposure is of high priority in understanding normal and abnormal eSF function in women. Effects of E2 and P4 alone on the hormone-epigenome interplay has largely been studied in breast cancer cell lines, providing key insights into hormone receptor topology, epigenetic genomic alterations, transcriptional regulation, and chromatin dynamics [28,29]. As these complex interactions are cell- and tissue-specific, extrapolating their properties to normal endometrium is limited, although a few studies have investigated the effects of E2 plus progestins, such as medroxyprogesterone acetate, in the presence of cAMP on eSF for 72 hrs [10] or longer with or without glucose in the culture medium on chromatin accessibility or histone marks [9,11]. These studies show altered chromatin accessibility in eSF decidualized with E2+MPA+cAMP [11] and provide significant insights into PGR binding across the genome and association with open chromatin [10].

The current study investigated the effects of estradiol, progesterone (individually) and their combination (without cAMP or other progestins) on the DNA methylome and transcriptome and their interplay in normal eSF at 14 days, mimicking in vivo exposure times. Moreover, we compared eSF from the inflammatory disorder, endometriosis, in the setting of lesser and more advanced stage disease to the normal eSF. The data herein revealed, for the first time, that E2 and P4 individually and together promote unique patterns and profiles in the normal DNA methylome of this cell type. E2 alone elicited broad changes, blunted by P4, and mostly result in open chromatin by inducing more loss of methylation and increased H3K27ac histone mark. Progesterone alone had a limited effect on the DNA methylome, and unlike E2, elicited loss and gain of methylation equally. E2+P4 affected the epigenome less robustly than E2 alone, but showed more loss than gain of methylation. In support of our observation Vrljicak et. al. using transposase accessible chromatin followed by sequencing (ATAC-Seq) found altered chromatin accessibility with more open than closed chromatin loci after 4days treatment with MPA and cAMP [11]. These data suggest that E2 and P4 interact differently with the epigenome when combined than when individually administered suggestive of different mechanisms involved in the response of eSF to E2 and to P4 individually and in combination. Hormone-specific patterns and profiles were abnormal in both disease stages with more severe abnormalities associated with stage IV disease. The range of differences in individual loci with differential methylation and the number of enriched clusters and gene functions induced by each hormonal treatment in disease versus normal suggest inherent differences in disease and disease stages. However, in disease, as in normal, E2 induced more extensive alteration than E2+P4 followed by P4. Despite these differences, hormone-induced changes overall mainly involved CpG sites at the 3’end, intergenic regions, and enhancers, limited involvement of 5’end and 1st exons and rarely involved CpG sites in close proximity to transcription start sites (TSS200) or CpG islands. Notably, E2 treatment of MCF-7 breast cancer cells also demonstrated minimal binding of ER to proximal promoter regions (up to 5kb) [30], despite their containing the majority of known EREs.

Whereas CpG sites in CpG islands (CGI) were minimally affected, CGI shores and shelves were more involved, regardless of methylation loss or gain or the type of hormone treatment, indicating a specific genome landscape interaction of hormones in this endometrial cell type. Whether the lack of involvement of CGIs reflects regulation of genes whose functions are not regulated by direct or indirect hormone-targeted mechanisms, or whether hormone response elements are not affected in CGIs is yet to be determined. In breast cancer cell lines gene expression [31] and DNA methylation profiles [32] as well as DNA methylation at several candidate genes at their CGIs [33] depend on their ER and PR status. This observation further highlights findings herein that the majority of differentially methylated loci in eSFnormal are located in the intergenic regions, 3’UTRs and enhancers, and do not involve regions in close proximity to TSS, 5’UTR and 1st exons, where most CpG islands are located.

Epigenetic signatures in endometriosis: Hormone response, disease subtypes, pre-existing abnormalities

Women with endometriosis have high prevalence of infertility with otherwise unknown etiologies and lower implantation, clinical pregnancy and live birth rates compared to those without disease [34]. Studies in humans [35] and animal models [36] suggest compromised implantation attributed, in part, to an abnormal response to P4 and the inflammatory milieu of the endometrium. The current study confirmed abnormal P4-regulated decidualization marker expression in eSFendo, largely attributed to “P4-resistance”, although P4 did have effects across the eSFendo genome and PGR targets. However, eSFendo additionally had different responses to E2 compared to eSFnormal, which likely also contributes to abnormal endometrial function in women with disease (as described below). Of note, aberrant lack of ERα down-regulation at the time of implantation in endometrium of endometriosis women is considered key in implantation failure in women with disease [37]. However, endometrial-based infertility and effects on pregnancy outcomes are controversial, as large studies on IVF/ICSI outcomes in women with endometriosis and ovarian endometriomas revealed no differences in pregnancy rates [38,39] or a significant difference in endometrial receptivity array test in women with endometriosis versus controls [40]. How the aberrancies observed herein in P4-, or E2-induced epigenetic signatures are linked with implantation outcomes in women with endometriosis warrants further investigation.

Whether stage I and stage IV endometriosis are distinct disease sub-types has been the subject of debate. That eSFstage-I differ greatly from eSFstage-IV in hormone response supports distinct disease subtypes. Also, disparities in the DNA methylomes of eSFstage-I and eSFstage-IV before hormonal treatments further support distinct disease subtypes. The latter observation underscores pre-existing abnormalities in the eSF epigenome in the setting of endometriosis, and the data showed eSFstage-IV with more extensive pre-existing differences affecting its responses, compared to eSFstage-I.

This is consistent with previous findings that the endometrial bulk tissue transcriptome differs between the two stages [18,26] and several endometriosis genome-wide analysis studies suggestive of a stronger genetically driven component for stage IV than stage I disease [41]. Of note, clinically, women with stage IV versus stage I endometriosis have significantly lower implantation rates (13.7% vs. 28.3%, respectively), pregnancy rates (22.6% vs. 40.0%, respectively) [42], and lower IVF pregnancy rates (13.84% vs. 21.12% respectively) [43]- believed due to endometrial abnormalities that reflect distinct subtypes of disease. Mapping hormone-genome interactions of these subtypes holds promise for innovative, targeted therapies to modify pre-existing and stage-specific abnormalities in endometrium of women with endometriosis and optimize endometrial receptivity for implantation and pregnancy success of women with endometriosis.

Dyson et. al. [44] have observed aberrant DNA methylation in the stromal fibroblast isolated from the endometriotic ectopic lesion. It remains to be determined whether and to what extent the ectopic lesion aberrances stem from the eutopic endometrial stromal fibroblasts.

Recently, Maekawa et. al. assessed the genome-wide methylome changes during decidualization and in contrast to our data reported no changes in the DNA methylome [45]. DNA methylation distribution follows a bimodal distribution with the majority of CpG sites either hypomethylated or hypermethylated [46], as also reported in their study as well as in the current study and as we have previously observed in normal and endometriosis endometrium [14,15] Furthermore, we have also observed that the majority of CpG sites remain unchanged in decidualized versus non-decidualized eSF. The differences in our observation could be due to different analyses methods, where Maekawa et. al limited the definition of differential methylation to >Δβ of 0.3, which would not detect smaller changes. In our analyses we considered smaller changes in the DNA methylome but with the stringency that they were observed in at least 75% of each sample group. Another reason could be due to differences in the samples, where we used normal controls while patients with myoma or cervical cancer were used in their study, or it is likely that E2+MPA used in that study affects the epigenome differently than E2+P4 in our study.

Progesterone “resistance”

Pursuing bulk tissue transcriptomic analysis, we first described “P4 resistance” in endometrium from women with endometriosis [8,35] a phenomenon also observed by others [47-49]. In samples obtained in the implantation window and timed to the LH surge, there was evidence for impaired expression of key epithelial and stromal fibroblast markers of embryo receptivity and decidualization, respectively [35]. Analysis of endometrium across the menstrual cycle from women with severe disease strikingly revealed persistent E2-regulated genes in the early secretory phase, consistent with impaired P4 action [8]. Moreover, these data were substantiated in a larger cohort [26], that also revealed a marked pro-inflammatory phenotype within the endometrium of women with disease. Inflammation can cause epigenetic changes in endometrium as demonstrated in an animal model of the disease [50]. We and others found P4-resistance in eSF [19,51]. Notably, inflammatory cytokines (e.g., IL-1β and TNFα) epigenetically silence the eSF PR, promoting P4-resistance with diminished expression of decidualization markers IGFBP1 and prolactin [50,52] and enhanced secretion of matrix metalloproteinases, which are normally suppressed in eSF by P4 [53]. Epigenetic mechanisms underlying P4-resistance in endometriosis have mostly focused on the disease itself (as opposed to the eutopic endometrium studied herein) which exhibits P4 and progestin resistance for pain relief initially or acquired over time [54].

We suggest that the nomenclature of “P4 resistance” be re-evaluated, since the data herein show eSFendo respond to P4 with regard to epigenetic marks, PGR target sites, and gene expression, albeit differently from eSFnormal, although they do not fully decidualize. Notably, endometrium of women with disease does not retain a proliferative phenotype throughout the cycle [6,8], although there is compromised implantation in women with disease [35,52]. Importantly, “P4 -resistance” was observed in endometrium of non-pregnant women who previously had severe pre-eclampsia, and this was also found in the decidua at delivery of women with this disorder [55], underscoring the need to understand this process for normal pregnancy. Thus, P4 signaling in the endometrium and aberrancies in decidualization therein in vivo are likely influenced by other cell types in the tissue, including the inflammatory status of the individual and warrant further investigation.

Abnormal response to E2

Herein, for the first time the observations have been made that in addition to aberrant eSFendo P4 response, eSF from women with endometriosis show vastly aberrant responses to E2. Specifically, E2-induced eSF DNA methylation changes blunted in stage IV disease and were more extensive in stage I. Loci with strong associations of DNA methylation and gene expression had distinct enrichment in gene functions in stage I and stage IV, including ion channels, ATP and nucleotide binding in stage I and plasma membrane and signaling in stage IV, suggesting functional impairment of eSF from women with versus without endometriosis. Moreover, they underscore that not only is the eSF response to P4 abnormal in women with disease, but also their response to E is abnormal. The latter has received little attention in the endometriosis literature, which is surprising, as the disorder is estrogen-dependent [12,21]. As eSF normally require E2 priming prior to the full decidualization P4 response, it is not unanticipated that with abnormal E2 signaling in eSF, P4 signaling would also be disrupted. Aromatase, essential for E2 production, as well as E2 levels are highly expressed in endometriotic tissue [56-59] with an increased COX2 expression in turn resulting in increased E2 production in a positive feedback loop in ectopic and eutopic tissue of endometriosis patients [57]. Whether the aberrant response to E2 observed herein could be affected by these aberrancies in E2 production in endometriosis remains to be determined.

Potential mechanisms of E2-epigenome interactions

The binding of some hormone NRs commonly occurs at accessible regions of the chromatin before hormone induction [60] or their recruitment occurs almost equally at the nucleosome-occupied and nucleosome-free states before hormone induction [61]. E2 (biological active estrogen) enters cells, binds to subtypes of ERα and β that have high affinity for E2 and are encoded by different genes [62]. While both ER subtypes are expressed in human endometrium, ERα is the primary mediator of E2 action in this tissue [63]. ERα recruitment is complex involving multiple mechanisms depending on cell type and culture conditions [64]. ERα can bind to compact chromatin while there are abundant accessible regions before E2 induction that will further recruit ERα [65]. The DNA methylation and histone modification findings herein suggest that E2 can increase open chromatin. Chromatin accessibility can be induced by ERα binding, as these accessible clusters are found near estrogen-target genes [66]. While increased open chromatin was found in eSFnormal in response to E2, the opposite was found in eSFendo. This could be due to pre-existing abnormalities in disease affecting chromatin structure, a combination of transcriptional machinery preloaded across the genome, or, as found in disease, up to 50% of loci displaying pre-existing differences in epigenetic signatures influencing this response. Furthermore, the state of chromatin compaction may play an important role. About half of EREs are in regions of DNA with open chromatin prior to estrogen induction [67], but many ERα binding sites in open chromatin are associated with differentially expressed genes after estrogen induction. These data indicate that chromatin compaction can directly affect ERα recruitment and subsequently target gene transcription. Thus, pre-existing and distinct differential methylation observed in stage I and stage IV can potentially affect chromatin compaction in patients with endometriosis. This is currently under investigation in our laboratory.

ERα can be activated by phosphorylation by growth factors binding to tyrosine kinase receptors such as EGFR [68], which were dysregulated in the current study. Genes targeted by phosphorylated ERα are distinct from those targeted by estrogen-induced ERα activation [69]. Signaling through EGFR is a key pathway in eSF response to E2, and constitutive activation of EGFR in eSF from women with endometriosis has been reported [70]. Inhibition of EGFR in eSF from women with disease restores decidualization markers [71], underscoring the complexity of the interplay between E2 and P4 signaling in eSF in endometrium of women with endometriosis. The data overall support phosphorylation of ERα in eSF treated with E2 may contribute, in part, to differential DNA methylation signatures and gene expression profiles observed in E2 versus E2+P4 in women without and with disease, which remains to be proven experimentally. Interestingly, E2 and EGF can induce ERα recruitment at three classes of enhancers [72], bound only with EGF stimulation, only with E2 stimulation, or either. Herein, enrichment of enhancer involvement upon E2 stimulation and with E2+P4 was observed in normal eSF as well as eSF from women with both stages of disease, but with different effects on downstream target genes. We propose that even small differences in EGFR signaling pathways could greatly alter the eSF responses to hormones, as observed herein.

Study strengths and limitations and future directions

In this study, effects of E2, P4 and their combination were elucidated on genome-wide DNA methylation marks of the endometrial stromal fibroblast, the predominant cell type in human endometrium essential for establishing and continuing pregnancy. The clinical phenotyping of truly normal controls and specific, well-phenotyped disease stages is a great strength of this study. Also, using the same cells for DNA methylation and gene expression analyses also added to the robustness of the data. Moreover, comparisons of the data herein with published gene expression and DNA methylation data in bulk tissue underscore signatures in the latter due to this predominant cell type. Single cell analysis of eSF from bulk tissue by FACS further underscores the signature of this cell type in overall bulk tissue analyses and opens the door for single cell RNA-Seq and DNA-me analyses in the future.

While use of an in vitro system can address whether/how steroid hormones directly affect the eSF epigenome, an in vivo model using freshly isolated, sorted endometrial cells would offer an opportunity to assess functionality of the ER and PR landscape in human eSF and other endometrial cell types. Primary epithelial cells were not viable in culture when treated with hormones and as such we did not expand the current study involving epithelial cells, but organoid systems may offer a tool for this investigation. Further analysis using chromatin immunoprecipitation followed by deep sequencing of ER and expanding on the PR binding sites by Mazur et.al [10], identification of ERE and PRE-specific to endometrium and the status important to pioneer and co-activators and in a larger sample size are required for detailed mapping of the steroid hormone landscapes and hormone-epigenome interplay in normal human endometrial cells and in disease. Transcriptome data from the same cells demonstrate extensive overlap with previously identified differentially expressed genes in whole tissue in the corresponding hormone milieu. We note that utilizing microarray instead of a more comprehensive transcriptome analysis such as RNA-Seq limited the number and the type of transcripts investigated herein. Furthermore, protein data will enable full assessment of epigenomic and transcriptomic effects of hormones in endometrial function normally and in women with endometriosis. An important limitation of this study is the small sample size in each group. We had used very strict criteria for sample selection, both in identifying normal samples without any gynecological and pelvic disease/disorder and for endometriosis to not have any other disease no matter how benign, such as uterine fibroids. A follow up study with a much larger sample size is required to confirm our observations.

Overall this study has elucidated the array of responses of eSF in health and disease in hormone milieu encountered in cycling women that can also serve for comparisons with actions of pharmaceutical steroids used clinically and potentially environmental estrogens that can compromise reproductive function. Moreover, the data reveal unique responses and pre-existing epigenetic abnormalities in women with endometriosis that can benefit endometrial-based diagnostic development and novel targeted therapies for endometrial dysfunction in women with this disorder.

Materials and methods

Ethics statement

This study was approved by the Committee on Human Research of the University of California, San Francisco (UCSF) (IRB# 10–02786). All samples were collected after written informed consent was obtained from all subjects.

Samples

Eutopic (within the uterus) endometrial tissue samples were collected through the UCSF/NIH Human Endometrial Tissue/DNA Bank. Stringent inclusion criteria were applied as follows: I) for normal controls, samples were collected from oocyte donor volunteers with no uterine or pelvic pathology (NUP, normal controls); endometriosis samples were collected from stage I and stage IV endometriosis patients (). Oocyte donor volunteers (controls) were extensively screened, had no gynecologic disorders, and donated endometrial samples six months post oocyte retrieval. Endometriosis patients (stage I and IV) were surgically confirmed and had no other gynecologic abnormalities. II) All samples were collected in the proliferative phase and matched for age, BMI, no smoking history (one exception), no contraceptive steroid use three months prior to sample collection, and endometrial stromal fibroblast (eSF) passage number. Menstrual cycle phase was determined by histological evaluation [73] as well as serum levels of E2 and P4. Disease stage was determined by ASRM criteria [74].

Stromal cell isolation and hormone treatment

Primary eSF were isolated from endometrial biopsies by digestion with collagenase and size fractionation as described [75] and cultured as monolayers in stromal cell medium (SCM, [18,19,76]). To ensure the purity of stromal cells in culture, after digestion of endometrial sample biopsies, the digested tissue was size fractionated using a 40μM filter to separate epithelial glands, followed by stromal cells selective attachment and growth in stromal cell medium. The purity of primary eSF was monitored morphologically during the culture and the homogeneity was verified by immunocytochemical localization of vimentin for eSF, keratin for epithelial cells and actin for vascular cells [76] before further hormone experiments. Only pure primary eSF (with <0.1% other cells) were used for this study. After 24 hours serum starvation, eSF from normal women (n = 7, controls) and endometriosis women (n = 6 stage I, n = 9 stage IV) were treated with four different hormonal treatments of 10 nM E2, 1μM P4, 10 nM E2 +1μM P4, or vehicle (0.1% ethanol) control for 14 days [76]) after which conditioned media and cells were collected for further analysis. Decidualization was assessed in E2+P4 treated eSF from normal, stage I and stage IV disease (see below). As eSF from women without endometriosis mostly have a robust decidualization response to E2+P4, our controls were eSF that fully decidualized (n = 4) by the decidualization marker IGFBP1 by ELISA and morphologically. As rarely do eSF from endometriosis women decidualize in vitro in response to E2+P4, eSF from stage I (n = 4) and stage IV (n = 4) with non-detectable decidualization (the most common phenotype) by morphology and IGFBP1 marker by ELISA were used for further analysis.

Decidualization assessment

Insulin-like growth factor binding protein-1 (IGFBP1), a P4-induced decidualization marker [77], was measured in media conditioned by 14 day E2+P4 treated cultures, by ELISA (Alpha Diagnostic International Inc., San Antonio, TX) as a marker for decidualization. Concentrations were measured in duplicate, averaged and normalized to cell number. eSF were assessed by microscopy for morphological changes corresponding to decidualization.

DNA and RNA extraction

After treatments cells were harvested, pelleted and frozen at -80C for DNA and RNA extraction as previously described [14,15]. Genomic DNA was extracted using QIAGEN (QIAamp DNA Tissue Kit, QIAGEN, Germantown, MD) and RNA was extracted using the Macherey-Nagel NuceloSpin Tissue Kit with DNase treatment (Macherey-Nagel Inc., Bethlehem, PA) according to manufacturers’ recommendations and stored at -80C.

DNA methylation

Genomic DNA was bisulfite converted at the University of Southern California (USC) Epigenome Center using the EZ-96 DNA Methylation Kit (Zymo Research, Irvine, CA), according to the manufacturer’s protocol, and as previously described [14,15]. Quality, completeness of bisulfite conversion and amount of bisulfite-converted DNA were assessed by a panel of MethyLight reactions [78]. All samples passed all quality controls (QCs) and were further assayed by the Illumina Infinium HumanMethylation450K DNA methylation platform (HM450) based on Illumina’s specifications. (All data files are submitted to GEO, under SuperSeries accession number GSE145702).

Gene expression microarray analysis

RNA quality was assessed by Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA). RNA samples were prepared for microarray analysis according to Affymetrix (Affymetrix, Inc., Santa Clara, CA) specifications [15]. cDNA sample quality was assessed by Bioanalyzer, and samples passing QCs were hybridized to Affymetrix HU133 Plus 2.0 gene array, interrogating >38,500 genes at the UCSF Genomics Core.

DNA methylation data analysis

HM450 interrogates 485,577 methylation targets across the genome. The ratio of methylated signal over total fluorescent signal was used to calculate β values, ranging from 0 (no) to 1 (complete) methylation. 850 control bead types were used as positive and negative controls and to calculate a detection P value to assess DNA methylation measurement quality for each probe of each sample [79]. The passing threshold of P value was set at P<0.05, and probes with P>0.05 were indicated as “missing” (no statistically significant differences from background). Probe dropout-rates (percent probes with missing values versus total number of platform probes) were calculated to exclude samples with dropout-rates >1%. All 12 samples passed these criteria. Probes with a “missing” value in >1 sample were removed. Differential DNA methylation in each hormone treatment in each group (control, stage I, stage IV) was identified compared to non-treated cells (vehicle control) in the same group. For each probe, β values of hormone treatment in an individual sample were compared to its corresponding vehicle treated β value (Δβ). Changes in β values (Δβ) < 10% were not considered as differentially methylated for each sample. Median or average changes between hormone treatment vs vehicle control were not used, to obviate limited numbers of strong signals affecting selecting differentially methylated loci. Instead, probes were considered differentially methylated if they exhibited >10% change in β value (Δβ >10%) in hormone treatment versus vehicle and in at least 3 of 4 samples within each group, and with the same direction of methylation change (gain or loss). To assess if the aberrant signatures observed in hormone induced changes in disease compared to normal were due to pre-existing aberrations in the non-treated cells, we compared the methylation signatures of the non-treated (vehicle) eSF in each disease stage to non-treated (vehicle) eSF in normal for each aberrantly methylated locus in response to each hormone in disease. Pre-existing differential methylation, loci differentially methylated in normal vs each stage of disease in untreated cells (vehicle) were determined. The percentage of the differences observed in hormonal treatment of disease were assessed. Association of CpG islands with enhancers, distribution across the genome, association with CGIs and CGI shores and shelves were extracted from Illumina Infinium HumanMethylation450K manifest.

Gene expression analysis

The raw.CEL gene expression data files were RMA normalized using GeneSpring (GX13.1 version, Agilent Technologies). Loci were considered differentially expressed with Benjamini-Hochberg corrected ANOVA p<0.05 and fold change (FC)≥1.5 in each comparison.

DNA methylation association with changes in gene expression

Differentially methylated loci and normalized gene expression were imported into R, and corresponding probes from each platform were matched using the transcript identifier. Every DNA methylation probe for a given locus identifier was compared to all corresponding transcripts of that locus using the non-parametric Spearman’s rank-order correlation method, as bivariate normality could not be assumed (DNA methylation data are not normally distributed). Spearman's rank correlation coefficients (ρ) on gene expression and DNA methylation were computed for each probe, along with a p value testing against the null hypothesis that ρ equals zero.

Pathways, biofunctions and genomic distribution analyses

Ingenuity Pathway Analysis (IPA, QIAGEN) software was used to determine pathways, upstream regulators and biofunctions of differentially expressed genes, as described [14,15,26]. Pathways with Z-scores ≥|2| were considered significantly enriched. For differentially methylated loci, DAVID and KEGG databases were used to identify functional classification, functional enrichment and pathways. For pathway selection an enrichment score ≥|2| with a Benjamini-Hochberg corrected p<0.05 was considered. For genomic distribution, element enrichment analyses, the null hypothesis that the observed proportions in two groups are the same, a test of proportions was performed in R using the prop.test() function.

Association of in vitro and in vivo differential gene expression

eSF treatment with E2 and E2+P4 used herein, mimicked in vivo proliferative endometrium (PE, E2-dominant phase) and mid-secretory endometrium (MSE (E2 and P4- dominant phase). To assess commonalities, the current data were compared to previously published [26] whole endometrial gene expression in PE and MSE in normal versus disease. E2-treated differential expression in eSF from stage I versus control and stage IV versus control were compared to endometriosis (all stages) versus control in bulk tissue PE, and E2+P4- treated compared to MSE. FACS-sorted eSF gene expression data [20] in disease versus control were also compared. As FACS-sorted eSFendo and eSFnormal included a mixture of phases and endometriosis stages, gene expression signatures of FACS-sorted eSF were compared to stage I, stage IV eSF treated with E2 and E2+P4.

ChIP-Seq for Histone H3K27me3 and H3K27ac in response to E2

We found that E2 affected the methylome more robustly than P4 or E2+P4 and that, unexpectedly, along with aberrant P4 response in disease, E2 response was also aberrant in both stages of disease. Therefore, we sought to investigate further the effect of E2 on two repressive and open chromatin histone marks, H3K27me3 and H3K27ac. eSF cells from two independent control participant were isolated from endometrial biopsies by digestion followed by size fractionation, and primary eSF were cultured in SCM and purity of cultured eSF was assessed as described above. eSF was passaged with trypsin and 1x105 cell/well were seeded. Confluent cultures were serum starved for 24hrs and treated with E2 or vehicle for 14 days. Cells were cross-linked by a final concentration of 1% formaldehyde and terminated after 10 minutes by 0.125 M final concentration glycine. Chromatin was extracted using Chromatin Extraction Kit according to manufacturer’s recommendation (ab117152, Abcam, Cambridge, UK) sonicated by Diagenode Bioruptor and the size of sheared chromatin was visualized on agarose gel (100-600bp). Chromatin Immunoprecipitation was done using Abcam ChIP Kit (ab117138, Abcam, Cambridge, UK) with antibodies for H3K27me3 (ab6002, Abcam, Cambridge, UK) or H3K27ac (ab4729, Abcam, Cambridge, UK). Input control and immunoprecipitated DNA were paired-end sequenced using Illumina NextSeq 500 after library preparation according to the manufacturer’s instructions. Data were analyzed by removing adapter sequences, then aligned to reference human genome. Peaks called using Macs2 callpeaks and were selected with q-value <0.05. Differential peaks were identified using Macs2 bdgdiff and log likelihood ratio >3.

IGFBP1 ELISA assay of normal, stage I and IV eSFs used in the study.

(PDF)

Click here for additional data file.

Genomic distribution of differentially methylated CpG sites in each hormonal treatment, by gain or loss of methylation, in normal (NUP), stage I (Endo I) and stage IV (Endo IV) eSF.

(PDF)

Click here for additional data file.

Heat-map of pre-existing aberrancies prior to hormone treatments and percentage of contribution to each hormone induced methylation in stage I and stage IV eSF.

(PDF)

Click here for additional data file.

Differentially methylated CpG sites in response to E2, P4, E2+P4 (vs. vehicle) in Normal eSF, based on loss and gain of methylation.

(XLSX)

Click here for additional data file.

Distribution of differentially methylated loci (in %) based on their location across the genome and in comparison to those interrogated within the platform in normal (NUP), stage I (Endo I) and stage IV (Endo IV).

(XLSX)

Click here for additional data file.

Distribution of differentially methylated loci (in %) based on their association with CpG islands (CGI), and CGI north/south shelf, CGI north/south shore and not associated with CGIs (open sea), and in comparison to the interrogated CGI in the platform in normal (NUP), stage I (Endo I) and stage IV (Endo IV).

(XLSX)

Click here for additional data file.

Differentially expressed genes in hormone treatments (E2, P4, E2+P4) vs. vehicle in normal eSF (NUP).

(XLSX)

Click here for additional data file.

Enriched pathways and the associated genes and upstream regulators of differentially expressed gene in response to hormones (E2, P4, E2+P4) vs. vehicle, in normal (NUP), stage I (Endo I), and stage IV (Endo IV) disease.

(XLSX)

Click here for additional data file.

Differentially methylated CpG sites in response to E2, P4, E2+P4 (vs. vehicle) in endometriosis stage I eSF, based on loss and gain of methylation.

(XLSX)

Click here for additional data file.

Differentially methylated CpG sites in response to E2, P4, E2+P4 (vs. vehicle) in endometriosis stage IV eSF, based on loss and gain of methylation.

(XLSX)

Click here for additional data file.

Differentially expressed genes in hormone treatments (E2, P4, E2+P4) vs. vehicle in stage I disease (Endo I).

(XLSX)

Click here for additional data file.

Differentially expressed genes in hormone treatments (E2, P4, E2+P4) vs. vehicle in eSF from stage IV disease (Endo IV).

(XLSX)

Click here for additional data file.

Unique and common up- and down-regulated genes in each hormone treatment and in normal (NUP), stage I (Endo I) and stage IV Endo (IV).

(XLSX)

Click here for additional data file.

Association of differentially methylated (DM) loci with differentially expressed (DE) genes in normal eSF (NUP) for all hormone treatments (E2, P4, E2+P4) vs. vehicle, based on positive or negative association (positive/negative rho) and gain and loss of methylation (gain, loss).

(XLSX)

Click here for additional data file.

Association of differentially methylated (DM) loci with differentially expressed (DE) genes in stage I eSF (Endo I) for all hormone treatments (E2, P4, E2+P4) vs. vehicle, based on positive or negative association (positive/negative rho) and gain and loss of methylation (gain/loss).

(XLSX)

Click here for additional data file.

Association of differentially methylated (DM) loci with differentially expressed (DE) genes in stage IV eSF (Endo IV) for all hormone treatments (E2, P4, E2+P4) vs. vehicle, based on positive or negative association (positive/negative rho) and gain and loss of methylation (gain/loss).

(XLSX)

Click here for additional data file.

Differentially methylated (DM) loci associated with differentially expressed (DE) assessed by Spearman rho, with either positive association (pos) or negative association (neg).

Common/unique genes are shown in columns and in Venn diagrams, in normal (NUP), stage I (Endo I), and stage IV (Endo IV).

(XLSX)

Click here for additional data file.

Overlap of in vivo with in vitro genes and biofunctions from hormonal treatments in culture and endometrial tissue cycle phases in endometriosis stage I or stage IV versus normal.

(XLSX)

Click here for additional data file.

Genes commonly differentially expressed in eSF from disease vs. normal, treated or untreated in vitro with hormones and in eSF FACS-isolated in disease vs. normal.

(XLSX)

Click here for additional data file.

Peaks and GO biofunctions enriched in histone marks H3K27me3 and H3K27ac in normal eSF induced by E2 versus vehicle.

(XLSX)

Click here for additional data file.

Sample information.

(XLSX)

Click here for additional data file.

Differentially methylated loci and the associated pathways/biofunctions that are affected in hormonal treatments of eSFnormal, eSFstage I and eSFstage IV with known roles/importance in normal endometrial function and dysfunction in endometriosis.

(PDF)

Click here for additional data file.

Table 1

Pathways associated with differentially methylated (DM) transcribed loci in each hormonal treatment (E2, P4, E2+P4) vs. vehicle in normal (NUP), stage I (Endo I) and stage IV (Endo IV).

Pathways in each hormone treatment	NUP (P-Values <0.05, Enriched but P>0.05)	Endo I (P-Values <0.05, Enriched but P>0.05)	Endo IV (P-Values <0.05, Enriched but P>0.05)
E2 vs. Veh	Gap junction (0.015), Long-term potentiation (0.016), Long-term depression (0.025), Melanogenesis (0.035), Glutamatergic synapse (0.038), Dopaminergic synapse (0.041), Sphingolipid signaling pathway (0.053), Ubiquitin mediated proteolysis (0.064), cGMP-PKG signaling pathway (0.07), Thyroid hormone signaling pathway (0.077), Retrograde endocannabinoid signaling (0.078), Cell adhesion molecules (CAMs) (0.08), mRNA surveillance pathway (0.092), Vascular smooth muscle contraction (0.097)	MAPK signaling pathway (0.00052), cGMP-PKG signaling pathway (0.0014), PI3K-Akt signaling pathway (0.0019), Focal adhesion (0.0025), Axon guidance (0.0049), ECM-receptor interaction (0.0058), Oxytocin signaling pathway (0.012), Melanoma 0.015 Platelet activation (0.015), ErbB signaling pathway (0.017), Gap junction (0.018), Long-term depression (0.018), Amoebiasis (0.022), AMPK signaling pathway (0.022), Olfactory transduction (0.024), Renin secretion (0.025), Choline metabolism in cancer (0.039), Arrhythmogenic right ventricular cardiomyopathy (ARVC) (0.041), T cell receptor signaling pathway (0.044), Protein digestion and absorption (0.046), VEGF signaling pathway (0.057), Cell adhesion molecules (CAMs) (0.058), Inflammatory bowel disease (IBD) (0.069), Proteoglycans in cancer (0.071), Renal cell carcinoma (0.074), HIF-1 signaling pathway (0.077), Ubiquitin mediated proteolysis (0.093), Regulation of actin cytoskeleton (0.098)	Term P-Value Melanoma (0.00043), Signaling pathways regulating pluripotency of stem cells (0.00069), Proteoglycans in cancer (0.0018), Estrogen signaling pathway (0.012), Glioma (0.012), Choline metabolism in cancer (0.013), Olfactory transduction (0.022), ErbB signaling pathway (0.032), Prostate cancer (0.034), Endometrial cancer (0.035), Ras signaling pathway (0.037), GnRH signaling pathway (0.037), FoxO signaling pathway (0.038), Non-small cell lung cancer (0.042), MAPK signaling pathway (0.064), Endocytosis (0.067), Cholinergic synapse (0.068), Fc epsilon RI signaling pathway (0.068), Oxytocin signaling pathway (0.068), Rap1 signaling pathway (0.071), Neurotrophic signaling pathway (0.085), RNA degradation (0.091), Pathways in cancer (0.095)
P4 vs. Veh	Staphylococcus aureus infection (0.071), MAPK signaling pathway (0.084), Tuberculosis (0.095), Tight junction (0.097)	Pathways in cancer (0.015), Cell adhesion molecules (CAMs) (0.02), Glutamatergic synapse (0.028), Toxoplasmosis (0.032), Cocaine addiction (0.037), Steroid hormone biosynthesis (0.056), Phototransduction (0.062), Inflammatory bowel disease (IBD) (0.071), Serotonergic synapse (0.086), Chronic myeloid leukemia (0.094)	Non-small cell lung cancer (0.00031), Prostate cancer (0.00033), Glioma (0.00062), Endometrial cancer (0.0022), Thyroid hormone synthesis (0.0065), Melanoma (0.0069), Glutamatergic synapse (0.0073), Calcium signaling pathway (0.012), ErbB signaling pathway (0.014), MAPK signaling pathway (0.019), Proteoglycans in cancer (0.02), Estrogen signaling pathway (0.021), Acute myeloid leukemia (0.022), Pathways in cancer (0.024), Long-term depression (0.026), Cholinergic synapse (0.031), Pancreatic cancer (0.032), Amphetamine addiction (0.033), Chronic myeloid leukemia (0.042), Dopaminergic synapse (0.048), FoxO signaling pathway (0.055), cAMP signaling pathway (0.06), Gap junction (0.068), Viral carcinogenesis (0.068), Hepatitis B (0.069) Adrenergic signaling in cardiomyocytes (0.071), Bladder cancer (0.071), GnRH signaling pathway (0.073), Rap1 signaling pathway (0.074), Oxytocin signaling pathway (0.089), Choline metabolism in cancer (0.094), Cocaine addiction (0.096)
E2+ P4 vs Veh	Morphine addiction (0.022), Steroid hormone biosynthesis (0.082), Long-term depression (0.089)	Focal adhesion (0.00097), ECM-receptor interaction (0.0019), PI3K-Akt signaling pathway (0.0087), Glutamatergic synapse (0.0099), Toxoplasmosis (0.012), FoxO signaling pathway (0.024), Dorso-ventral axis formation (0.033), Arrhythmogenic right ventricular cardiomyopathy (ARVC) (0.038), Hypertrophic cardiomyopathy (HCM) (0.053), Dilated cardiomyopathy (0.069), NF-kappa B signaling pathway (0.077)	Transcriptional mis-regulation in cancer (0.012), PI3K-Akt signaling pathway (0.043), AMPK signaling pathway (0.051), Epstein-Barr virus infection (0.066), Prostate cancer (0.079), Type I diabetes mellitus (0.083), Taste transduction (0.093)

Table 2

Comparison of Functional Enrichment clusters in normal, stage I and stage IV in each hormone treatment (E2, P4, E2+P4) vs. vehicle.

	Normal	Stage I	Stage IV
E2 Treatment	4 clusters (Enrichment Score ≥2)1) 255 loci: zinc/zinc ion binding, metal binding2) 576 loci: disulfide bond, signal peptide, Cell membrane, Signal, Extracellular, Glycoprotein3) 14 loci: TK4) 19 loci: Spectrins, Spectrin/alpha-actinin, Calponin homology domain, CHs, Actin-binding, Actinin-type, DNA replication, recombination, and repair	6 clusters (Enrichment Score ≥2)1) 449 loci: Cell membrane, cytoplasmic, glycoprotein, membrane, transmembrane2) 13 loci: SAM3) 48 loci: synapse/cell junction4) 20 loci: Fibronectin types5) 301 loci: disulfide bond, signal peptide, Signal, Glycoprotein6) 16 loci: PDZ	Only 1 cluster (Enrichment Score ≥2)1) 135 loci: signal peptide, Receptor, Signal, disulfide bond, Cell membrane, Extracellular, Cytoplasmic, plasma membrane, glycosylation site: N-linked, Glycoprotein, integral component of membrane, transmembrane region, Transmembrane helix, Membrane
P4 Treatment	2 (+ 1) clusters (Enrichment Score ≥2)1) 18 loci: Pleckstrin homology domain, PH2) 31 loci: PH1, PH2, guanyl-nucleotide exchange factor activity, regulation of small GTPase mediated signal transduction, positive regulation of apoptotic process, Dbl homology (DH) domain, RhoGEF, Rho guanyl-nucleotide exchange factor activity, regulation of Rho protein signal transduction, DH, Guanine-nucleotide releasing factor, positive regulation of GTPase activity, Pleckstrin homology-like domain3) (ES = 1.99) 100 loci: topological domain: Cytoplasmic, Cell membrane, plasma membrane	4 clusters (Enrichment Score ≥2)1) 25 loci: Synapse, cell junction2) 34: Cadherins, cell adhesion3) 149: in Extracellular, Cell membrane, plasma membrane, Cytoplasmic, glycosylation site: N-linked, transmembrane region, integral component of membrane, Transmembrane, membrane4) 7: SAM	2 clusters (ES 2.9 and ES 1.99)1) 24 loci: Cadherins, calcium ion binding, homophilic cell adhesion via plasma membrane adhesion molecules, Calcium, Cell adhesion2) 14 loci: postsynaptic density, Postsynaptic cell membrane, cell junction, Synapse, Cell junction
E2+P4 treatment	NO clusters of ES≥2; 2 clusters of 1.611) 126 loci: transcription factor activity, sequence-specific DNA binding, DNA binding, regulation of transcription, DNA-templated, nucleoplasm, Transcription/T regulation, Nucleus2) 3 loci: Diacylglycerol kinase, catalytic domain, ATP-NAD kinase-like domain, DAGKc	6 clusters (Enrichment Score ≥2)1) 14: PDZ, PDZ domain2) 49: EGF, EGF-like, EGF-like Ca binding3) 21: ECM-receptor interaction; Focal adhesion; PI3K-Akt signaling pathway4) 43: postsynaptic density, Postsynaptic cell membrane, Synapse, cell junction,5) 24: Spectrins, Actin/actin-binding, Spectrin/alpha-actinin6) 264:: Extracellular, cell membrane, Glycoprotein, transmembrane/ transmembrane region	E2/P4: 1 cluster (Enrichment Score ≥2)1) 24 loci: Fibronectins, Immunoglobulin-like fold

Table 3

Top up/down regulated differentially expressed genes in response to each hormone treatment (E2, P4, E2+P4) vs. vehicle in normal (NUP), stage I (Endo I) and stage IV (Endo IV).

Transcripts Cluster Id	Gene symbol	Fold change	Regulation	Chr.	Entrez gene	Gene Description
NUP_E2 vs Veh
7965873	IGF1	21.46	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8040292	GREB1	8.19	up	chr2	9687	growthregulationbyestrogeninbreastcancer1
7951165	PGR	4.78	up	chr11	5241	Progesterone receptor
7971461	LCP1	4.70	up	chr13	3936	lymphocytecytosolicprotein1(L-plastin)
8101659	SPARCL1	4.60	up	chr4	8404	SPARC-like1(hevin)
8088560	ADAMTS9	3.81	up	chr3	56999	ADAMmetallopeptidasewiththrombospondintype1motif,9
7987315	ACTC1	3.79	up	chr15	70	actin, alpha, cardiacmuscle1
8107823	ADAMTS19	3.13	up	chr5	171019	ADAMmetallopeptidasewiththrombospondintype1motif,19
8086352	ULK4	3.00	up	chr3	54986	unc-51likekinase4
8080562	IL17RB	2.91	up	chr3	55540	interleukin17receptorB
7965335	DUSP6	-2.08	down	chr12	1848	dualspecificityphosphatase6
7917561	GBP4	-2.09	down	chr1	115361	guanylatebindingprotein4
8146863	SULF1	-2.11	down	chr8	23213	sulfatase1
8097692	EDNRA	-2.18	down	chr4	1909	Endothelin receptor type A
7921916	RGS5	-2.25	down	chr1	8490	regulatorofG-proteinsignaling5
7968417	FRY	-2.32	down	chr13	10129	Furry homolog (Drosophila)
7909503	SERTAD4	-2.41	down	chr1	56256	SERTAdomaincontaining4
8067969	CHODL	-2.46	down	chr21	140578	chondrolectin
8006433	CCL2	-2.84	down	chr17	6347	chemokine(C-C motif) ligand2
7906919	RGS4	-3.58	down	chr1	5999	regulatorofG-proteinsignaling4
NUP_P4 vs Veh
8089145	ABI3BP	14.27	up	chr3	25890	ABI family, member3(NESH)binding protein
8101659	SPARCL1	13.19	up	chr4	8404	SPARC-like1(hevin)
8155864	RORB	5.17	up	chr9	6096	RAR-related orphan receptor B
7965873	IGF1	4.61	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
7977933	SLC7A8	4.31	up	chr14	23428	solutecarrierfamily7(amino acid transporter light chain, L system), member 8
8138231	THSD7A	4.11	up	chr7	221981	thrombospondin, type I, domaincontaining7A
8132694	IGFBP1	3.87	up	chr7	3484	insulin-likegrowthfactorbindingprotein1
7971461	LCP1	3.70	up	chr13	3936	lymphocytecytosolicprotein1(L-plastin)
8122099	ENPP1	3.22	up	chr6	5167	Ectonucleotide pyrophosphatase/phosphodiesterase1
8043995	IL1R1	3.21	up	chr2	3554	interleukin1receptor, type I
8046048	CSRNP3	-2.07	down	chr2	80034	cysteine-serine-richnuclearprotein3
7923978	CD34	-2.09	down	chr1	947
8166747	SYTL5	-2.10	down	chrX	94122	synaptotagmin-like5
7945680	H19	-2.17	down	chr11	283120///100033819///6206	H19, imprintedmaternallyexpressedtranscript(non-proteincoding)|microRNA675|ribosomalproteinS12
7917561	GBP4	-2.25	down	chr1	115361	guanylatebindingprotein4
8150428	SFRP1	-2.29	down	chr8	6422	secretedfrizzled-relatedprotein1
8146863	SULF1	-2.30	down	chr8	23213	sulfatase1
8102587	NDNF	-2.31	down	chr4	79625	neuron-derived neurotrophic factor
8138289	ETV1	-2.66	down	chr7	2115	etsvariant1
8006433	CCL2	-2.73	down	chr17	6347	Chemokine (C-C motif) ligand2
NUP_E2+P4 vs Veh
8101659	SPARCL1	53.10	up	chr4	8404	SPARC-like1(hevin)
8089145	ABI3BP	37.91	up	chr3	25890	ABI family, member3 (NESH) binding protein
8132694	IGFBP1	24.11	up	chr7	3484	insulin-likegrowthfactorbindingprotein1
8155864	RORB	17.01	up	chr9	6096	RAR-related orphan receptor B
8040292	GREB1	15.10	up	chr2	9687	growthregulationbyestrogeninbreastcancer1
7977933	SLC7A8	13.95	up	chr14	23428	solutecarrierfamily7(amino acid transporter light chain, L system), member 8
7965873	IGF1	13.61	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8138231	THSD7A	10.44	up	chr7	221981	thrombospondin, type I,domaincontaining7A
8144917	LPL	8.89	up	chr8	4023	Lipoprotein lipase
7971461	LCP1	8.88	up	chr13	3936	lymphocytecytosolicprotein1(L-plastin)
8055323	NCKAP5	-3.31	down	chr2	344148	NCK-associatedprotein5
8121916	RSPO3	-3.39	down	chr6	84870	R-spondin3
8129573	MOXD1	-3.49	down	chr6	26002	monooxygenase, DBH-like1
7917561	GBP4	-3.80	down	chr1	115361	guanylatebindingprotein4
7945680	H19	-3.95	down	chr11	283120///100033819///6206	H19, imprinted maternally expressed transcript (non-protein coding) |microRNA675|ribosomalproteinS12
8138289	ETV1	-4.12	down	chr7	2115	etsvariant1
8131803	IL6	-4.34	down	chr7	3569	interleukin6
8006433	CCL2	-5.52	down	chr17	6347	Chemokine (C-C motif) ligand2
8150428	SFRP1	-6.05	down	chr8	6422	secretedfrizzled-relatedprotein1
7933194	CXCL12	-6.08	down	chr10	6387	chemokine(C-X-C motif) ligand12
Endo I_E2 vs Veh
7965873	IGF1	11.37	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
7951165	PGR	3.82	up	chr11	5241	Progesterone receptor
8117020	MYLIP	3.34	up	chr6	29116	Myosin regulatory light chain interacting protein
8145766	NRG1	2.94	up	chr8	3084	neuregulin1
8111490	PRLR	2.82	up	chr5	5618	Prolactin receptor
8102950	INPP4B	2.77	up	chr4	8821	inositolpolyphosphate-4-phosphatase, type II,105kDa
7942674	TSKU	2.67	up	chr11	25987	tsukushi, small leucine rich proteoglycan
8043995	IL1R1	2.55	up	chr2	3554	interleukin1receptor, type I
8088560	ADAMTS9	2.45	up	chr3	56999	ADAMmetallopeptidasewiththrombospondintype1motif,9
8106516	JMY	2.26	up	chr5	133746	Junction mediating and regulatory protein, p53 cofactor
8135069	SERPINE1	-1.61	down	chr7	5054	Serpin peptidase inhibitor, clade E (nexin,plasminogenactivatorinhibitortype1),member1
7902565	LPHN2	-1.66	down	chr1	23266///101927458	latrophilin2|uncharacterizedLOC101927458
7921916	RGS5	-1.74	down	chr1	8490	regulatorofG-proteinsignaling5
7909503	SERTAD4	-1.74	down	chr1	56256	SERTAdomaincontaining4
7968417	FRY	-1.75	down	chr13	10129	Furry homolog(Drosophila)
7922343	TNFSF4	-1.79	down	chr1	7292	Tumor necrosis factor(ligand)superfamily,member4
8006433	CCL2	-1.90	down	chr17	6347	chemokine(C-C motif) ligand2
7906919	RGS4	-2.02	down	chr1	5999	regulatorofG-proteinsignaling4
8021081	SLC14A1	-2.37	down	chr18	6563	solutecarrierfamily14(urea transporter), member1(Kidd blood group)
7922337	TNFSF18	-2.53	down	chr1	8995	Tumor necrosis factor(ligand) superfamily, member18
Endo I_P4 vs Veh (All DE loci)
7964834	CPM	2.06	up	chr12	1368	Carboxypeptidase M
8089145	ABI3BP	1.90	up	chr3	25890	ABI family, member3(NESH) binding protein
8043995	IL1R1	1.89	up	chr2	3554	interleukin1receptor, type I
8157524	TLR4	1.88	up	chr9	7099	toll-likereceptor4
8101659	SPARCL1	1.83	up	chr4	8404	SPARC-like1(hevin)
7907271	FMO2	1.83	up	chr1	2327	flavincontainingmonooxygenase2(non-functional)
7933204	C10orf10	1.65	up	chr10	11067	chromosome10openreadingframe10
7965873	IGF1	1.65	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8122660	UST	1.59	up	chr6	10090	uronyl-2-sulfotransferase
8052355	EFEMP1	1.55	up	chr2	2202	EGFcontainingfibulin-likeextracellularmatrixprotein1
7969861	ITGBL1	1.51	up	chr13	9358	integrin, beta-like1(with EGF-like repeat domains)
8100154	CORIN	1.50	up	chr4	10699	corin, serine peptidase
7922337	TNFSF18	-1.56	down	chr1	8995	Tumor necrosis factor (ligand)superfamily, member18
Endo I_E2+P4 vs Veh
8101659	SPARCL1	23.12	up	chr4	8404	SPARC-like1(hevin)
7965873	IGF1	12.49	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8089145	ABI3BP	10.15	up	chr3	25890	ABI family, member3 (NESH) binding protein
8155864	RORB	8.92	up	chr9	6096	RAR-related orphan receptor B
8147030	STMN2	7.62	up	chr8	11075	stathmin2
8043995	IL1R1	6.73	up	chr2	3554	interleukin1receptor, type I
8125919	FKBP5	5.05	up	chr6	2289	FK506bindingprotein5
7977933	SLC7A8	4.04	up	chr14	23428	solutecarrierfamily7(amino acid transporter light chain, L system), member 8
8122099	ENPP1	4.04	up	chr6	5167	Ectonucleotide pyrophosphatase/phosphodiesterase1
7964834	CPM	3.96	up	chr12	1368	Carboxypeptidase M
8138289	ETV1	-2.32	down	chr7	2115	etsvariant1
7921916	RGS5	-2.34	down	chr1	8490	regulatorofG-proteinsignaling5
8081001	ROBO2	-2.34	down	chr3	6092	roundabout,axonguidancereceptor,homolog2(Drosophila)
8136248	MEST	-2.44	down	chr7	4232	Mesoderm specific transcript
7922343	TNFSF4	-2.51	down	chr1	7292	Tumor necrosis factor (ligand)superfamily, member 4
7906919	RGS4	-2.76	down	chr1	5999	regulatorofG-proteinsignaling4
8006433	CCL2	-3.14	down	chr17	6347	chemokine(C-C motif) ligand 2
7922337	TNFSF18	-3.29	down	chr1	8995	Tumor necrosis factor (ligand) superfamily, member 18
8021081	SLC14A1	-3.33	down	chr18	6563	solutecarrierfamily14(urea transporter), member 1(Kidd blood group)
7945680	H19	-3.62	down	chr11	283120///100033819///6206	H19,imprintedmaternallyexpressedtranscript(non-proteincoding)|microRNA675|ribosomalproteinS12
Endo IV_E2 vs Veh
7965873	IGF1	12.71	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
7951165	PGR	5.19	up	chr11	5241	Progesterone receptor
8117020	MYLIP	4.13	up	chr6	29116	Myosin regulatory light chain interacting protein
8111490	PRLR	3.44	up	chr5	5618	Prolactin receptor
8102950	INPP4B	2.88	up	chr4	8821	inositolpolyphosphate-4-phosphatase, typeII,105kDa
8040292	GREB1	2.77	up	chr2	9687	growthregulationbyestrogeninbreastcancer1
8107823	ADAMTS19	2.71	up	chr5	171019	ADAMmetallopeptidasewiththrombospondintype1motif,19
8088560	ADAMTS9	2.51	up	chr3	56999	ADAMmetallopeptidasewiththrombospondintype1motif,9
8145361	NEFM	2.48	up	chr8	4741	neurofilament, medium polypeptide
8106516	JMY	2.43	up	chr5	133746	Junction mediating and regulatory protein, p53 cofactor
7920123	S100A10	-1.94	down	chr1	6281	S100calciumbindingproteinA10
8131803	IL6	-2.03	down	chr7	3569	interleukin6
8135218	LRRC17	-2.11	down	chr7	10234	leucinerichrepeatcontaining17
7906919	RGS4	-2.15	down	chr1	5999	regulatorofG-proteinsignaling4
7921916	RGS5	-2.15	down	chr1	8490	regulatorofG-proteinsignaling5
7997139	CALB2	-2.20	down	chr16	794	calbindin2
8021081	SLC14A1	-2.23	down	chr18	6563	solutecarrierfamily14(urea transporter), member 1(Kidd blood group)
7922343	TNFSF4	-2.51	down	chr1	7292	Tumor necrosis factor (ligand) superfamily, member4
7922337	TNFSF18	-2.58	down	chr1	8995	Tumor necrosis factor(ligand) superfamily, member 18
8006433	CCL2	-3.30	down	chr17	6347	chemokine(C-C motif) ligand 2
Endo IV_P4 vs Veh (All DE loci)
7907271	FMO2	1.99	up	chr1	2327	flavincontainingmonooxygenase2(non-functional)
7908459	CFH	1.93	up	chr1	3075	Complement factor H
7965873	IGF1	1.85	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8043995	IL1R1	1.82	up	chr2	3554	interleukin1receptor, type I
8157524	TLR4	1.61	up	chr9	7099	toll-likereceptor4
7961514	MGP	1.59	up	chr12	4256	Matrix Gla protein
8089145	ABI3BP	1.58	up	chr3	25890	ABI family, member3 (NESH)binding protein
7933204	C10orf10	1.58	up	chr10	11067	chromosome10openreadingframe10
8101659	SPARCL1	1.51	up	chr4	8404	SPARC-like1(hevin)
8111490	PRLR	1.50	up	chr5	5618	Prolactin receptor
8111941	HMGCS1	-1.49	down	chr5	3157	3-hydroxy-3-methylglutaryl-CoAsynthase1(soluble)
Endo IV_E2+P4 vs Veh
8101659	SPARCL1	29.40	up	chr4	8404	SPARC-like1(hevin)
8089145	ABI3BP	17.62	up	chr3	25890	ABI family, member3 (NESH)binding protein
7965873	IGF1	17.30	up	chr12	3479	insulin-likegrowthfactor1(somatomedin C)
8122099	ENPP1	10.98	up	chr6	5167	Ectonucleotide pyrophosphatase/phosphodiesterase1
8155864	RORB	10.33	up	chr9	6096	RAR-related orphan receptor B
7977933	SLC7A8	8.57	up	chr14	23428	solutecarrierfamily7(amino acid transporter light chain, L system), member8
8043995	IL1R1	8.54	up	chr2	3554	interleukin1receptor, type I
7908459	CFH	6.96	up	chr1	3075	Complement factor H
8138231	THSD7A	6.18	up	chr7	221981	thrombospondin, type I, domaincontaining7A
8125919	FKBP5	6.12	up	chr6	2289	FK506bindingprotein5
7997139	CALB2	-3.00	down	chr16	794	calbindin2
8150428	SFRP1	-3.12	down	chr8	6422	secretedfrizzled-relatedprotein1
7922337	TNFSF18	-3.26	down	chr1	8995	Tumor necrosis factor (ligand) superfamily, member 18
7917561	GBP4	-3.27	down	chr1	115361	guanylatebindingprotein4
8055323	NCKAP5	-3.36	down	chr2	344148	NCK-associatedprotein5
8138289	ETV1	-3.59	down	chr7	2115	etsvariant1
8021081	SLC14A1	-3.81	down	chr18	6563	solutecarrierfamily14 (urea transporter), member1(Kidd blood group)
7922343	TNFSF4	-4.00	down	chr1	7292	Tumor necrosis factor (ligand) superfamily, member4
7906919	RGS4	-4.21	down	chr1	5999	Regulator of G-proteinsignaling4
8006433	CCL2	-6.47	down	chr17	6347	Chemokine (C-C motif)ligand2

Table 4

Differentially expressed loci common in each hormone treatment across normal (NUP), stage I (Endo I) and stage IV (Endo IV).

E2-induced upregulated loci common in NUP_Endo I_ Endo IV	E2-induced downregulated loci common in NUP_Endo I_ Endo IV	P4-induced upregulated loci common in NUP_Endo I_ Endo IV	E2+P4-induced upregulated loci common in NUP_Endo I_Endo IV		E2+P4-induced downregulated loci common in NUP_Endo I_ Endo IV
IGF1*^	NUAK1	ABI3BP	SPARCL1	SLC29A1	LIMCH1
GREB1^	NCAM1	SPARCL1	ABI3BP	ZBTB16	ROBO2
PGR^	SERPINE1	IGF1*	RORB	RAB31^	KRT18
LCP1^	TNFSF18	IL1R1*	GREB1^	KLF6	PIK3R3
ADAMTS9	SLC14A1U***	C10orf10	SLC7A8	PEMT	PTCHD4
ADAMTS19	FLT1	TLR4	IGF1^	RGS9	RASSF2
IL17RB	TNFSF4		THSD7A	PCBP3	GALNT5
SNCA^	MEST***		LCP1^	ITPR1	FRY***
ISOC1	LMO7***		FBXO32	AHNAK2	RNA5SP104
MYLIP	EDNRA***		MAOB	MAP3K4	FAM46C
TSKU^	RGS5***		FKBP5	SOD2	TYMS
NRG1^	FRY***		ENPP1	CNTN3	LYPD1
SEMA6A	SERTAD4***		CRISPLD2	ADAMTS1	SLC14A1***
GUCY1A2	CCL2***		ULK4	SPTSSA	ARHGAP18
PRLR^	RGS4***		IL1R1^	SORBS1	PLK2
PRICKLE2^			IMPA2	TSC22D3	ENC1
JMY			MUM1L1	LAMA3	SULF1
MIR503			GPX3	ITGBL1	MYO1B
KLF4			STMN2	SYTL4	SYTL5
IL1R1*^			CRYAB	PMP22	ATP8B1
AFAP1L2			SPSB1	SH3PXD2B	AMIGO2
ASPN^			TMEM37	MFGE8	UGCG
INPP4B			CPM	PIK3R1	KRTAP1-5
PXK			GALNT15	NID1	PLEKHG1
TMEM120B			MYOCD	EVA1C	GPR39
MIR503HG			ST6GALNAC2	APOD	LOXL4**
CPXM1			THBD	ABCC9	ATP2B1
LIN7A			FAM134B	CD151	TMEM130
NEFM			LPAR1	EPS8	MGAT5
BMP2			GPRC5B	SEPP1	GLT8D2
SLC35F6^			PRLR^	LAMA2	CCDC14
FAM102A			CERS6	FOXO1	MEST***
CPZ			MGST1	ARHGAP20	PPP1R3C
OSBPL3			ABHD5	PPAP2B	RARRES2
SNCAIP			C10orf10	NRG1^	RASGRF2
LOXL4**			CD68	TXNIP	PGRMC1
RNU2-6P			CXCR4	RAP2A	ACKR4
ACER3			GLULP4	FBLN2	LMO7***
RAB31^			EFEMP1	ETS2	GBP2
PAPSS2			SESTD1	TMOD1^	ARHGAP29
ASCC3			NLGN4X	HAND2-AS1	SERTAD4***
TMOD1^			PRUNE2	PTGER2	CSRNP3
C16orf45			OLFML2B	LHFP	GABRE
CDH4			SLC40A1	PPP1R14A	FJX1
			HAND2	LAMB1	RGS5***
			ADAMTS2	TLR4	FHOD3
			TSKU^	LMCD1	CD200
			GLUL	RAP1B	TNFRSF19
			ABLIM3	CFH	DUSP6
			CORIN	YBX3	PTN
			ACSL1	PGR^	TGFBI
			IRS2	RPS6KA2	KRT19
			MEDAG	SLC35F6^	F2RL2
			HOMER1	PPP1R3B	DACH1
			ADRA2C	NPC1	NCAM2
			DPT	PLIN2	MXRA5
			ASPN^	PRICKLE2^	RGS4***
			MOB3B	HSPA2	EDNRA***
			NFIL3	DPP4	NCKAP5
			INSR	ITGB1BP1	MOXD1
			ARRDC4	PQLC3	GBP4
			DKK1	APCDD1	H19
			SNCA^	ENDOD1	ETV1
			SRPX	ANTXR2	IL6
			UST	ANG	CCL2***
			ELMO1	GRIA1	SFRP1

* in both common E2 and common P4 up-regulated genes

^ in both E2 and E2+P4 up-regulated; BOLD, in E2, P4 and E2+P4 common loci

** Up regulated in E2 and down regulated in E2+P4

*** Down-regulated in E2 and E+P. NOTE: There are no P4-induced downregulated loci in common in NUP, Endo I, Endo IV.

Table 5

Functional enrichment cluster analysis for loci with strong association of differential methylation and gene expression changes in each hormone treatment in normal (NUP), stage I (Endo I) and stage IV (Endo IV) disease.

Enrichment score>2 and P<0.05	E2 vs. Veh	P4 vs. Veh	E2+P4 vs. Veh
NUP, pos/neg rho (gain and loss)	Cluster 1: Spectrins, Spectrin/alpha-actinin, Actinin-type, actin-binding, Calponin homology (CH) domain, Cell division and chromosome partitioning, DNA replication, recombination, and repair. Cluster 2: ATP-binding, Protein kinase-like domain, protein kinase activity, Nucleotide-binding, protein phosphorylation, intracellular signal transduction, Serine/threonine-protein kinase. Cluster 3: Pleckstrin homology-like domain, PH. Cluster 5: Metal-binding, Zinc, Zinc-finger, zinc ion binding. Cluster 6: Cell adhesion, Cadherins. Cluster 7: Biological rhythms, regulation of circadian rhythm, rhythmic process. Cluster 8: postsynaptic membrane, Cell junction, Synapse. Cluster 9: GTPase activator activity, Rho GTPase activation protein, regulation of small GTPase mediated signal transduction. Cluster 10: Fibronectins. Cluster 11: zinc finger region: ZZ-type.	Cluster 1: Dynein heavy chain/domain, microtubule motor activity. cluster 2; Zinc, zinc ion binding, Zinc-finger, Metal-binding.	No functional cluster enrichment with significant p-value and enrichment score>2.
Endo I, pos/neg rho, (gain and loss)	Cluster1: cell adhesion, cadherins, calcium ion channel. Cluster 2: a number of Sushi domains. Cluster 3: SAM. Cluster 4: ATP-binding, nucleotide phosphate-binding region: ATP, Nucleotide-binding, Kinase.	cluster 1: homophilic cell adhesion via plasma membrane adhesion molecules, Cadherins, Cell adhesion, calcium ion binding, CA. Cluster 2: Synapse, Cell junction, Postsynaptic cell membrane. Cluster 3: regulation of small GTPase mediated signal transduction, positive regulation of GTPase activity, Rho GTPase activation protein, signal transduction. Cluster 4: Sterile alpha motif domain (SAM). Cluster 5: regulation of small GTPase mediated signal transduction, Dbl homology (DH) domain, Rho guanyl-nucleotide exchange factor activity, Pleckstrin homology-like domain, positive regulation of apoptotic process, Src homology-3 domain, intracellular signal transduction.	Cluster 1: PDZ/PDZ domain. Cluster2: EGF, EGF-like, EGF-like calcium-binding, Laminin, Insulin-like growth factor binding protein, TB domain, Extracellular matrix, Secreted. Cluster 3: Calponin homology (CH) domain. Cluster 4: Cell junction, Synapse, Postsynaptic cell membrane, neuron projection.
Endo IV, pos/neg rho (gain and loss)	Cluster 1: plasma membrane, topological domain: Extracellular, Glycoprotein, Signal, Disulfide bond, Transmembrane. Cluster 2: Neurexin/syndecan/glycophorin C, cell adhesion.	No enrichment in functional clusters (total probes = 106) with our stat cut off, but some enrichment without ES>2 include: cluster 1 at 1.69: Fibronectins, Laminin, cell adhesion, extracellular matrix, secreted. Cluster 2: 1.13: cell adhesion, glycoprotein, cell membrane, signal.	cluster 1: repeat I, II, III, IV (this includes these: calcium voltage-gated channel subunit alpha1 C(CACNA1C), calcium voltage-gated channel subunit alpha1 A(CACNA1A), integrin subunit beta like 1(ITGBL1). (cluster 2: 1.8: Phosphotyrosine interaction domain. Cluster 3 (1.56):Pleckstrin homology-like domain. Cluster 4: (1.4): synapse, cell junction).

Table 6

Genes commonly regulated in eSF treated with E2 and E2+P4 and in whole endometrial tissue in PE and MSE or in FACS sorted eSF.

eSF E2 treated in common with WHOLE Tissue PE in Disease vs. Control
Gene Symbol	In Vitro array.Probe.Set. ID	FC ([Endo IV-E] vs [NUP-E])	FC ([Endo I-E] vs [NUP-E])	FC ([PE.Endo] vs [PE.Control])	Tissue array. Probe.Set.ID	Gene Title
SPARCL1	8101659	-2.86	-4.51	-2.77	200795_at	SPARC-like 1 (hevin)
CXCL12	7933194	-2.79	-2.02	-2.20	203666_at	chemokine (C-X-C motif) ligand 12
GREB1	8040292	-2.77	-3.49	-1.66	205862_at	growth regulation by estrogen in breast cancer 1
IL1R1	8043995	-2.53	-1.76	-2.30	202948_at	interleukin 1 receptor, type I
COLEC12	8021946	-2.23	-2.02	-2.33	221019_s_at	collectin sub-family member 12
LCP1	7971461	-2.12	-2.66	-1.73	208885_at	lymphocyte cytosolic protein 1 (L-plastin)
FAM134B	8111136	-2.09	-1.82	-2.82	218532_s_at	family with sequence similarity 134, member B
IGF1	7965873	-2.07	-2.76	-2.68	209540_at	insulin-like growth factor 1 (somatomedin C)
ADAMTS9	8088560	-2.02	-2.06	-3.77	226814_at	ADAM metallopeptidase with thrombospondin type 1 motif, 9
SLC40A1	8057677	-2.01	-2.25	-1.98	223044_at	solute carrier family 40 (iron-regulated transporter), member 1
PDE7B	8122222	-1.97	-2.00	-2.44	230109_at	phosphodiesterase 7B
PLIN2	8160297	-1.94	-1.65	-2.18	209122_at	perilipin 2
IL17RB	8080562	-1.93	-2.36	-2.27	224156_x_at	interleukin 17 receptor B
SPTLC2	7980438	-1.91	-1.72	-5.38	203128_at	serine palmitoyl transferase, long chain base subunit 2
CD109	8120719	-1.88	-1.49	-4.35	226545_at	CD109 molecule
TMX4	8064939	-1.85	-1.61	-2.67	201580_s_at	thioredoxin-related transmembrane protein 4
CCBE1	8023575	-1.85	-1.59	-1.83	243805_at	collagen and calcium binding EGF domains 1
FAM46A	8127778	-1.84	-1.68	-1.90	224973_at	family with sequence similarity 46, member A
EPHA5	8100578	-1.83	-1.86	-2.42	237939_at	EPH receptor A5
SPON1	7938608	-1.82	-1.84	-2.13	213994_s_at	spondin 1, extracellular matrix protein
ABCG2	8101675	-1.81	-1.94	-3.01	209735_at	ATP-binding cassette, sub-family G (WHITE), member 2
RSPO3	8121916	-1.75	-1.92	-2.21	228186_s_at	R-spondin 3 homolog (Xenopus laevis)
CLIC2	8176234	-1.74	-1.54	-2.68	213415_at	chloride intracellular channel 2
SLC18A2	7930837	-1.68	-1.78	-2.58	205857_at	solute carrier family 18 (vesicular monoamine), member 2
PRICKLE2	8088550	-1.67	-1.48	-2.83	225968_at	prickle homolog 2 (Drosophila)
KLHL13	8174654	-1.61	-1.51	-2.31	227875_at	kelch-like 13 (Drosophila)
PLA2R1	8056151	-1.59	-1.51	-1.51	207415_at	phospholipase A2 receptor 1, 180kDa
PLD1	8092134	-1.57	-1.47	-1.59	226636_at	phospholipase D1, phosphatidylcholine-specific
WT1	7947363	-1.53	-1.82	-2.49	216953_s_at	Wilms tumor 1
C17orf58	8017831	-1.48	-1.45	-1.87	226901_at	chromosome 17 open reading frame 58
TWISTNB	8138454	-1.46	-1.65	-2.81	226784_at	TWIST neighbor
LONRF2	8054281	-1.46	-2.44	-3.55	225996_at	LON peptidase N-terminal domain and ring finger 2
SERPINE1	8135069	1.45	1.65	3.39	1568765_at	serpin peptidase inhibitor, clade E (nexin, plasminogen activator inhibitor type 1), member 1
IER3	8124848	1.48		5.13	201631_s_at	immediate early response 3
EIF4A1	8004506	2.01	1.60	1.89	214805_at	Eukaryotic translation initiation factor 4A1
TBRG4	8139482	2.12	1.73	2.93	220789_s_at	transforming growth factor beta regulator 4
IGFBP5	8058857	3.28	3.24	2.45	203425_s_at	insulin-like growth factor binding protein 5
IER3	8178435		1.47	5.13	201631_s_at	immediate early response 3
eSF E2+P4 treated in common with WHOLE Tissue MSE in Disease vs Control
	In Vitro array.Probe.Set. ID	FC ([Endo IV-EP] vs [NUP-EP])	FC ([Endo I-EP] vs [NUP-EP])	FC ([MSE.Endo] vs [MSE.Control])	Tissue array. Probe.Set.ID	Gene Title
LPL	8144917	-4.54	-5.48	-1.62	203548_s_at	lipoprotein lipase
LGR4	7947199	-3.35	-3.41	-3.15	218326_s_at	leucine-rich repeat-containing G protein-coupled receptor 4
GREB1	8040292	-3.19	-4.97	-2.10	205862_at	growth regulation by estrogen in breast cancer 1
F3	7917875	-2.92	-2.11	-2.50	204363_at	coagulation factor III (thromboplastin, tissue factor)
CAB39L	7971590	-2.64	-2.33	-2.62	225915_at	calcium binding protein 39-like
SPARCL1	8101659	-2.28	-2.92	-1.91	200795_at	SPARC-like 1 (hevin)
GPR155	8056837	-2.24	-2.50	-1.63	239533_at	G protein-coupled receptor 155
FAM134B	8111136	-2.23	-2.76	-2.14	218532_s_at	family with sequence similarity 134, member B
RORB	8155864	-2.19	-2.73	-2.93	242385_at	RAR-related orphan receptor B
THSD7A	8138231	-2.10	-4.77	-2.28	214920_at	thrombospondin, type I, domain containing 7A
ADAMTS9	8088560	-2.09	-2.27	-2.74	226814_at	ADAM metallopeptidase with thrombospondin type 1 motif, 9
PRLR	8111490	-2.03	-2.67	-2.91	206346_at	prolactin receptor
LPAR1	8163257	-2.02	-2.14	-1.76	204037_at	lysophosphatidic acid receptor 1
ATP13A3	8092849	-2.02	-2.22	-1.58	212297_at	ATPase type 13A3
ADAMTS5	8069689	-1.99	-2.11	-3.14	229357_at	ADAM metallopeptidase with thrombospondin type 1 motif, 5
SAT1	8166469	-1.95	-2.20	-2.66	213988_s_at	spermidine/spermine N1-acetyltransferase 1
APOO	8171823	-1.94	-1.78	-2.29	221620_s_at	apolipoprotein O
CNKSR2	8166355	-1.91	-1.97	-1.79	229116_at	connector enhancer of kinase suppressor of Ras 2
ABHD5	8079153	-1.74	-1.91	-2.59	213935_at	abhydrolase domain containing 5
CTSS	7919800	-1.67	-1.56	-1.79	202902_s_at	cathepsin S
ANGPT2	8149071	-1.66	-1.65	-1.86	205572_at	angiopoietin 2
CADM1	7951807	-1.64	-1.87	-3.47	209032_s_at	cell adhesion molecule 1
ZCCHC6	8162147	-1.63	-1.87	-2.25	220933_s_at	zinc finger, CCHC domain containing 6
SESTD1	8057394	-1.63	-1.95	-2.31	226763_at	SEC14 and spectrin domains 1
KAL1	8171248	-1.62	-1.91	-1.97	205206_at	Kallmann syndrome 1 sequence
ETS2	8068593	-1.60	-1.46	-1.56	201328_at	v-ets erythroblastosis virus E26 oncogene homolog 2 (avian)
TSPAN12	8142524	-1.54	-1.52	-2.14	219274_at	tetraspanin 12
MDM1	7964810	-1.51	-1.58	-3.31	213761_at	Mdm1 nuclear protein homolog (mouse)
TSPAN13	8131600	-1.51	-1.52	-1.69	217979_at	tetraspanin 13
TMEM133	7943369	-1.48	-1.52	-1.77	223595_at	transmembrane protein 133
ARRDC4	7986350	-1.47	-1.90	-2.41	225283_at	arrestin domain containing 4
FERMT2	7979204	-1.45	-1.89	-2.29	209209_s_at	fermitin family member 2
TBRG4	8139482	1.99	1.57	2.35	220789_s_at	transforming growth factor beta regulator 4
eSF signature in Disease vs Control in E2, or E2+P4 treated eSF common with FACS sorted eSF in Disease vs Control
Compared to E2 treatment
	In Vitro array.Probe.Set. ID	FC ([Endo IV-E] vs [NUP-E])	FC ([Endo I-E] vs [NUP-E])	FC eSF FACS Endo vs eSF FACS control	Gene Title
TNXB|TNXA	8179935	-2.53	-1.83	-2.32	tenascin XB | tenascin XA pseudogene
LCP1	7971461	-2.12	-2.66	-2.12	lymphocyte cytosolic protein 1 (L-plastin)
FAM134B	8111136	-2.09	-1.82	-1.64	family with sequence similarity 134, member B
IL17RB	8080562	-1.93	-2.36	-2.06	interleukin 17 receptor B
TGFBR3	7917649	-1.80	-1.57	-1.97	transforming growth factor, beta receptor III
ABI3BP	8089145	-1.76	-2.10	-1.91	ABI family, member 3 (NESH) binding protein
KLHL13	8174654	-1.61	-1.51	-1.55	kelch-like 13 (Drosophila)
SNORD46	7901048	1.90	1.74	1.76	small nucleolar RNA, C/D box 46
SERTAD4	7909503	2.28	2.04	1.54	SERTA domain containing 4
RGS4	7906919	2.36	2.89	1.90	regulator of G-protein signaling 4
SULF1	8146863	3.79	2.31	1.61	sulfatase 1
Compared to E2+P4 Treatment
	In Vitro array.Probe.Set. ID	FC ([Endo IV-EP] vs [NUP-EP])	FC ([Endo I-EP] vs [NUP-EP])	FC eSF FACS Endo vs eSF FACS control	Gene Title
ABI3BP	8089145	-2.36	-6.50	-1.91	ABI family, member 3 (NESH) binding protein
ADAMTS5	8069689	-1.99	-2.11	-1.87	ADAM metallopeptidase with thrombospondin type 1 motif, 5
CRISPLD2	7997642	-1.58	-1.88	-1.55	cysteine-rich secretory protein LCCL domain containing 2
FAM134B	8111136	-2.23	-2.76	-1.64	family with sequence similarity 134, member B
KIAA0040	7922474	-2.46	-2.33	-1.50
LCP1	7971461	-2.12	-2.90	-2.12	lymphocyte cytosolic protein 1 (L-plastin)
MAOB	8172204	-2.24	-3.23	-1.54	monoamine oxidase B
PARM1	8095751	-1.90	-1.56	-2.30	prostate androgen-regulated mucin-like protein 1
SERTAD4	7909503	1.67	1.71	1.54	SERTA domain containing 4
SNORD46	7901048	2.16	1.60	1.76	small nucleolar RNA, C/D box 46
SULF1	8146863	2.12	1.67	1.61	sulfatase 1
TMEM37	8044813	-1.61	-2.23	-1.54	transmembrane protein 37

Table 7

Differentially methylated (DM) genes with PGR binding sites affected by E2+P4 treatment in normal eSF (NUP), stage I eSF (Endo I) and stage IV eSF (Endo IV) and the associated pathways by GO analysis in DAVID.

DM genes with PGR binding sites in E2+P4 treated NUP	DM genes with PGR binding sites in E2+P4 treated Endo I	DM genes with PGR binding sites in E2+P4 treated Endo IV
AFF1, ASH1L, C3orf21, CAP2, ELK3, CLSTN2, LOC442459, MYH13, MYLK4, NSMCE2, NUBPL, PHACTR1, RBM46, SAV1, SOX5, SPAG16, TBC1D1, TMEM232, TOP2B, WDR27, WWC3, NLK***, TTL***ATP8A2, C1QTNF9B, C8orf44, CACNA2D3, CDC73, CSRNP3, FMN1, GPM6B, HSD17B2, INPP5A, JPH1, KIAA1033, MCF2, MGAT4C, NEK10, ODZ2, OXR1, PIK3C3, PLXNC1, PTPRQ, RAB3C, SEMA5A, SGIP1, SNX25, SRGAP1, TBC1D4, TNFRSF9, VENTXP1, VTA1, WBP4, YEATS4, RBM20*, SRBD1*	ACTG2, ANKRD11, AP4E1, ARHGAP12, ARHGAP15, ATF6, BCAS1, C13orf26, C15orf61, CDH6, CLASP2, COL1A2, CPEB4,CROCC, CSTF3, DIRC3, DOPEY1, DSEL, EBF1, EYS, FARP1, FBN1, FBN3, FBXL2, FIGN, FLJ43860, FRAS1, FSIP1, GLI3, GLIS1, GNG7, GPR143, GPR176, GRID2, HERC2, HHLA2, HTRA3, ICA1L, INTS3, IQGAP2, IQGAP3, ITGA8, ITGB5, KCTD9, KDM4B, KIF26B, LOC285692, LRBA, LRRC8C, LTBP2, MAGI2, MED13L, METAP2, MYH8, NFIB, NHS, NOTCH4, NTM, PACSIN3, PCNX, PCNXL2, PDE3B, PLCG2, PPM1H, PRIM2, PTBP2, PTPN4, PTPRG, RGS7, RHOQ, ROBO1, RPS6KA2, RTKN2, SHROOM3, SLC44A1, SLIT3, SNTB1, SYNPR, TBX4, TGFBR3, THSD4, THSD7B, TNC, TRIO, UTRN,VWA3B, ZDHHC6, RBM20*, SRBD1*, NLK***, TTL***, CNTNAP5**, ITGBL1**, ZNF275**	AFF1, ASH1L, C3orf21, CAP2, ELK3, CLSTN2, LOC442459, MYH13, MYLK4, NSMCE2, NUBPL, PHACTR1, RBM46, SAV1, SOX5, SPAG16, TBC1D1, TMEM232, TOP2B, WDR27, WWC3, NLK***, TTL***ACER3, ARL15, ATP8B4, C10orf11, C16orf62, C17orf104, C17orf67, COL14A1, CREB5, DCBLD1, DENND2A, DLC1, DOCK4, EFR3A, FAM155A, FAM188A, FAM73A, FKBP5, FYCO1, GAL3ST2, GALNTL6, GULP1, HPS3, IGF1R, IL31RA, KIAA0564, KLHL29, LRRC16A, MOXD1, MYO1E, NMBR, NXN, PDE4DIP, PPFIA2, PSMD14, RECK, RGL1, SCHIP1, SIK3, SLC2A13, SNX29, SPSB1, ST7, TEAD1, TNFRSF1B, UST, WWTR1, YWHAQ, ITGBL1**, CNTNAP5**, ZNF275**
Pathways of DM genes with PGR binding sites in E2+P4 treated NUP (P<0.05)	Pathways of DM genes with PGR binding sites in E2+P4 treated Endo I (P<0.05)	Pathways of DM genes with PGR binding sites in E2+P4 treated Endo IV (P<0.05)
Cellular component morphogenesis; cell morphogenesis; cell part morphogenesis;single-organism organelle organization; cytoskeleton organization; endosome organization; neuron development; cell morphogenesis involved in differentiation; cell projection morphogenesis; negative regulation of developmental growth; neuron projection development; regulation of axon extension; regulation of extent of cell growth	Tissue morphogenesis; cellular response to transforming growth factor beta stimulus; organ morphogenesis; morphogenesis of an epithelium; metanephros development; extracellular matrix organization; animal organ development; system development; cell morphogenesis involved in differentiation; response to growth factor; epithelium development; cell development; urogenital system development; neuron projection development; cell morphogenesis; neuron development; anatomical structure formation involved in morphogenesis; cellular component morphogenesis; cellular component morphogenesis; kidney development; cell surface receptor signaling pathway; renal system development; cellular response to organic substance; regulation of cellular response to transforming growth factor beta stimulus; negative regulation of transmembrane receptor protein serine/threonine kinase signaling pathway; gland development; odontogenesis; embryonic morphogenesis; regulation of transmembrane receptor protein serine/threonine kinase signaling pathway; regulation of cell projection organization; neuron projection guidance; signal transduction; tube morphogenesis; mammary gland development; regulation of signal transduction; telencephalon cell migration; circulatory system development; negative regulation of cellular component movement; forebrain cell migration; negative regulation of cellular response to transforming growth factor beta stimulus; forebrain generation of neurons; regulation of neurogenesis; brain development; tangential migration from the subventricular zone to the olfactory bulb; regulation of cell development	Negative regulation of Wnt signaling pathway; intracellular signal transduction; single-organism organelle organization

* Common in NUP and Endo I

** Common in Endo I and Endo IV

*** Common in NUP, Endo I and Endo IV; , common in NUP and Endo IV