Abstract
Free full text
Mapping the temporal and spatial dynamics of the human endometrium in vivo and in vitro
Abstract
The endometrium, the mucosal lining of the uterus, undergoes dynamic changes throughout the menstrual cycle in response to ovarian hormones. We have generated dense single-cell and spatial reference maps of the human uterus and three-dimensional endometrial organoid cultures. We dissect the signaling pathways that determine cell fate of the epithelial lineages in the lumenal and glandular microenvironments. Our benchmark of the endometrial organoids reveals the pathways and cell states regulating differentiation of the secretory and ciliated lineages both in vivo and in vitro. In vitro downregulation of WNT or NOTCH pathways increases the differentiation efficiency along the secretory and ciliated lineages, respectively. We utilize our cellular maps to deconvolute bulk data from endometrial cancers and endometriotic lesions, illuminating the cell types dominating in each of these disorders. These mechanistic insights provide a platform for future development of treatments for common conditions including endometriosis and endometrial carcinoma.
Main
The human endometrium is the site of implantation that provides nutritional support to the placenta throughout pregnancy. Unlike other mucosal tissues, it undergoes dynamic, cyclical changes of shedding, regeneration and differentiation throughout reproductive life coordinated by the hypothalamic–pituitary–ovarian axis. Endometrial dysfunction underpins many common disorders, including abnormal uterine bleeding, infertility, miscarriage, pre-eclampsia, endometriosis and endometrial carcinoma, that collectively affect many women across the world1–5. Throughout reproductive years, the functional upper layer of the endometrium, the stratum functionalis, is shed at menstruation. The subsequent tissue repair and proliferation are driven by rising levels of estrogen, secreted by the ovarian follicle during the first half of the menstrual cycle (proliferative phase). Following ovulation, progesterone, produced by the corpus luteum, induces the secretory phase during which the initial changes of decidualization occur. Menstruation and spontaneous decidualization are unique to higher simian primates6–8. Thus, dissecting the mechanisms that regulate cellular differentiation across the menstrual cycle in humans is crucial for understanding how normal endometrium is regulated.
Essential to endometrial function are the lumenal and glandular epithelia, composed of a mixture of ciliated and secretory cells. The lumenal epithelium is the site of embryo attachment covering the endometrial surface. Long tubular glands open into the lumenal epithelium and produce secretions that are rich in growth factors and lipids necessary for placental growth9. We and others have established a three-dimensional in vitro organoid culture model of human endometrial epithelium10,11. These organoids are generated from dissociated endometrial tissue and menstrual fluid samples12; they retain the morphology, function and gene signature of the tissue in vivo and respond functionally to ovarian hormones with differentiation into ciliated and secretory cells. They are therefore powerful platforms to investigate endometrial disorders and study mechanisms regulating endometrial differentiation in humans. However, a systematic, quantitative comparison of endometrial organoids at a single-cell level with epithelial cell states in vivo is lacking. This is required to confirm their suitability for exploring the cellular pathways and processes involved in normal and pathological endometrial function.
The explosion in spatial transcriptomics technologies13–16 provides a unique opportunity to resolve tissue architecture in conjunction with underlying cellular interactions. The spatial arrangement of cells is key to understanding a morphologically complex tissue such as the endometrium, where a cell’s function may differ depending on signals it receives from neighboring cells17. Many spatially resolved transcriptomics methods are not quite at single-cell resolution and rely on the computational integration of coupled single-cell (or single-nucleus) transcriptomes to achieve this level of detail18–20. These genomic technologies are the basis of the Human Cell Atlas initiative, which aims to map all cells in the human body21.
In this study, by using single-cell and spatial transcriptional profiling, we interrogate the cellular states and spatial localization of human endometrial cells during the proliferative and secretory phases of the menstrual cycle in women of reproductive age. We develop CellPhoneDB v.3.0 to measure intercellular communication taking into account spatial coordinates of cells and use this tool to define cell signaling in both lumenal and glandular epithelial microenvironments. We define a complementary role for WNT and NOTCH signaling in regulating differentiation toward the two main epithelial lineages (ciliated and secretory). We profile three-dimensional endometrial organoids at single-cell resolution to characterize their hormonal responses in vitro and design a computational toolkit to compare the results with those observed in vivo to benchmark this model system. Finally, by modulating WNT and NOTCH pathways in the organoid cultures, we develop lineage-specific endometrial epithelial cells and define the molecular events involved in their response to ovarian hormones.
Results
A single-cell map of the full-thickness human uterus
To generate a cellular map of the human endometrium that accounts for the temporal and spatial changes across the menstrual cycle, uterine samples were analyzed by single-cell transcriptomics—single-cell RNA sequencing (scRNA-seq) and single-nucleus RNA sequencing (snRNA-seq)—alongside spatial transcriptomics methods (10x Genomics Visium slides and high-resolution microscopy) (Fig. (Fig.1a1a and Extended Data Fig. Fig.1a).1a). Two different types of samples from women of reproductive age were integrated in our analysis: endometrial biopsies from live donors screened for potential endometrial disorders (n=3) and the whole endometrium with attached subjacent myometrium from the uteri of donors who died of nongynecological causes (n=6) (Supplementary Table 1). This latter approach allows sampling of the endometrial basal layer and myometrium, which are absent from endometrial biopsies.
We were able to generate a uterine map of 98,568 cells from 15 individuals by integrating our dataset with previous scRNA-seq data obtained from superficial endometrial biopsies22 (Fig. (Fig.1b,1b, Extended Data Fig. 1b–f and Supplementary Table 2). We identify 14 clusters that were assigned cell identity based on their expression of known markers (Extended Data Fig. Fig.1g).1g). These clusters can be grouped into five main cellular categories: (1) immune (lymphoid and myeloid); (2) epithelial (SOX9+, lumenal, glandular and ciliated); (3) endothelial (arterial and venous); (4) supporting—perivascular cells (PV STEAP4 and PV MYH11), smooth muscle cells and fibroblasts expressing C7 (fibroblasts C7); and (5) stromal—nondecidualized endometrial (eS) and decidualized endometrial (dS). SOX9+ epithelial cells and eS are characteristic of the regenerating proliferative phase (Fig. (Fig.1c).1c). Subclustering of immune cells resolves their heterogeneity, including identification of the three uterine natural killer cell subsets we have previously defined in the early pregnant uterus23 (Fig. (Fig.1d1d and Extended Data Fig. Fig.1h).1h). snRNA-seq data from four additional full-thickness uterine samples of proliferative endometrium confirm the populations found in the proliferative phase (Fig. 1e,f and Extended Data Fig. 2a–d). We additionally found lymphatic endothelial and mast cells enriched in these samples; these are likely to originate from the myometrium (Fig. (Fig.1e1e and Extended Data Fig. Fig.2d2d).
To systematically map the location of the cell types identified by scRNA-seq within the endometrium and myometrium, we used Visium Spatial Transcriptomics technology. We examined four full-thickness uterine samples in the proliferative and secretory phases from two individuals (Extended Data Fig. 3a–e). After integrating single-cell transcriptomics and Visium data using our recently developed cell2location algorithm18, cell states were mapped to the endometrium and/or the myometrium. We identify specific perivascular cells: PV MYH11 are characteristic of myometrium while PV STEAP4 are only present in the endometrium (Fig. (Fig.1g1g and Extended Data Fig. Fig.3f).3f). In addition, we find that fibroblasts C7 are enriched in the basal layer of the endometrium in both the proliferative and secretory phases (Fig. 1h,i and Extended Data Fig. Fig.3f3f).
Altogether, our analysis yields a comprehensive catalog of the major subsets of uterine cells together with their cellular position in endometrium and myometrium. We have made an open-source web server available at www.reproductivecellatlas.org.
Spatiotemporal characterization of proliferative epithelium
We next focused on the two main lineages of endometrial epithelial cells, secretory and ciliated, across the menstrual cycle and analyzed these subsets individually (Fig. (Fig.2a).2a). Epithelial cells are classified into four main groups based on their marker expression: (1) SOX9 populations, enriched in the proliferative phase and expressing genes characteristic of rising estrogen levels (MMP7, ESR1); (2) ciliated cells (PIFO, TPPP3); (3) lumenal cells (LGR5); and (4) glandular cells (SCGB2A2) (Fig. 2b–d and Extended Data Fig. 4a,b). A fraction of the glandular epithelial cells express molecules characteristic of uterine milk in the secretory stage (PAEP, CXCL8).
With our integrative scRNA-seq maps, we can now resolve three clusters within the SOX9 population: (1) SOX9+LGR5+ cells, expressing KRT17 and WNT7A; (2) SOX9+LGR5− cells, expressing IHH; and (3) proliferative SOX9+ cells, including both LGR5+ and LGR5− cycling cells (that is, cells in G2M/S phase). By integrating scRNA-seq and Visium data (Fig. (Fig.2e2e and Extended Data Fig. Fig.4c),4c), we define specific spatial coordinates for these SOX9 subsets: (1) noncycling SOX9+LGR5+ cells are enriched in the surface epithelium; (2) noncycling SOX9+LGR5− cells are located in the basal glands; and (3) cycling SOX9+ cells map to glands in the regenerating superficial layer. Single-molecule fluorescence in situ hybridization (smFISH) with RNAscope probes shows higher expression of the proliferative marker MKI67 in the superficial layer of the endometrium during the proliferative phase, validating Visium data (Extended Data Fig. Fig.5a).5a). We also validated the presence of the markers LGR5 and WNT7A, characteristic of the SOX9+LGR5+ cells present in surface epithelium during the proliferative phase (Fig. (Fig.2f2f and Extended Data Fig. 5b,c). SOX9 and LGR5 are expressed in stem cells/progenitors in several tissues including gut, kidney, skin and ovaries and may label similar populations in the human endometrium24–29.
Ciliated cells are present in both the proliferative and secretory phases, but, as expected, PAEP secretory cells are only present following ovulation (Fig. 2b,c and Extended Data Fig. Fig.4b).4b). This indicates that estrogen alone can induce ciliary differentiation, while secretory differentiation depends on the addition of progesterone. During the proliferative phase, in addition to the FOXJ1, PIFO, TP73 ciliated population, we define a distinct subset of preciliated cells that express MUC12, HES6 and cell cycle genes (CDC20B, CCNO) (Fig. (Fig.2d2d).
After ovulation, secretion of progesterone induces the differentiation of SOX9+ cells into specialized secretory cells found at the surface and in glands (Extended Data Fig. Fig.5d).5d). To date, it has not been possible to distinguish between differentiated lumenal and glandular subsets using specific transcriptomic signatures and markers. We can now identify these subsets in our scRNA-seq data, confirmed by integration with spatial transcriptomic data (Fig. (Fig.2g2g and Extended Data Fig. 4d,e) and revealing markers that are validated at the protein level (Fig. (Fig.2h).2h). There is enriched expression of both COX1 (encoded by PTGS1) and KRT5 in the lumenal epithelium and SCGB2A2 in the glandular epithelium.
SOX9+ epithelial cells in endometrial disorders
Disorders in endometrial function have a profound impact on women’s health and reproductive outcomes. There has been limited progress in the study of these disorders over the past decade, partly due to the challenges in analyzing this highly dynamic and complex tissue. To identify the cells involved in these disorders, we looked for transcriptomic signatures in bulk RNA data from endometrial cancers and endometriotic lesions.
We first deconvoluted the transcriptomes from bulk RNA sequencing data from serous and endometrioid endometrial adenocarcinomas available in The Cancer Genome Atlas (TCGA), using our single-cell atlas as a reference. Deconvolution reveals that the signature of the SOX9+ epithelial population is dominant in tumors (Fig. (Fig.3a).3a). Based on expression signals in epithelial cells, endometrial adenocarcinomas resolve into three patterns. The first, characterized by SOX9+LGR5+, occurs in 63% of serous and 24% of endometrioid adenocarcinomas. The second, SOX9+LGR5−, is found in 33% of endometrioid adenocarcinomas but is absent from serous adenocarcinomas. The third pattern, with expression of differentiated ciliated or glandular signals, is present in <10% of endometrial carcinomas. In line with these results, markers characteristic of the SOX9+LGR5+ subset (for example, WNT7A, MSLN, KRT17, PTGS1 and VTCN1) are all upregulated in tumors with a SOX9+LGR5+ transcriptomic signature compared to tumors with a SOX9+LGR5− signature (Extended Data Fig. Fig.6a6a).
We then correlated the clinical stages of endometrial adenocarcinomas with our cell signals. The more advanced stages of endometrial adenocarcinomas (stages III and IV) have a greater SOX9+LGR5+ signal (Wilcoxon test, stages I+II against stages III+IV, P value 7.54×10−6) (Fig. (Fig.3b3b and Supplementary Table 3). We also linked the SOX9+LGR5+ signals to the four molecular subtypes of endometrial cancer defined by TCGA30. The SOX9+LGR5+ signature is stronger in tumors characterized by high copy number (copy-number-high) alterations (Wilcoxon text, P value 6.72×10−13), typical of serous endometrial adenocarcinomas and linked with a worse prognosis (Supplementary Table 3). Indeed, all the copy-number-high tumors considered in our analysis (31 serous and 7 serous-like endometrioid adenocarcinomas) were classified as SOX9+LGR5+. In contrast, copy-number-low tumors have a lower SOX9+LGR5+ signal (Wilcoxon text, P value 5.91×10−6) (Supplementary Table 3).
We then explored the expression of specific epithelial markers in microarray expression data available from peritoneal biopsies from donors with endometriosis31 (Supplementary Table 3). As expected, endometriotic peritoneal lesions upregulate markers characteristic of proliferative endometrium (SOX9+ and preciliated markers) compared with normal peritoneum (Fig. (Fig.3c).3c). In particular, peritoneal lesions upregulate markers specific for the SOX9+LGR5+ subset (such as WNT7A and KRT17) with expression levels similar to those in proliferative endometrium (Fig. (Fig.3c).3c). In contrast, markers of secretory epithelial cells, PAEP and SCGB2A2, or ciliated cells, PIFO and TP73, are expressed at similar levels as in normal peritoneum (Extended Data Fig. Fig.6b6b).
In summary, our analysis suggests that dysfunctional epithelium is a major driver of endometrial disease. By defining the transcriptomes of our two SOX9 populations, we can show the specific cell signals that dominate in endometrial carcinomas and endometriosis.
Effect of microenvironments on epithelial identity
Having defined these distinct epithelial cell states and their potential role in pathology, we focused next on the transcription factors (TFs) that regulate epithelial differentiation throughout the menstrual cycle by analyzing expression of TFs and their consensus target genes32 (Fig. (Fig.4a4a and Supplementary Table 4). There is high activity of WNT targets (for example, FOXJ1) in the ciliated epithelium. In contrast, the glandular subsets show high expression for TFs induced by WNT inhibition (for example, CSRNP1 and FOXO1) and NOTCH activation (for example, HES1 and HEY1). This suggests different roles for NOTCH and WNT in shaping the identity and function of ciliated versus secretory cells.
To investigate the cell signals operating in the lumenal and glandular microenvironments that could influence differentiation into ciliated and secretory lineages, we used spatial transcriptomics and performed clustering on the 10x Genomics Visium spots assigned to epithelial subsets. We resolve five clusters corresponding to cells in the lumenal (one cluster), functional (two clusters) and basal (two clusters) layers (Fig. (Fig.4b).4b). In addition to cell-type-specific markers, signatures of WNT and NOTCH signaling pathways are present in distinct endometrial regions (Fig. (Fig.4c4c and Supplementary Table 5). Genes involved in the WNT pathway, FOXJ1 and LGR5, are highly expressed at the lumenal surface while NOTCH2 is enriched in glands in the functional layer. To validate expression of NOTCH2 and WNT7A in these compartments, we stained uterine tissue with smFISH probes for both genes alongside EPCAM using immunohistochemistry (IHC) (Fig. 4d,e). Lumenal and glandular epithelial cells were classified automatically based on EPCAM expression, and the distance of the signal from the lumen was then measured (Methods). Our results show that NOTCH2 expression increases in glands moving away from the lumen while WNT7A expression is higher in the lumenal epithelium compared with glands (Fig. (Fig.4e).4e). By contrast, the noncanonical WNT molecule WNT5A was mainly expressed in stromal cells surrounding the glands (Extended Data Fig. Fig.7a).7a). These findings suggest that canonical WNT is downregulated in the glandular microenvironment where noncanonical WNT pathways dominate.
To investigate how surrounding cells may shape signaling in the surface and glandular compartments, we developed CellPhoneDB v.3.0, an updated version of our cell–cell communication pipeline that takes into account spatial cellular colocalization when mapping ligand–receptor pairs33 (Methods and Fig. Fig.5a).5a). CellPhoneDB considers the multimeric composition of the majority of ligands and receptors, which is highly relevant for the complex regulation of WNT signaling (Fig. (Fig.5b).5b). We define three endometrial microenvironments centered on epithelial cells based on the cellular coordinates provided by cell2location: (1) lumenal—preciliated, ciliated and SOX9+LGR5+ epithelium (proliferative phase) and ciliated and lumenal (secretory phase); (2) functional—SOX9+ proliferative epithelium, immune and eS (proliferative phase) and immune, glandular and dS (secretory phase); and (3) basal—SOX9+LGR5− and fibroblasts C7. We ran CellPhoneDB on each of the microenvironments (Supplementary Table 6) and found significant epithelial–stromal interactions (log fold change>0.02 and false discovery rate (FDR)<0.005).
NOTCH interactions are mainly mediated by the epithelial compartment (Fig. (Fig.5c),5c), with the NOTCH ligand JAG1 more highly expressed in the lumen than in glands, as shown by smFISH (Extended Data Fig. 7b,c). As well as JAG1 being coexpressed with HEY1, sparse JAG1high epithelial cells are often adjacent to JAG1low cells expressing HEY1 (Extended Data Fig. Fig.7c).7c). A similar lateral NOTCH inhibition model has been described in the gut34.
WNT ligands are expressed by both epithelial and stromal cells (Fig. (Fig.5c);5c); the latter express WNT agonists that can potentially bind the cognate WNT receptors expressed by all epithelial subsets during the proliferative phase (Fig. (Fig.5c).5c). Focusing on genes that are differentially expressed following ovulation, glandular secretory subsets show a dramatic decline in WNT receptor expression, potentially limiting the activity of this pathway (Fig. (Fig.5c).5c). They also show a decrease in the WNT target AXIN2 when compared with their lumenal counterparts (Extended Data Fig. Fig.7d).7d). In addition, decidualized stromal cells express significantly higher levels of DKK1, a potent inhibitor of the WNT pathway, than their nondecidualized counterparts. Expression of DKK1 surrounding the glands of the secretory endothelium is also found in spatial transcriptomics (Fig. (Fig.5d).5d). Overall, these findings strongly suggest that WNT signaling is inhibited in the secretory cell lineages, meaning that NOTCH signaling will then dominate (Fig. (Fig.5e5e).
Response of endometrial organoids to ovarian hormones
To test our predictions on the potential roles of WNT and NOTCH signaling pathways on endometrial epithelium in vitro, we first profiled endometrial organoids at a single-cell level to benchmark this model system against our in vivo data. We derived organoids10 from three different donors and primed them as previously described with estrogen for 48h, followed by stimulation with progesterone, prolactin and cyclic AMP (cAMP) in the presence of estrogen for 4d (ref. 10) (Fig. (Fig.6a6a and Supplementary Tables 7 and 8).
To identify epithelial cells in the organoids, we looked for markers specific to the clusters. Before hormonal treatment, the majority of cells within the organoids are proliferative (expressing TOP2A, PCNA) and all express the estrogen receptor ESR1 (Fig. 6b–d and Extended Data Fig. 8a,b). Two additional populations emerge when the organoids are treated with estrogen: (1) an estrogen-induced population expressing the progesterone receptor (PGR), a target gene of estrogen, and (2) a preciliated population sharing markers with the equivalent cluster defined in vivo. Upon further stimulation with progesterone, markers of more advanced stages of differentiation emerge in both secretory and ciliated populations. Markers of secreted products (PAEP, DEFB1) and glands (SCGB2A2) are both seen. Ciliated cells express typical markers (FOXJ1, TP73) and closely match their in vivo counterparts. IHC of endometrial organoids confirmed the expression of glandular and ciliary markers after hormonal stimulation (Fig. (Fig.6e6e).
We next assessed how closely the hormonal responses of the organoids that generate glandular and ciliated cells relate to their in vivo counterparts. To focus on the hormonal responses to estrogen and progesterone, we performed a quantitative assessment of organoids exposed to estrogen and progesterone while maintaining them in their expansion medium (ExM)10. Next, we projected the epithelial in vivo reference data onto the hormone-treated in vitro epithelial subsets. Assignments were made based on logistic regression predictions. In the absence of progesterone, the SOX9+ populations are the best match for estrogen-induced cells while preciliated organoid cells align with their preciliated in vivo counterparts (Fig. (Fig.6f,6f, Extended Data Fig. Fig.8c8c and Supplementary Table 9). In response to hormones, a large fraction (over a quarter of secretory cells) correspond to glandular epithelium, and all ciliated cells match perfectly with their in vivo counterparts, indicating that organoids respond similarly to hormones as in vivo (Fig. (Fig.6f6f and Supplementary Table 9).
To determine whether the pathways driving hormonally induced differentiation toward the ciliated and secretory lineages are similar in vitro and in vivo, we looked at the differential expression and activity of lineage-specific TFs (Fig. (Fig.6g6g and Supplementary Table 10). We calculated differentially expressed or active TFs in the hormonal subsets in vivo and in vitro. WNT-activated TFs (FOXJ1) are present in the ciliated lineage while WNT-inhibitory TFs (CSRNP1, FOXO1) are in the secretory lineage. NOTCH-induced TFs, HEY1 and HES1, are activated in the secretory lineage. These results indicate that ovarian hormones activate similar pathways both in vivo and in vitro.
These results mean that it is possible to reconstruct pseudotime and recapitulate cell fate decisions of epithelial cells in response to hormones using endometrial organoids. We performed scRNA-seq on two clonal organoids from the same individual (Extended Data Fig. Fig.8d).8d). Both clones showed similar behavior and were integrated under the same manifold (Extended Data Fig. 8e,f). Annotation of clusters was performed based on known markers, and then clusters expressing the progesterone receptor were selected to reconstruct epithelial differentiation in response to hormones (Fig. (Fig.6h6h and Extended Data Fig. Fig.8g).8g). A subset of cells emerging from the estrogen-induced population differentiate into preciliated cells in response to estrogen and, following additional progesterone, into ciliated cells. The secretory lineage also emerges from the estrogen-induced population. Thus, there is a common progenitor for both lineages.
WNT and NOTCH inhibition mediates epithelial differentiation
To test the roles of the WNT and NOTCH pathways in ciliated and secretory differentiation, we cultured organoids in the presence of inhibitors of either NOTCH (DBZ or DAPT) or WNT (IWP-2 or XAV939). We used functional, histological and single-cell transcriptomic assays to assess outcomes (Fig. (Fig.7a).7a). Organoid viability is high under all conditions (Extended Data Fig. Fig.9a).9a). scRNA-seq analysis reveals a higher proportion of preciliated and ciliated cells and a lower proportion of secretory cells in the presence of DBZ (NOTCH inhibitor) (Fig. 7b–d, Extended Data Fig. 9b,c and Supplementary Table 11). This finding was validated by IHC and quantitative PCR with reverse transcription (RT–qPCR) using two NOTCH inhibitors (DBZ and DAPT) (Fig. (Fig.7e7e and Extended Data Fig. Fig.9d).9d). Ciliated cells are virtually absent when WNT is inhibited (presence of XAV939), highlighting the strong dependence on this pathway for ciliary commitment, while the proportion of secretory cells is increased under these conditions (Fig. 7b–d, Extended Data Fig. 9b,c and Supplementary Table 11). The drive toward the secretory lineage under WNT-inhibitory conditions (IWP-2 and XAV939) was further validated by RT–qPCR and IHC (Fig. (Fig.7e7e and Extended Data Fig. Fig.9d).9d). These results demonstrate the fine balance between NOTCH and WNT signaling to regulate commitment to endometrial secretory or ciliary lineages.
Hormonal stimulation in the presence of WNT and NOTCH inhibitors modified the secretory cell transcriptome (Fig. (Fig.7b).7b). To quantify the similarity of the secretory populations that emerged in the presence of inhibitors with their in vivo counterparts, we measured the expression levels of genes encoding secretory products. Expression levels for PAEP and DEFB1 are higher when WNT is inhibited and downregulated with NOTCH inhibition (Fig. (Fig.7f).7f). We validated the increase of one secretory product, PAEP, under WNT-inhibitory conditions by ELISA (Extended Data Fig. Fig.9e).9e). Although PAEP expression levels increase slightly in the presence of WNT inhibitors even in the absence of hormones, although this is not significant (Fig. (Fig.7f).7f). In addition, we built a logistic regression model training the data on the organoids stimulated with distinct hormonal and inhibitory conditions. Computational projection of the in vitro dataset onto the in vivo dataset shows similar overlaps of the secretory populations emerging from WNT-inhibitory conditions and controls with their in vivo counterparts (Fig. (Fig.7g7g and Supplementary Table 12). This shows that WNT inhibitors target a specific gene module relevant for glandular secretions.
We next dissected the regulatory programs in the secretory lineage by comparing expression and activity of TFs between populations emerging after treatment with and without hormones and WNT inhibition. The expression of NOTCH-regulated TFs (HEY1) is only upregulated in the presence of WNT inhibitors when hormones are present. This probably explains why WNT inhibitors are not sufficient to induce the secretory lineage on their own (Fig. 7h,i and Supplementary Table 13). In the presence of hormones, WNT inhibition represses TFs characteristic of the ciliated lineage (FOXJ1, TP73, RFX2) (Fig. 7h,i and Supplementary Table 13). Switching off these genes therefore drives secretory lineage differentiation.
To summarize, we show that by inhibiting NOTCH and WNT pathways it is possible to influence cell fate decisions between ciliary and secretory differentiation. This effect depends on the presence of estrogen and progesterone, and we can dissect the dialog between signaling pathways and hormonal stimulation in endometrial differentiation.
Discussion
Profiling the uterus in space and time is essential to define the cell states and signaling pathways of normal human endometrium. This spatiotemporal reference, using samples from healthy women, will improve understanding of the molecular and cellular aberrations occurring in common conditions including infertility, endometriosis and endometrial carcinoma. The uterine lining in women of reproductive age is a challenging tissue to study due to difficulty in accessing samples covering the dynamic changes occurring across all stages of the menstrual cycle. Here, we have used single-cell expression analysis, spatial transcriptomics and high-resolution quantitative multiplex imaging to generate profiles of uterine cell states throughout the cycle. We focused on epithelial populations as they are major players in endometrial function and pathology, as suggested by our comparison of gene signatures with bulk RNA data from endometrial diseases. We use and develop computational tools to integrate and analyze scRNA-seq and spatial data and investigate the molecular mechanisms driving epithelial differentiation in the glandular and lumenal microenvironments. We utilize our reference atlas to benchmark endometrial organoids and engineer lineage-specific organoids informed by signaling factors predicted by our in vitro–in vivo comparisons. Our work shows the potential for using human cell atlases as blueprints for tissue engineering experiments.
There are two advances from our spatiotemporal uterine cell reference atlas. First, we have profiled uterine cells from transplant donors, allowing us to include the endometrial basal layer and the myometrium, which are absent from endometrial biopsies. To improve the temporal resolution of the functional endometrial layer, we have combined our dataset with another recent single-cell atlas of endometrial biopsies22. By integrating these two datasets, we have revealed other cell states, including a population of fibroblasts (fibroblasts C7) restricted to the basal layer. Second, our strategy of spatial mapping with 10x Genomics Visium and quantitative multiplexed smFISH techniques allows us to determine three-dimensional cellular arrangements described in transcriptomic analysis of cell isolates. By mapping cells into tissues with our deconvolution method18, we allocated epithelial cells into the three main endometrial layers: lumenal, functional and basal. Our expanded CellPhoneDB v.3.0 analysis framework dissects the cell signaling mechanisms in the lumenal and glandular endometrial microenvironments, which revealed that NOTCH and WNT pathways control ciliated and secretory epithelial cell commitment.
Using our single-cell transcriptomics data as a reference to deconvolute bulk data, we show that endometrial carcinomas and endometriotic lesions have a less differentiated epithelial phenotype than normal endometrium. Endometrial adenocarcinomas have two main signatures, SOX9+LGR5+ and SOX9+LGR5−, indicating that these distinct transcriptomic signatures reflect differences in pathogenesis and disease progression. SOX9+LGR5+ is the dominant signature in serous endometrial and some endometrioid adenocarcinomas and is positively associated with the ‘Copy-Number high’ molecular subtype from TCGA, as well as the clinically more aggressive stage III and IV adenocarcinomas. Our study is consistent with previous analyses showing that serous adenocarcinomas share molecular features with a subset of endometrioid tumors30. We also show an enrichment of markers from SOX9+LGR5+ epithelial cells in endometriotic lesions, in line with findings that organoids derived from higher stages (III–IV) of endometriosis have higher expression of SOX9 (ref. 35). These organoid models will allow investigation of the potential roles of these SOX9 populations in driving endometrial cancer, endometriosis and other disorders35.
Endometrial organoids are a powerful model to study normal endometrial epithelium10,11,36. We have systematically benchmarked the cellular composition of organoids relative to our in vivo reference map. Machine learning approaches, such as logistic regression scoring of expression profiles as well as correlation analysis, have been previously used to compare in vitro datasets with their in vivo counterparts37,38. Using a logistic regression approach, we demonstrate that endometrial organoids recapitulate the in vivo response to hormones, with the ciliary lineage becoming fully differentiated while the secretory lineage maintains progenitors in addition to differentiated secretory cells. This provides compelling evidence for the validity of our model. The organoids continue to expand during exposure to ovarian hormones because they are cultured in ExM that promotes proliferation. We also compared TFs operating in vivo and in vitro and showed that similar programs are induced. Informed by these findings, further optimization of the culture conditions to achieve more complete secretory differentiation is underway. Our computational kit for in vivo and in vitro comparisons will be of general utility for tissue engineering experiments using Human Cell Atlas data as a blueprint.
Our extensive validation assessing the hormonal responses of the endometrial organoids means that we could use them to test the effects of NOTCH and WNT signaling on epithelial cell fate10. Inhibition of WNT signaling, which mimics the low-WNT microenvironment in differentiated glands, inhibits ciliary commitment and induces secretory cells. In the presence of hormones, these cells produce more secretions, probably through stronger silencing of ciliary genes. WNT inhibition alone does not result in secretory differentiation, as the NOTCH pathway is not induced without hormonal stimulation. These results reinforce previous findings suggesting tight coordination between these signaling pathways and ovarian hormones39,40. NOTCH inhibitors were found to promote the generation of ciliated cells in fallopian tube41 and endometrial organoids42,43. We also show that NOTCH-inhibited secretory cells express lower amounts of uterine milk proteins. Using single-cell mapping, we pinpoint the effect of NOTCH to early ciliary differentiation, as suggested by the strong effect NOTCH inhibition has on preciliated cells. Thus, we demonstrate opposing roles of WNT and NOTCH in shaping distinct endometrial epithelial lineages. In vivo, this is regulated by the boundaries set by the localization of distinct cellular populations in the lumenal versus glandular microenvironments.
Our integrative map of cellular profiles of the normal endometrium will serve as an essential reference for the study of many neglected endometrial disorders. Organoids, which can be biobanked, have been established from samples of endometriosis and endometrial adenocarcinomas that resemble the original tumors10,11,35,44. Our study shows that the combination of genomics, imaging and organoids can create a robust platform for studying endometrial physiology. This will have a wide-ranging impact on women’s health and reproductive medicine.
Methods
Uterine tissue retrieval
Full-thickness uterine wall samples were obtained from deceased transplant organ donors (A13, A30) after ethical approval (reference 15/EE/0152, East of England–Cambridge South Research Ethics Committee) and informed consent from the donor families. Uterus was removed within 1h of circulatory arrest.
Full-thickness uterine wall samples were collected from four women during autopsy (Trv2, Trv3, Trv4 and Trv5). All four individuals died of noncancer-related causes, either traumatic injuries (Trv2, Trv3 and Trv4) or brain edema (Trv5). Samples were collected within less than 10h of death (postmortem interval was 5, 2, 4 and 6h, respectively). Once collected, all tissue biopsies were snap-frozen in liquid nitrogen and subsequently stored at −80°C. The use of these tissues was approved by the London, Surrey Research Ethics Committee (REC reference 17/LO/1801, 26/10/2017).
Endometrial biopsies were obtained from live donors with written, informed consent from all participants from multiple centers.
Endometrial biopsies for sequencing were obtained from individuals recruited from Newcastle Upon Tyne Hospitals after ethical approval (reference 16/NE/0167, North East–Newcastle & North Tyneside 1 Research Ethics Committee).
Proliferative endometrial biopsies for deriving organoids were obtained from Addenbrooke’s Hospital under ethical approval from the East of England–Cambridge South Research Ethics Committee (08/H0305/40).
Endometrial scratch samples from secretory-phase endometrium for deriving organoids were obtained from Bourne Hall Clinic under ethical approval from the East of England–Central Research Ethics Committee for the ‘Biology of the Human Uterus in Pregnancy and Disease Tissue Bank’ run by the Centre for Trophoblast Research (17/EE/0151).
Endometrial biopsies were obtained using a disposable endometrial cell sampler, starting from the uterine fundus and moving downward to the internal cervical ostium. None of the participants were on hormonal treatments for at least 3months before the procedure.
Endometrial tissues were staged based on standard histological criteria.
Tissue dissociation for all fresh tissues was conducted within 24h of tissue retrieval in a two-step digestion protocol. Briefly, endometrial tissue was treated with collagenase to retrieve the stromal fraction (step 1)45. After filtering, pieces of tissue retained on the 100-µm filter were treated with trypsin to enrich for the epithelial fraction (step 2)46. Nuclei were released via Dounce homogenization as described previously47.
Endometrial organoid cultures
Endometrial organoids were grown as previously described10. Briefly, organoids were cultured at 37°C in a humidified atmosphere of 5% CO2. The medium was refreshed every 2–3d, and the organoids were passaged at an average ratio of 1:3 every 5–7d. The organoid suspension was centrifuged for 6min at 600g between passaging steps. The passaged organoid pellet was resuspended in 25-μl ice-cold Matrigel (Corning, 356231) droplets, plated in a 48-well plate (Costar, 3548), allowed to solidify at 37°C for 15–30min and covered with 250μl of endometrial organoid ExM. Components of ExM for culturing human endometrial organoids are available in Supplementary Table 14.
Hormonal stimulation and inhibition experiment of endometrial organoids
Endometrial organoids were stimulated with hormones and treated with NOTCH γ-secretase inhibitors (DBZ, Tocris 4489 and DAPT, Tocris 2634) as well as WNT inhibitors (tankyrase inhibitor XAV939, Tocris 3748 and porcupine inhibitor IWP-2, Tocris 3533) for 6d. First, 10,000 single cells were plated per 25-μl Matrigel droplet into a 48-well plate with ExM supplemented with Rho kinase inhibitor (Y-27632-CAS 146986-50-7) and CHIR 99021. At 10d after plating, organoids were primed with 10nM estrogen (E2) and treated with NOTCH and WNT inhibitors (20μM DAPT, 1μM DBZ, 2μM IWP-2, 2μM XAV939 in the ExM). R-spondin-1, a WNT signaling activator, was depleted from the ExM in the conditions where WNT inhibitors were used. After 48h, they were stimulated with 10nM E2, 1μM progesterone (P4), 100μgml−1 cAMP and 20ngml−1 prolactin while still being treated with NOTCH and WNT inhibitors (20μM DAPT, 1μM DBZ, 2μM IWP-2, 2μM XAV939). Conditions in which (1) inhibitors but no hormones, (2) no inhibitors but hormones and (3) no inhibitors and no hormones were added were used as controls.
Doublet detection, alignment of data across different batches and clustering of scRNA-seq data
The 10x Genomics scRNA-seq data were analyzed with Scanpy48, with the pipeline following their recommended standard practices. In addition, we implemented a number of enhancements described below.
Individual samples of single cells or single nuclei were initially analyzed separately before being batch corrected into an integrated dataset and had two-step diffusion doublet identification performed49,50. The first step was performed with Scrublet51 (v.0.2.1) on a per-sample basis, with the scores diffused by overclustering the cells and reporting each cluster’s median value. Doublets were identified from a distribution of these scores centered at the median and using a mean absolute standard deviation estimate, with statistically significant cells after FDR correction flagged as doublets. The second diffusion step takes place in a joint multi-sample manifold, with the frequency of identified doublets in granular (Leiden resolution 10) clusters serving as the basis for the distribution and the statistical significance analysis described in ref. 50, with Bonferroni for FDR correction and a significance threshold of 0.01. In the organoid samples, additional cells were identified as doublets by Souporcell52 when multiple genotypes were found in single droplets.
After filtering out cells with fewer than 500 genes and more than 15% mitochondrial reads (20% for organoids), samples were integrated using scVI53 (v.0.6.5). While the raw count matrices were used for single cells, the counts from single nuclei were denoised from ambient RNA before the manifold identification. For that task, decontX54 (under the ‘celda’ R package v.1.5.11) was used on each sample separately. We then excluded cell cycling genes from G2/M and S phases listed inside the Seurat package55. For the in vivo dataset, the expression levels of the 5,000 genes that were identified by the scVI native method were modeled by its generative model with 64 latent variables for 500 iterations. Epithelial, endothelial and immune cells were subsequently reanalyzed separately by using 16 latent variables. The resulting latent variables were used for neighbor identification for Leiden clustering56 and uniform manifold approximation and projection (UMAP) visualization. The resulting clusters specific to a single donor, or that had lower numbers of genes expressed or a lower percentage of mitochondrial expression, were excluded. In the epithelial and immune cell reanalysis, cells exhibiting high doublet scores were excluded, and epithelial cells were subsampled to balance donor contribution (Annotation of organoid data using scRNA-seq in vivo dataset reference).
Annotation of scRNA-seq datasets
Identification and labeling of the major cell types in the in vivo dataset was by manual inspection of marker genes and interpretation of these based on the literature. Cluster-specific marker genes were defined using two approaches. To account for the donor effect, we first used DEseq2 (ref. 57), where cells were aggregated into in silico mini-bulks by summing the raw expression of single cells separately according to their donor origin. Every cell-type-specific mini-bulk was compared against a matched mini-bulk that corresponded to all other cells from the same donor; however, such mini-bulk pairs defined by aggregating less than ten cells were excluded (for example, clusters associated with a specific menstrual stage). Secondly, the Wilcoxon test was used to report genes that were differentially expressed. To account for change of sequencing depth, cells were partitioned into four groups corresponding to the quartiles of the sequencing depth of cells considered (independent of donors), and the four resulting Z-scores were combined. P values were adjusted with the Benjamini–Hochberg method.
For each cell, we estimated the cell cycle phase (G1, S or G2/M) based on its expression of G2/M and S phase markers, following the method described in ref. 58 and implemented in Scanpy score_genes_cell_cycle function.
Efficiency of organoid differentiation
To identify clusters predominantly appearing in organoid cultures upon treatment with WNT or NOTCH inhibitors, we evaluated the proportion of cells in each cluster coming from the organoids with and without inhibitors using Fisher’s exact test. Odds ratios were computed against the control organoids (grown with no inhibitors) at the matched time points on a per-genotype basis. In addition, the robustness of the observed effect was evaluated using a paired t-test that compares these fractions of cells from organoids that were treated with inhibitor with their respective controls.
Annotation of organoid data using scRNA-seq in vivo dataset reference
To identify transcriptomic similarities with the in vivo scRNA-seq epithelial subset, we used a regularized logistic regression approach. To best link variation in gene expression in organoids to changes in donors with known stages of the menstrual cycle, we subsampled the in vivo epithelial cells up to 1,000 cells per donor, which balances the cell number representing each stage. To limit the influence of cell cycle, we excluded the SOX9 proliferative cluster composed of a majority of cells at G2/M or S phase, and also excluded the G2/M and S genes from Seurat. Subsequently, we subset both datasets to their shared highly variable genes and further prune half of remaining genes that are not cell-type-specific according to three heuristic measures (fold-increase, fold-increase×fraction-positive, fold-increase×fraction-positive0.5)59. Gene expression was log-transformed and normalized by the maximum RNA expression for both datasets. The model was trained with the in vivo epithelial identities, with 10,000 iterations, and used to classify organoid cells. To visualize results as radial projection we followed La Manno et al.37, where the overlaid position of each cell corresponds to the weighted average of radially balanced unitary vectors (each pointing toward a different corner of a regular polygon), where the weight is each posterior probability of a cell to belong to a given cell type. The latter transformation is a softmax function where each component is multiplied by 15 before being exponentiated.
Annotation of snRNA-seq using scRNA-seq in vivo dataset as a reference
A quantitative measure of the similarity was reported by evaluating the cosine distance from in vivo single nuclei to the centroids of expression defined by the cell-type clusters identified in in vivo single cells, where each vector compared either contains expression of all selected genes in a single cell or their respective mean expression for a given cell-type cluster. Results were visualized with a radial projection as described previously by La Manno et al.37.
Trajectory analysis of organoid data
To identify the branching point where cells commit to a ciliated or secretory fate, we first identified cell clusters in the organoid cultures. Then, we used Palantir60 (v.1.0.0) on cells from the clusters that did exhibit a mean expression level of the progesterone receptor higher than 0.2. A randomly selected cell corresponding to proliferating cells at day 2 was selected as a cell of origin.
Cellular signal analysis
Tumor bulk transcriptomes for endometrioid (430 samples) and serous (112 samples) endometrial adenocarcinoma were downloaded from TCGA. Cellular signal analysis was then applied61 to identify the major transcriptional programs used by tumor cells based on our single-cell endometrial atlas. This method fits the raw bulk messenger RNA counts to a weighted linear combination of transcriptomic signals derived from reference single-cell data. To limit the effect of the cell cycle, we only included cells in the G1 phase and excluded a proliferative SOX9 cluster. Proportions of epithelial-derived signals in the bulk samples were computed as the fraction of samples for which signals (exposures) derived from an epithelial cell cluster exceed the intercept term of the model. Clinical data associated with the samples that exhibited exposure for SOX9+LGR5+ or SOX9+LGR5− above the intercept value were further investigated; we used the Kruskal–Wallis test to confirm that SOX9+LGR5+ signature contribution (exposure) in tumors differs in different cancer stages (as defined in TCGA clinical data), and used the Wilcoxon test and t-test to assess the significance of the increase in exposure noted in later stages of cancer, as well as several clinical data (Supplementary Table 3).
Location of cell types in Visium data
To spatially map cell types defined by scRNA-seq analysis within the Visium spatial transcriptomics data, we used cell2location (ref. 18) (v.0.5-alpha). Briefly, cell2location decomposes multi-cell spatial transcriptomics data into cell-type abundance estimates in a spatially resolved manner. First, the model derives expression signatures of cell types by calculating average expression counts of each gene in each cell type in the raw count scRNA-seq data, selecting genes expressed in at least three cells. Next, to obtain cell-type locations, the model performs a hierarchical non-negative decomposition of the gene expression profiles at spatial locations (spots with multiple cells) into the reference signatures. Each Visium section was analyzed separately with parameters set to default values, except train_args=‘n_iter’: 30,000; posterior_args=‘n_samples’: 1,000; model_kwargs=‘cell_number_prior’: {‘cells_per_spot’: 8, ‘factors_per_spot’: 4}; and ‘gene_level_prior’: {‘mean’: 1/2, ‘sd’: 1/4, ‘mean_var_ratio’: 1}. We visualize the absolute amount of mRNA contributed by each cell population to each spot. We used a 5% percentile of the posterior distribution of this parameter (mRNA counts), representing the number of mRNA molecules confidently assigned to each cell type.
Clustering of spots
Visium data were processed using Scanpy48 following the recommended tutorial with normalization using a scaling factor of 10,000; log transformation; variable gene detection with ‘Seurat’ flavor; principal-component analysis (PCA); neighborhood graph building; and UMAP calculation. Each sample was analyzed independently.
Clusters were defined by the Louvain algorithm and assigned as myometrium or endometrium based on visual inspection of the H&E image of the tissue aligned with each spot. The cluster of spots corresponding to epithelial cells in the endometrium for sample A30, 152807 slide, was further clustered using the same approach. One of the spot subpopulations was excluded due to the low percentage of epithelial cells in the spot after visual inspection. The other ‘epithelial’ spot subpopulations were labeled based on the endometrial layer in which they were found—basal, lumenal and glandular.
Differentially expressed genes in each subpopulation of ‘epithelial’ spots were calculated using the limma R package62.
Calculating TF activities
TF activities were estimated via the combined expression levels of their targets. Target genes were retrieved from Dorothea32, where TF–target relationships are scored from A to E, with decreasing confidence. Here we updated Dorothea regulons as follows: first, synonymous gene names were corrected; second, bona fide TF–target relationships manually curated from Uniprot were added as a new curated source; third, signed and curated interactions were upgraded to score B; fourth, TRRUST63 curated interactions were updated to v.2_20180416 version and signed interactions supported by more than one PubMed source were upgraded to A; and, finally, we created a new category (AA) for the most trustable TF–target interactions that were either detected by all approaches or in more than two curated resources. For each TF, we used the highest scored set of targets with at least ten target genes, as in the original publication32.
Next, we estimated TF activities by performing a Gene Set Enrichment Analysis (GSEA)-like analysis of the gene expression signatures of each cluster resulting from the Wilcoxon test (Annotation of scRNA-seq datasets) with the msVIPER function in the Viper R package64 v.1.22.0. New regulons are available in Supplementary Table 15.
CellPhoneDB v.3.0
To study the interactions between epithelial and other cell populations identified in our endometrial samples, we updated our CellPhoneDB approach to v.3.0 (ref. 33). First, we retrieved the interacting pairs of ligands and receptors satisfying the following criteria: (1) all the members were expressed in at least 10% of the cells in the cluster under consideration and (2) at least one of the members in the ligand or the receptor is a differentially expressed gene (Wilcoxon test; Annotation of scRNA-seq datasets). To account for the distinct temporal and spatial location of cells (that is, microenvironment), we further classified the epithelial interactions based on (1) the phase of the menstrual cycle where cell subsets coexist and (2) their location in the three main endometrial layers (luminal, glandular and basal) according to cell2location (ref. 18).
To account for the complexity of WNT cell–cell signaling, several ligands and functional heteromeric receptors were further curated manually and re-annotated in the CellPhoneDB database (Supplementary Table 16).
Image stitching and manual annotation of selected glands
Confocal image stacks were stitched as two-dimensional maximum intensity projections using the BIOP Perkin Elmer Acapella Stitcher (EPFL, Lausanne; https://www.perkinelmer.com/PDFs/downloads/TCH-Workflows-In-Depth-High-Content-Analysis-Operetta.pdf).
smFISH quantification
RNA spot quantification of smFISH targeting WNT7A and NOTCH2 was achieved with a three-step process:
Step 1. Segmenting glands within the tissue. Ilastik65 was used to train a random-forest-based pixel classifier to detect valid gland areas based on the Nuclear-DAPI channel and the Gland-EPCAM (IHC) channel. Three rounds of ilastik classification were used to achieve adequate rejection of off-target signal to segment only the glands.
Step 2. Segmenting RNA spots within the glands. Another ilastik pixel classifier was used on the spot channels to segment areas that corresponded to genuine spots that were situated in the glands segmented in step 1. The spots were verified by only including spots visible in one channel only. This was done to remove blood inclusions, which gave a confounding signal across multiple channels.
Step 3. The edge of the lumen was manually annotated using napari66. The distance of each pixel on the image was calculated to the nearest point on the lumen edge. Then the total fluorescence intensity was measured for spots in glands and binned into intervals of distance away from the lumen. The gland area was also calculated for each distance interval. The spot fluorescence was divided by the gland interval to give a value of spot intensity that was normalized by area.
Reporting Summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article.
Online content
Any methods, additional references, Nature Research reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at 10.1038/s41588-021-00972-2.
Supplementary information
Supplementary Methods, statistics and reproducibility, and supplementary references.
Acknowledgements
This publication is part of the Human Cell Atlas; www.humancellatlas.org/publications. We thank the Sanger Cellular Generation and Phenotyping (CGaP) Core Facility and Sanger Core Sequencing pipeline for support with sample processing and sequencing library preparation; A. Wilbrey-Clark for coordinating sequencing experiments; A. Garcia for graphical images; A. Kimber and J. Riches for patient consent; K. Elder and T. Cindrova-Davies for endometrial samples and derivation of organoid cultures; I. Martincorena, G. Burton, F. Vieira Braga and A. Goncalves for discussions; N. Huang and A. Aivazidis for data analysis; and S. Williams and the 10x R&D team for help with Visium analysis. The material from the deceased organ donors was provided by the Cambridge Biorepository for Translational Medicine. The endometrial biopsies were obtained from the Newcastle Uteroplacental Bank. We thank the donors and donor families for granting access to the tissue samples as well as the Specialist Nurses on Organ Donation (SNODs). This work was supported by the Wellcome Trust grant ‘Wellcome Strategic Support Science award’ (grant no. 211276/Z/18/Z); the Wellcome Sanger core funding (grant no. WT206194); the Wellcome Trust joint investigator award (grant no. 200841/Z/16/Z); the MRC Human Cell Atlas (grant no. MR/S036350/1); the Royal Society Dorothy Hodgkin Fellowship (grant no. DH160216); the Royal Society Research Grant (grant no. RG93116); the Centre for Trophoblast Research; the European Union’s Horizon 2020 research and innovation program under grant agreement no. 874867; the European Research Council grant no. 853546, and the Chan Zuckerberg Initiative DAF, an advised fund of Silicon Valley Community Foundation (grant no. 2018-190766/RG98793). L.M. is a recipient of the Jean Shank/Pathological Society of Great Britain and Ireland Intermediate Research Fellowship (grant reference no. 1175).
Extended data
Source data
Source Data Fig. 4(25K, csv)
Statistical source data.
Source Data Extended Data Fig. 9(43K, xlsx)Statistical source data.
Author contributions
R.V.-T. and M.Y.T. designed and supervised the experiments and analysis with contributions from S.A.T. and O.A.B.; R.V.-T. and C.S.-S. performed sampling and library prep with help from L.M., R.H., T.P., E.P., M.D., K.T.M., P.A., M.R.S. and K.S.-P. L.G.-A., L.-F.H. and R.V.-T. performed computational analyses with contribution from V.L., H.M., K.P., J.P., S.v.D., M.P., V.Y.K., V.K. and A.A. K.R., B.W., E.T., T.L., K.K., A.T., C.I.M., V.V. and O.A.B. performed and analyzed the smFISH; L.G. and C.L. performed immunohistochemistry; K.N., R.C.F. and M.Y.T. performed the organoid experiments; L.M., M.J.-L., L.C. and A.M. provided pathological expertise; R.V.-T., A.M. and M.Y.T. wrote the manuscript with contributions from L.G.-A., K.R., K.N., L.M., O.A.B. and S.A.T. All authors read and approved the manuscript.
Data availability
Datasets were uploaded into ArrayExpress under accession numbers E-MTAB-10287 (scRNA-seq in vivo), E-MTAB-9260 (Visium in vivo) and E-MTAB-10283 (scRNA-seq in vitro). snRNA-seq in vivo data were uploaded into the European Genome-phenome Archive; EGAD00001007909. Requests for data access should be sent directly to the Data Access Committee of this work: datasharing@sanger.ac.uk (https://ega-archive.org/datasets/EGAD00001007909). Tumor bulk transcriptomes for endometrioid and serous endometrial adenocarcinomas were downloaded from TCGA. Additional single-cell transcriptomes of ten endometrial biopsies were downloaded from the Gene Expression Omnibus with GSE111976 (ref. 22). Processed matrices can be accessed and downloaded from www.reproductivecellatlas.org. Image datasets are available at the EMBL-EBI BioImage Archive under accession number S-BIAD190. Data availability includes all datasets in the manuscript. Source data are provided with this paper.
Code availability
Code is available at https://github.com/Ventolab/UHCA. Source code and tutorials for CellPhoneDB v.3 are available at https://github.com/Ventolab/CellphoneDB.
Competing interests
In the past 3years, S.A.T. has worked as a consultant for Genentech, Roche and Transition Bio, is a remunerated member of the scientific advisory boards of Biogen, GlaxoSmithKline, Foresite Labs and Qiagen and is an equity holder of Transition Bio. The remaining authors declare no competing interests.
Footnotes
Peer review information Nature Genetics thanks Shalev Itzkovitz and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
These authors contributed equally: Luz Garcia-Alonso, Louis-François Handfield, Kenny Roberts, Konstantina Nikolakopoulou.
These authors jointly supervised this work: Omer A. Bayraktar, Sarah A. Teichmann, Margherita Y. Turco, Roser Vento-Tormo.
Change history
12/20/2022
A Correction to this paper has been published: 10.1038/s41588-022-01287-6
Contributor Information
Omer A. Bayraktar, Email: ku.ca.regnas@5bo.
Sarah A. Teichmann, Email: ku.ca.regnas@9ts.
Margherita Y. Turco, Email: hc.imf@ocrut.atirehgram.
Roser Vento-Tormo, Email: ku.ca.regnas@4vr.
Extended data
is available for this paper at 10.1038/s41588-021-00972-2.
Supplementary information
The online version contains supplementary material available at 10.1038/s41588-021-00972-2.
References
Full text links
Read article at publisher's site: https://doi.org/10.1038/s41588-021-00972-2
Read article for free, from open access legal sources, via Unpaywall: https://www.nature.com/articles/s41588-021-00972-2.pdf
Citations & impact
Impact metrics
Citations of article over time
Alternative metrics
Discover the attention surrounding your research
https://www.altmetric.com/details/117985990
Smart citations by scite.ai
Explore citation contexts and check if this article has been
supported or disputed.
https://scite.ai/reports/10.1038/s41588-021-00972-2
Article citations
Cellular heterogeneity and dynamics of the human uterus in healthy premenopausal women.
Proc Natl Acad Sci U S A, 121(45):e2404775121, 29 Oct 2024
Cited by: 1 article | PMID: 39471215
Current advances in understanding endometrial epithelial cell biology and therapeutic applications for intrauterine adhesion.
Stem Cell Res Ther, 15(1):379, 25 Oct 2024
Cited by: 0 articles | PMID: 39456113 | PMCID: PMC11515228
Review Free full text in Europe PMC
Engineering human immune organoids for translational immunology.
Bioact Mater, 44:164-183, 18 Oct 2024
Cited by: 0 articles | PMID: 39474315 | PMCID: PMC11513670
Review Free full text in Europe PMC
Organoid modeling meets cancers of female reproductive tract.
Cell Death Discov, 10(1):410, 27 Sep 2024
Cited by: 0 articles | PMID: 39333482 | PMCID: PMC11437045
Review Free full text in Europe PMC
Combined assembloid modeling and 3D whole-organ mapping captures the microanatomy and function of the human fallopian tube.
Sci Adv, 10(39):eadp6285, 27 Sep 2024
Cited by: 0 articles | PMID: 39331707 | PMCID: PMC11430475
Go to all (190) article citations
Data
Data behind the article
This data has been text mined from the article, or deposited into data resources.
BioStudies: supplemental material and supporting data
GEO - Gene Expression Omnibus
- (3 citations) GEO - GSE141549
Similar Articles
To arrive at the top five similar articles we use a word-weighted algorithm to compare words from the Title and Abstract of each citation.
Estrogen Signaling Drives Ciliogenesis in Human Endometrial Organoids.
Endocrinology, 160(10):2282-2297, 01 Oct 2019
Cited by: 42 articles | PMID: 31290979
Long-term, hormone-responsive organoid cultures of human endometrium in a chemically defined medium.
Nat Cell Biol, 19(5):568-577, 10 Apr 2017
Cited by: 299 articles | PMID: 28394884 | PMCID: PMC5410172
Wnt/Β-catenin and sex hormone signaling in endometrial homeostasis and cancer.
Oncotarget, 1(7):674-684, 01 Nov 2010
Cited by: 64 articles | PMID: 21317462 | PMCID: PMC3248134
Review Free full text in Europe PMC
Endometrial Organoids: A Rising Star for Research on Endometrial Development and Associated Diseases.
Reprod Sci, 28(6):1626-1636, 02 Feb 2021
Cited by: 12 articles | PMID: 33533008
Review
Funding
Funders who supported this work.
European Research Council (1)
Grant ID: 853546
Medical Research Council (1)
Cell atlas of the human female reproductive system across the lifespan
Dr Sarah Teichmann, Wellcome Trust Sanger Institute
Grant ID: MR/S036350/1
Wellcome Trust (3)
Grant ID: WT206194
Grant ID: 211276/Z/18/Z
Maternal NK cell recognition of the placenta determines reproductive outcome
Dr Francesco Colucci, University of Cambridge
Grant ID: 200841/Z/16/Z