Single-cell analysis of endometriosis reveals a coordinated transcriptional programme driving immunotolerance and angiogenesis across eutopic and ectopic tissues

Abstract

Endometriosis is characterized by the growth of endometrial-like tissue outside the uterus. It affects many women during their reproductive age, causing years of pelvic pain and potential infertility. Its pathophysiology remains largely unknown, which limits early diagnosis and treatment. We characterized peritoneal and ovarian lesions at single-cell transcriptome resolution and compared them to matched eutopic endometrium, unaffected endometrium and organoids derived from these tissues, generating data on over 122,000 cells across 14 individuals. We spatially localized many of the cell types using imaging mass cytometry. We identify a perivascular mural cell specific to the peritoneal lesions, with dual roles in angiogenesis promotion and immune cell trafficking. We define an immunotolerant peritoneal niche, fundamental differences in eutopic endometrium and between lesion microenvironments and an unreported progenitor-like epithelial cell subpopulation. Altogether, this study provides a holistic view of the endometriosis microenvironment that represents a comprehensive cell atlas of the disease in individuals undergoing hormonal treatment, providing essential information for future therapeutics and diagnostics.

Extended Data Fig. 1. Overview of experiment design and comparison of bulk RNA-seq and scRNA-seq transcriptomic profiles from Ctrl and endometriosis tissues.

a, Experimental workflow. b, UMAP showing distribution of cell based on tissue types, PID, and endometriosis stage, before and after batch correction with Harmony. c, Box plot showing Spearman’s correlation rank (ρ) between bulkRNA-seq and pseudobulk from scRNA-seq in Eutopic (Ctrl & EuE, n = 144 ), Peritoneal (EcP & EcPA, n =90 ), or Ovary (EcO, n =24 ). Each dot represents a sample pair. The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. d, Scatterplot showing distribution of average gene expression between bulk RNA-seq and scRNA-seq (Spearman ρ). Each dot represents a gene. e, Volcano plots representing DEGs between scRNA-seq pseudo bulk (red) and bulk RNA-seq from undissociated tissue (blue) (edgeR, FDR < 0.001, LogFC > 3). The genes highlighted are exclusively expressed in bulk RNA-seq and associated with erythrocytes (orange), neuronal projections (green), adipocytes (brown), and muscle cells (purple). Related to Fig 1.

Extended Data Fig. 2. Proportion of major cell types in each replicate and IMC panel for spatial profiling of Ctrl and endometriosis tissues.

a, Major cell types were determined based on UMAP. The mean distribution for all 5 major cell populations is represented for each tissue type Ctrl, EuE, EcP, EcPA, and EcO (left of the line). Each pie chart represents major cell type proportions for each replicate (right of the line). b, Each antibody was selected according to the cell types identified by the scRNA-seq data analysis. Representative images show single channels for each metal-conjugated antibody in a EuE biopsy. A total of 26 antibodies was used to identify cellular heterogeneity within stromal, endothelial, epithelial, lymphocyte, and myeloid major cell types. Additional antibodies (in “Others”) were used to identify cell proliferation (Ki67), active metabolism (pS6), extracellular matrix (Collagen1), and nuclei (DNA). A complete list of cell subpopulations identified through this panel of markers is listed on Supplementary Table 8b. Related to Fig 1.

Extended Data Fig. 3. Stromal cell analysis across sample types.

a, Bar plot representing the proportion of stromal cell types in control endometrium and endometriosis lesions. Endometrial fibroblasts were found in all lesions. Fibroblast C7 is the predominant fibroblast type in EcO. b, Density plot showing distribution of mural cells for each tissue. Arrows points to Prv-CCL19. c, Heatmap of markers genes for mural cell subtypes. d, Track plot representing gene expression pattern for selected DEG in Prv-CCL19 subpopulations. GGT5 and ABCC9 are pan-markers for this cell subtype. e, Box plot showing the proportion of CCL19-expressing cells in Prv-CCL19 subpopulation within each tissue type. Each dot represents the percentage of CCL19+ cells in a tissue biopsy (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. Related to Fig. 2.

Extended Data Fig. 4. Characterization of endothelial cells (EC) across sample types.

a, Unique cell-to-cell interaction counts obtained from a modified CellPhoneDB procedure. To recover meaningful interactions, we analyzed ligand-receptor interaction in each sample independently. Unique interactions in each tissue type are counted as follows; each ligand-receptor pair observed in a specific cell type pair is counted as one interaction; this is tabulated for all possible pairwise cell type combinations (up to 58 subpopulations in this study) within a sample (n). The total count (Σ, n_celltype_pairs) represents the commonality of the ligand-receptor interaction of interest. The more common interactions (observed in multiple cell type pairs and in all individual samples) will have higher counts while restricted interactions (observed in specific cell type pairs) will have lower counts. We arbitrarily restricted our analysis to interactions observed fewer than 150 times to narrow the scope of analysis and focus on potentially uncovering unique cell-to-cell interactions. b, Box plot showing the proportion of DLL4-expressing cells in EC-tip subpopulation within each tissue type. c, Density plot showing distribution of endothelial cells for each tissue. d, EC proportions by sample type. e, AQP1+ cell abundance is substantially increased in peritoneal lesions (EcP and EcPA). f, (top) Proportion of aPCV among ECs across tissue types. (bottom) Swarm plot showing AQP1 expression per cell. Horizontal lines represent the median value. For box plots, each dot represents percentage of DLL4+ cells in EC-tip cluster (b) or AQP1+ cells in EC-aPCV cluster (e), in a tissue biopsy (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. Related to Fig. 3.

Extended Data Fig. 5. Myeloid cell diversity in control and endometriosis.

a, Heatmap representing marker genes for each myeloid subpopulation. b, Dendrogram showing the hierarchical clustering (Pearson correlation) for the myeloid cell clusters. c, Bar plot showing the representation of each myeloid subtype across tissue types. Related to Fig. 4.

Extended Data Fig. 6. DC subpopulations.

a, Bar plot represents the proportion of DCs among all myeloid cells for each patient (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). Patient-to-patient variability was observed in DC proportions within the myeloid population and across different tissue types. The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. b, PAGA and RNA velocity trajectory analyses suggest that pre-cDC2 differentiate towards cDC2 and DC3 in Ctrl and EuE. Red arrows indicate that some cDC2 and DC3 cells derive from a smaller intermediate cell population. c, Cell cycle analysis for pre-cDC2, cDC2 and DC3 populations. d, Expression of DC progenitor markers FLT3, AXL, and SIGLEC6. e, Phagocytosis pathway is enriched in cDC2 subpopulations of peritoneal lesions. Bar plot shows the Normalized Enrichment Score (NES) for the top 10- GSEA pathways in cDC2 cells of EuE and EcPA (FDR < 0.1). Related to Fig. 5.

Extended Data Fig. 7. Lymphocyte subpopulations in control and endometriosis tissues.

a, Density plot showing distribution of lymphocyte cells for each tissue. b, Dot plot representing marker genes for each lymphocyte subpopulation, including four natural killer cell (NK) clusters, innate lymphoid cells (ILCs), effector memory T-cells (T_EM), cytotoxic T-lymphocytes (CTL), naïve/central memory T-cells (T_N/T_CM), T regulatory cells (T_Reg), CD4- and CD8- tissue resident T cell (CD4-T_RM and CD8-T_RM, respectively), CD8− mucosal-associated invariant T cells (CD8-MAIT), plasma cells, and B cells. c, Representative IMC images showing the presence and proximity of myeloid cells labelled with CD68 (yellow) with T cells labelled with CD3 (cyan), and T_Reg labelled with FOXP3 (magenta) in EcO (n = 5); nuclei are marked with DNA intercalation (blue). Scale bar = 100 μm. d, Proportion bar plot of CTLA4 expressing cells from the total T_Reg subpopulation. e, Proportion box plot of BCL6, SEMA4A, CXCR5 expressing cells from the total B cells within each sample type. For box plots, each dot represents a unique patient (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. Related to Fig. 6.

Extended Data Fig. 8. Characterization of in vivo epithelial and in vitro endometrial epithelial organoid (EEO) cells.

a, Proportions of epithelial subpopulations per sample type. b, Representative IMC images of MUC5B+ epithelial cells in eutopic endometrium (Ctrl: C07, EuE: E12, E06) from multiple tissues. Epithelial cells are marked with PanCK, EpCAM, E-cadherin (green); MUC5B (magenta); nuclei (white). Scale bar = 100 μm. c, Proportion box plot of SAA1 expressing cells from the total MUC5B+ cells within each sample type. Each dot represents a unique patient (Ctrl n = 3, EuE n=9, EcP n=8, EcPA n=6, EcO n=4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. d, Sequencing metrics from EEO scRNA-seq; UMIs and unique genes counts are shown for Control (C) and endometriosis (E) patients and across tissue type. Undetermined (UD) group represents single cells which could not be assigned due to the lack of multiplexing hashtag but otherwise passed QC. e, Density plot showing distribution of EEO cells derived from Ctrl, EuE, EcP, and EcPA (UD cells were not included). f, UMAP showing the co-expression of MUC5B, SAA1, TFF3, and RUNX3 in the MUC5B+ population comprising in vivo epithelial cells and EEO. Related to Fig. 7.

Extended Data Fig. 9. MUC5B+ cells display a progenitor-like capacity in in vitro organoid culture.

A, Schematic and FACS sorting gating strategy to isolate MUC5B+ and MUC5B- epithelial cells from eutopic tissue for organoid generation. B, Representative brightfield images showing the progression of organoid generation from sorted single cells at day 2, 6 and 10. MUC5B+ cells formed EEO faster than MUC5B- cells. Each panel shows a whole Matrigel dome and magnified organoids are shown in the inset. Inset scale bar = 100 μm. c, Line graph showing area (top) and number (bottom) of EEO generated from MUC5B+ (dark blue, n=1) and MUC5B- (sky blue, n=1) cells over time. Area and Count of EEO is significantly higher in MUC5B+ compared to MUC5B- (paired t-test, two-tailed p < 0.0001). d, IF staining of EEO generated from MUC5B+ (n=1) and MUC5B- (n=1) sorted cells showing the co-localization of endometrial epithelial (E-Cadherin, in green) and MUC5B+ (magenta) staining. Nuclei were counterstained with DAPI (gray). Scale bar = 100 μm. Related to Fig. 7.

Extended Data Fig. 10. Schematic illustrating the proposed microenvironment alterations for ectopic peritoneal and ovary lesions.

In peritoneal lesion (left), the proportion of myeloid and endothelial is increased, and endometrial-like epithelial population is reduced. CCL19-expressing perivascular cells mediate immune cell recruitment, such as macrophages and T cells, which contributes to the immunomodulatory microenvironment. We observe the presence of MSR1-expressing dendritic cells contributing to immunomodulation. TLS is also observed in some lesions. In addition, Mø1-LYVE1 and perivascular cells contribute to angiogenesis by regulating endothelial tip proliferation. In contrast, ovarian ectopic lesions (right) show a striking increase in the proportion of stromal cell and a reduced endometrial-like-epithelial cell presence. The immunomodulatory microenvironment is mainly driven by Mø1-LYVE1 expressing IL10. In ovary lesions, the regulation of angiogenesis is marked by endothelial cell arrest, resulting in mature vasculature. Created with Biorender.com.

Fig. 1 |. scRNA-seq from control and endometriosis patient.

a, Schematic and photographic representation of collected tissue biopsies. Control (Ctrl) specimens were obtained from eutopic endometrium of women without endometriosis. Eutopic endometrium (EuE), ectopic peritoneal endometriosis (EcP), ectopic peritoneal adjacent lesion (EcPA), and ectopic ovarian endometriosis (EcO) were obtained from women with endometriosis. Peritoneal lesions were collected with a surrounding margin of up to 1 cm². The margin (EcPA) was separated upon macroscopic tissue assessment from the lesion (EcP) when possible and before single cell dissociation, as depicted in the representative image for one biopsies. b, Diagram showing scRNA-seq metrics per patient (left) and tissue type (right) after QC. These metrics indicate unique molecular identifier (UMIs) and total genes per cells across patients and tissue types. The cord diagram (center) indicates the representation of each patient (E: endometriosis, C: control) in each tissue type. c, Violin plot representing marker gene expression for each major cell type identified in the scRNA-seq dataset. d, UMAP plot showing the 108,497 single-cells from control and endometriosis tissues. Five major cell types are identified (center UMAP plot) and subsequently subclustered into 58 subpopulations (radial UMAP plots). Each subpopulation was identified using marker genes curated form the literature. The presence of basophils and neutrophils (arrows) indicate that the cell recovery workflow was well-suited to capture delicate cell types known to be easily lost during tissue dissociation. e, Diagram showing number of major cell types (bar plot) and the cell type proportion in each tissue type (heatmap plot). Cell proportions are indicated within each square. Unique combinations of cell markers from each major cell cluster were used to design an IMC panel. Assigned colors represent each major cell type identified in EcP (f) and EcO (g). White arrows indicate endometriotic epithelial glands. Scale bar = 100 μm.

Fig. 2. |. Cellular composition of Ctrl and endometriosis eutopic endometrium.

a, Representative H&E images of eutopic tissues from Ctrl (n=5) and EuE (n=8), before and after classification. b, Box plot showing the proportion of epithelial cells in endometrial tissue. Ctrl (n=5) tissues show significantly higher epithelial cell proportion compared to EuE (n=8), Welch’s T-test, two-sided, p-value = 0.013. Each dot represents a tissue section; See Source Data Fig2. Box represents the interquartile range, whiskers represent min and max, and box centerline represents median. c, UMAP representation for the expression of TOP2A in Ctrl and EuE in global clustering (top). Circle denotes the stromal cell population. Representation of TOP2A (proliferating cells) and MME (endometrial fibroblasts) expressing cells in Ctrl and EuE stromal subclusters (bottom). Arrows depict all endometrial fibroblast (eF) expressing TOP2A in EuE. d, Proportion of cells in G1, G2M, S cell cycle phases within all eF and in eF2 subpopulation. e, Representative IMC image showing the presence of proliferating cells labelled with KI67 (green) in Ctrl (n=4) and EuE (n=5); epithelial cells marked with PanCK, EpCAM, E-cadherin (magenta); stromal cells marked with COL1A1 and CD10 (yellow/orange); nuclei marked with DNA (blue). f, Matrix plot representing the overall similarity of endometrium biopsies from control and endometriosis (Pearson correlation based on gene expression from each patient). EuE clustered into two groups, each showing an enrichment of fibroblasts or immune cells (E: endometriosis, C: control). g, Violin plot showing significantly upregulated genes (OGN and NES) in EuE relative to Ctrl in decidualized stroma (dS2) subpopulation. h, Representative IMC image confirming increase of OGN (cyan) secretion within stroma (orange) in EuE (n=5) relative to Ctrl (n=4). a, e, h, Both Ctrl and EuE representative images are taken from patient receiving same hormonal treatment, Scale bar = 100 μm.

Fig. 3 |. Role of Stromal cell diversity in angiogenesis and immune trafficking in endometriosis lesions.

a, UMAP plot of the 12 identified stromal subpopulations and classified into three general cell subtypes: endometrial fibroblast (eF), C7 fibroblast (fib C7) and mural cell (n = 42,713 cells). b, Violin plot showing markers of mural cell subpopulations. c, UMAP plot of endothelial cells (EC), represented across 7 subclusters: lymphatic EC (LEC), high endothelial venule (EC-HEV), tip EC (EC-tip), capillary (EC-capillary), post-capillary vein (EC-PCV), activated PCV (EC-aPCV), and arterial (EC-artery). d, (top) Proportion of Prv-CCL19 within stromal cells. A major increase of Prv-CCL19 is observed in EcPA. Bars represent the mean value. (bottom) The swarm plot shows CCL19 expression in individual cells from each lesion. e, Dot plot showing significantly upregulated genes involved in angiogenesis and immune cell trafficking (edgeR, FDR < 0.05) in Prv-CCL19. f, Schematic of mural and EC localization. Larger arteries and veins are unsheathed by VSMC, while smaller vessels (e.g., capillaries) are unsheathed by perivascular cells. Lesions Prv cells increase expression of pro-angiogenic genes when compared to Ctrl. G, Dot plot showing significant representative DEGs involved in new vessel formation in tip EC (edgeR, FDR <0.05). h, Dot plot showing significant DEGs involved in cell adhesion and permeability in a-PCV (edgeR, FDR <0.05). I, Representative IMC image from a peritoneal lesion (n=7). CD3⁺ T-cells (cyan) and CD68⁺ myeloid cells (magenta) localize within and surrounding blood EC vasculature marked by CD31 and AQP1 (yellow). Nuclei counterstained by DNA labeling (blue). Scale bar = 100μm.

Fig. 4 |. Macrophage heterogeneity in control and endometriosis.

a, UMAP plot of myeloid cells, clustered into 15 different subtypes (n = 12,262 cells). b, Dot plot showing expressed marker genes for tissue-resident (TRM), blood-infiltrated, and activated macrophages across identified Mɸ subpopulations and tissues. c, Density plot showing macrophage distribution for in each tissue type. d, UMAP plot showing RNA velocity streamlines for monocytes and macrophages in Ctrl. Streamlines represent the predicted transition path of cells across subpopulations. e, Bar plot showing the proportion of LYVE1-expressing cells to all macrophages within each tissue type. Each dot represent percentage of LYVE1+ cells in a tissue biopsy (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. f, Dot plot showing DEG involved in immunotolerance in Mɸ1-LYVE1 population. g, IMC image from FFPE tissue section of a peritoneal lesion. Images depict LYVE1⁺ macrophages (LYVE1, CD68) localization near endothelial cells (CD31, AQP1) (white arrows). Scale bar =100μm. h, Matrix plot showing expression of pro-inflammatory and pro-tolerogenic related genes in Mɸ4 subpopulation in Ctrl and endometriosis.

Fig. 5 |. Immunomodulatory role of DC in peritoneal endometriosis.

a, Violin plot showing markers of dendritic cell (DC) subpopulations; CD1C expression is prevalent in three DC subpopulations: pre-cDC2, cDC2 and DC3. b, CD1C⁺cells represent majority of the DC population, accounting for more than 83% of DCs in all but ectopic ovary tissue. c, Density plot showing the increased cDC2 populations in peritoneal lesions compared to EuE. d, Expression of cDC2 markers CD207 and CD1A; and proliferation marker TOP2A. e, Proportion of CD207 expressing cells across all cDC2 populations. CD207+ cells were consistently observed in eutopic endometrium, but variable in peritoneal lesions and not observed in ovarian lesions. Each dot represents the percentage of CD207⁺ for each tissue biopsy (Ctrl n = 3, EuE n = 9, EcP n = 8, EcPA n = 6, EcO n = 4). The box represents the interquartile range with median and minimum/maximum represented by box centerline and whiskers, respectively. f, Track plot representing the expression of DEGs upregulated in cDC2-CD1A in EcPA (Wilcoxon, FDR < 0.05). Each bar represents a cell. Differential expression for the represented genes is detected in EcPA cells (black box). g, Density plot of cDC2 from EcP and EcPA showing the distribution of cDC2 on UMAP representing different cell states (left). Scatter plot showing CD207+/MSR1- (n = 237), CD207-/MSR1+ (n = 121), double positive (n = 82) and double negative (n = 141) cells (right). h, Top 12 DEGs between CD207+/MSR1- and CD207-/MSR1+ populations from cDC2 subpopulations in EcP and EcPA (Wilcoxon, FDR < 0.05, logFC > 1).

Fig. 6 |. TLS presence in peritoneal endometriosis.

a, UMAP plot of lymphocyte subpopulations. Represented clustering highlights 14 different subpopulations (n=22,225) based on known markers. b, Schematic showing CD86-CTLA4 ligand-receptor interaction between macrophages Mɸ1 and Tregs. The dot plot shows gene expression for this interacting pair in each tissue type. c, Dot plot showing DEGs associated with T_Reg self-tolerance maintenance (edgeR, FDR <0.05, ^# marks non-significant DEG). d, Violin plot representing ENTPD1 gene expression in tissue-resident NK1 cells across sample types. e, H/E staining from FFPE tissue section of a peritoneal lesion. This patient sample presented TLS-like formation highlighted in the white frame, detected in n=2 out of 7 EcP. f, IMC image from the same lesion showing endometrial fibroblasts (CD10, red), B-cells (CD20, yellow), epithelial cells (Pan-KRT, green), stroma (Col1A1, cyan) and antigen presenting cells (HLA-DR, magenta). TLS are primarily located through an accumulation of CD20⁺ cells forming GC-like structures in the periphery of the lesion (white frame, arrow). HLA-DR overlap with CD20 indicates an antigen presenting capacity within the GC. g, Magnified image showing GC structures with accumulation of B-cells (CD20, yellow) in the center surrounded by T-cells (CD3, cyan). KI67 labels proliferative B-cells within the GC (green, middle). CD31 and AQP1 label blood endothelial cells (green, right panel). PDPN marks follicular dendritic cells (magenta on middle or cyan on right image). h-i, H/E (left) and corresponding IMC (right) representative images from endometriotic lesions without TLS in EcP (n=5/7) (h) and EcO (n= 6/6) (i) for identical antibodies panels in (f). Scale bar =100 μm.

Fig. 7 |. Characterization of epithelial cell subpopulations in Ctrl and endometriosis patients.

a, Unsupervised clustering of epithelial cells led to 10 subpopulations (n=19,200) represented in the UMAP. b, Density plot showing the distribution of epithelial subtypes across tissues. c, Markers for each epithelial subtype and menstrual phase across each epithelial cell subtype. d, Immunohistochemistry (IHC) staining confirms the presence of MUC5B⁺ cells in EuE (Left panel, n = 1). Immunofluorescence (IF) showing co-localization of endometrial epithelial (E-Cadherin⁺, in green) and MUC5B⁺ cells (magenta). Nuclei were counterstained with DAPI (cyan) in EuE (n = 4). Scale bar = 100μm. e, Formyl Peptide Receptor 2 (FPR2) expression is specific to myeloid cells (left), and more precisely to monocytes and Mɸ4-infiltrated cells (right). f, Representative image of endometrial epithelial organoid (EEO) cultures derived from dissociated single cell of endometrium and endometriotic lesions. g, UMAP plot representing the merge dataset for in vivo (tissue derived) and in vitro (EEO) epithelial cells. Classification follows previously described subpopulations in vivo.