Epithelial populations of the mouse prostate

To generate a comprehensive aggregated single-cell dataset for the mouse prostate, we gathered publicly available scRNA-seq datasets generated from C57BL/6 and FVB mice and generated new data focusing on the proximal and periurethral regions of the prostate, which have been less studied (Additional file 2: Table S1). We analyzed these datasets using two independent computational approaches to confirm the reproducibility of our interpretations. In the first approach, we de-noised each dataset using random matrix theory (RMT), which improves the ability to separate and detect rare cell populations [32]. We then sequentially clustered in each dataset to identify cell populations, following the same strategy used previously [4]. This analysis of 11 datasets resulted in an aggregated dataset of 21,952 cells arranged in 30 prostate cell clusters (Fig. 1A, Additional file 1: Fig. S1A, Additional file 1: Fig. S2). As a second approach, we used a standard Seurat pipeline to generate an aggregated dataset from 13 datasets of sufficiently high quality, which was composed of 30,433 cells in 18 distinct clusters (Additional file 1: Fig. S3A).

These parallel approaches allowed for the identification and comparison of cell populations across datasets in a uniform manner, independent of differences in reporting and labeling between publications. To define robust cell populations, we required that the population be identified in at least three independent datasets and have nearly complete overlap in globally distinguishing gene expression. For rare cell populations, we only required that the population be present in at least two independent datasets. Of note, the majority of clusters were identical using both the RMT and Seurat approaches. The RMT approach handled sparse data differently, yielding a greater number of small clusters and providing better discrimination between populations with low cell numbers.

We found that epithelial populations were remarkably consistent across datasets and approaches. Interestingly, no distinct subclusters were formed based on mouse strain background, which did not significantly contribute to prostate epithelial heterogeneity. In particular, basal cells formed a single contiguous cluster in individual datasets (Fig. 1B–F; Additional file 1: Fig. S3B–F), as previously reported [4,5,6,7,8], and in our aggregated datasets (Fig. 1A, Additional file 1: Fig. S3A). We did not observe evidence of a distinct basal subcluster with expression of Zeb1 or other epithelial-mesenchymal transition (EMT) markers [41]. However, in the Seurat pipeline, we observed that a small subset of basal cells adjoins the periurethral (PrU) cluster, is proximally enriched, and expresses slightly more proximal and luminal markers (Additional file 1: Fig. S3A, G).

We identified multiple luminal epithelial clusters, which represent distinct cell types that are separated by prostate region (Fig. 1A–F, Additional file 1: Fig. S3A–F), as previously reported for individual datasets [4,5,6,7,8]. Since the nomenclature for these populations differs between laboratories (summarized in [3]), we follow a descriptive naming system [4] that denotes lobe-specific prostate populations (e.g., LumA for the distal anterior lobe) as well as proximal populations (LumP for proximal prostate). Notably, although the dorsal and lateral lobes have often been combined as a “dorsolateral lobe,” highly distinct dorsal (LumD) and lateral (LumL) populations were always found in each individual dataset as well as the aggregated datasets (Fig. 1A, B, D–F, Additional file 1: Fig. S3A, B, D–F). In contrast, the anterior (LumA) and dorsal (LumD) distal luminal populations consistently displayed the most transcriptomic overlap (Fig. 1A, B, D–F, Additional file 1: Fig. S3A, B, D–F).

Unlike distal luminal cells, which differ by lobe, proximal luminal cells (LumP) formed a single cluster without lobe-specific identity (Fig. 1A–F, Additional file 1: Fig. S3A–F) [4,5,6,7,8]. The vast majority of LumP cells are located in the proximal region of the prostate, though rare distal cells can be observed [4, 5, 8], and functional heterogeneity within the population has been reported [5]. In this regard, in the Seurat pipeline, we observed that a subset of LumP cells is adjacent to distal luminal cells (Additional file 1: Fig. S3A).

Neuroendocrine (NE) cells represent a rare and historically elusive epithelial population that could be detected in both analytical pipelines (Fig. 1, Additional file 1: Fig. S3). Ionocytes are another rare population that was recently described in the prostate [6], and our meta-analysis revealed their presence in additional datasets (Fig. 1A, B, D, Additional file 1: Fig. S3A, B, D) [4, 7]. Though ionocytes have some transcriptional similarities to NE cells, they express Foxi1 and Atp6v1g3 but not specific luminal or basal markers (Additional file 1: Fig. S2I, J). Both cell types were observed in higher proportions in the proximal dataset, and the PrU population is described in detail below.

Using the aggregated datasets, we generated reference gene expression signatures that are specific for each prostate epithelial cell type (Additional file 3: Table S2). In addition, we examined the Gene Set Enrichment Analysis Hallmark signatures and found increased expression of genes involved in protein secretion in seminal vesicle and distal luminal cells, and the lowest levels of Notch signaling genes in NE cells (Additional file 1: Fig. S4). Finally, we observed rare epithelial clusters in individual datasets that may represent cell states. In particular, a subset of LumA cells expresses both LumA and basal markers and may correspond to “intermediate” cells with hybrid luminal and basal features (Additional file 1: Fig. S1B, C).

Non-epithelial cell populations

Our scRNA-seq meta-analysis also provided consistent insights into non-epithelial cell types in the mouse prostate. The mesenchymal/stromal cells present in these datasets are predominantly fibroblasts and can be divided into several different clusters (Fig. 1, Additional file 1: Fig. S2E, F, Additional file 1: Fig. S3). The Mesenchyme 1 (Mes1) population is proximally enriched, lies adjacent to the epithelium, and expresses Srd5a2 as well as many Wnts and other signaling factors, whereas Mesenchyme 2 (Mes2) is enriched more distally, is located slightly farther from the epithelium, and expresses many chemokines and complement components [4, 10]. We also identified distinct myofibroblast and smooth muscle populations that express smooth muscle actin (Acta2), and observed that a subset of myofibroblasts expresses Lgr5 [4, 6, 10, 42]. Although a third fibroblast population has been reported [10], it did not appear as a distinct cell type in our analyses, but rather as a subset of Mes2 (Additional file 1: Fig. S2E, F). Interestingly, several mesenchymal cell types reported to exist in the prostate (e.g., telocytes) were not detected in any dataset, suggesting that the prostate stromal compartment is incompletely captured in existing scRNA-seq data.

Hematopoietic lineage populations (such as B and T lymphocytes, dendritic cells, and NK cells) were also detected across multiple datasets, with the immune compartment displaying a notable myeloid bias. In particular, macrophages divided into distinct subclusters along a continuous spectrum, which was most evident in the RMT pipeline. Since profiles for M1 and M2 macrophages could not be definitively identified, we have named these populations alphabetically (Fig. 1, Additional file 1: Fig. S2A–D). In addition to the macrophage populations, we detected a population with substantial overlap in gene expression to macrophages, which appeared to correspond to differentiating monocytes (Additional file 1: Fig. S2A–D).

Finally, we also observed contaminating seminal vesicle cells across multiple datasets. Seminal vesicle epithelial cells could be clustered into a single basal population as well as luminal populations with more proximal markers or more distal markers (Additional file 1: Fig. S2G, H), suggesting potential epithelial heterogeneity within this tissue.

The periurethral region

We define the periurethral (PrU) region as the most proximal extent of each prostate lobe nearest the junction with the urethra. PrU cells make up most of the epithelium in this region. Because this region is located exclusively within the rhabdosphincter and hence is more difficult to dissect, many prostate scRNA-seq samples have not captured the epithelial populations in this region. However, our meta-analysis detected PrU epithelial cells in several datasets [4, 7] as well as many in our proximal prostate scRNA-seq dataset (Fig. 1G, Additional file 1: Fig. S3G). Uniquely, PrU epithelial cells display hybrid luminal and basal features, similar to urothelial cells in the adjacent urethra [4]. However, PrU cells can be readily distinguished from urethral cells by lineage-tracing with an Nkx3.1-Cre driver [4].

To understand the unique morphological features of PrU cells, we imaged the periurethral region by electron microscopy and immunofluorescence staining (Fig. 2). At the ultrastructural level, PrU cells share some features with distal luminal (LumDist) cells, such as organelles involved in protein secretion, and many features with LumP cells, including a high density of mitochondria (Fig. 2A–F). Interestingly, several features of PrU cells also resemble urothelial cells of the urethra, including the nuclear orientation of more basally situated PrU cells, as well as the lumen-facing structures of apically situated cells, which may resemble the rigid, uroplakin-filled surface of urothelial cells. Thus, PrU cells share ultrastructural features of both the prostate and the urethral urothelium and may represent a physical transition between the two tissues.

Fig. 2

Imaging of mouse PrU cells reveals unique and shared features with prostatic and urethral cells. Scanning electron microscopy (EM) images of PrU cells show a focal region of cells where they appear to be multilayered (A), a region that is not multilayered and displays unique features (B), and a higher magnification of this region (C). The features of distal LumA cells (D), proximal LumP and basal cells (E), and LumP cells at higher magnification (F) are shown for comparison. Arrows indicate basal nuclear orientation (purple), mitochondrial density (red), apical membrane structures (orange), rough endoplasmic reticulum (green), and Golgi apparatus (blue). Scale bars in A–F indicate 5 µm. G–L Immunofluorescence staining show changes in basal and proximal keratin expression. G Overview of the periurethral region with neighboring urethral and proximal cells at low power. Insets show co-expression of basal keratins CK5 (red) and CK14 (green) in distal (H) and proximal (I) basal cells, and consistent CK5 but reduced CK14 in periurethral (J) and periurethral and urethral (K) basal cells. Proximal keratin CK4 (white) is maintained through the proximal and periurethral region (L). No superficial-like cells were observed in the periurethral region. Scale bars in G–L indicate 50 µm

At the level of gene expression, PrU cells uniquely express Lmo1, Anxa8, Dapl1, and Aqp3 and have higher Ly6d and Sca-1 expression than LumP cells [4] (Additional file 1: Fig. S5B). Although Krt5 and Krt14 expression overlaps in the basal layer throughout more distal regions of the prostate, Krt14 expression becomes intermittent in the PrU region and Krt5 is maintained, whereas basal cells of the urothelium rarely express Krt14 (Fig. 2G–L). Based on our re-analysis of a scRNA-seq dataset of the proximal prostate and urethra [7], we could define two distinct urethral populations, a luminal-intermediate urothelial cell group (which we term Urethral 2) with transcriptomic similarity with LumP cells and a basal-intermediate urothelial cell group (Urethral 1) with similarity with PrU cells (Additional file 1: Fig. S5). Notably, at homeostasis, PrU and LumP cells can be readily distinguished from urothelial cells by key markers, including several uroplakins (Additional file 1: Fig. S5). Thus, PrU cells also represent a transition population in terms of molecular features, such as gene expression.

The transcriptomic response to androgen deprivation and restoration

The prostate regresses in response to androgen-deprivation and regenerates after androgen restoration, which can be repeated through at least 30 cycles in the mouse [43, 44]. Following castration, the prostate undergoes rapid shrinkage and involution resulting in a stable regressed state, whereas restoration of androgen levels results in prostate regrowth to its former size [45]. To examine the response of individual cell populations to androgen-deprivation and restoration, we examined scRNA-seq data of mouse prostate through time courses of regression and regeneration [5, 6]. For this analysis, we defined a “cell type score” to represent the average of the most specific and differentially expressed genes for each cell type (“Methods”). In response to castration, every cell type except endothelial cells showed a significant decrease in its cell type score (Fig. 3A, B). Interestingly, the rates of transcriptomic change were different for each population, as distal luminal (LumDist) cells, myofibroblasts, and Mes1 cells rapidly lost almost all of their cell-type specific gene expression, whereas LumP, basal, smooth muscle, and Mes2 cells only lost approximately half of their cell-type specific gene expression, with Mes2 cells retaining their gene expression profile the longest (Fig. 3A, B).

Fig. 3

Time course of prostate regression and regeneration reveals androgen-dependent plasticity. A, B Meta-analyses of single-cell RNA-seq datasets for prostate regression and regeneration. A As described in [6], for regression time points, wild-type mice were castrated and prostate tissues from 2 biological replicates were collected at 1 day, 7 days, 14 days, and 28 days after castration; for regeneration time points, mice that had been castrated for at least 4 weeks were subcutaneously implanted with dihydrotestosterone (DHT) pellets, with 2 biological replicates collected at 1, 2, 3, 7, 14, and 28 days after pellet implantation. B As described in [5], prostate tissues were collected from wild-type mice at 7 and 28 days after castration. “Cell type score” is defined as the percentage of most specific differentially expressed genes for each population, averaged over the whole population (“Methods”). Changes in gene expression that are enriched in urethral but not PrU cells, such as Areg and Ociad2, in the LumA (C), basal (D), and LumP (E) populations, showing distinguishing genes for each population (left column), genes for general compartmental markers, and genes that are enriched for PrU and not co-expressed in LumP (right column), where the line indicates the average expression for each gene across the population and the bar indicates confidence interval (± 95%)

Interestingly, our analysis indicated that mouse prostate epithelial cells shift toward a PrU-like expression profile during regression. A detailed examination of gene expression patterns in LumA, basal, and LumP populations showed that each population lost expression of many specific genes but retained its distinctive expression of select distal luminal, basal, or proximal luminal markers during the regression-regeneration cycle (Fig. 3C–E). However, each epithelial population gained expression of multiple PrU markers following castration and lost this expression after androgen restoration; furthermore, the markers retained by LumP cells during regression were those that are co-expressed by PrU cells. The epithelial populations did not shift toward urethral gene expression profiles, as only rare LumP cells expressed any urothelial markers (Additional file 1: Fig. S6, Additional file 1: Fig. S7). Notably, while the normal PrU profile includes some genes that are co-expressed by either LumP or the urethral urothelium, the regressed epithelium expresses many PrU-specific genes that are distinct from both (Additional file 1: Fig. S6E). These findings highlight PrU-like transcriptomic profiles and provide a broader context for the previously reported shift from LumA toward LumP in the anterior prostate following androgen deprivation [6].

A transcriptomic shift was also observed in the prostate stroma during regression, as both the Mes1 and Mes2 fibroblast populations altered gene expression in response to androgen deprivation. Mes1 cells rapidly shifted toward a Mes2 expression profile and lost expression of several defining factors including Wnts, whereas Mes2 cells changed gene expression more slowly (Fig. 3A, B, Additional file 1: Fig. S6, Additional file 1: Fig. S7). Thus, we conclude that transcriptomic reprograming following androgen deprivation is not exclusive to the luminal or distal compartments, but instead represents a tissue-wide alteration of cell states.

Atlas of the human prostate

Next, we performed a meta-analysis of published scRNA-seq datasets to establish a corresponding reference atlas of the normal human prostate [4, 5, 11, 13] using the criteria described for the RMT pipeline (“Methods”). Despite differences in the relative proportions of cell populations between these datasets, the data were remarkably consistent. We found that the human prostate has a single basal epithelial population, two luminal populations corresponding to luminal acinar (LumAcinar) and luminal ductal (LumDuctal), and a periurethral-like (PrU) population (Fig. 4A–F, Additional file 1: Fig. S8A). The stromal populations were more variable and less well-represented across datasets, but corresponded to at least 1 endothelial population and 3 fibroblast-like populations (Fig. 4A, Additional file 1: Fig. S8A, B). Of the 3 fibroblast-like populations, the first expressed several classic fibroblast markers and did not subdivide readily (we denote these as general fibroblasts), the second corresponded to fibroblasts that express several muscle genes (myofibroblasts), and the third to fibroblast-like cells that express many contractile muscle genes (fibromyocytes) (Fig. 4F) [46]. Based on differential gene expression, we generated signatures for each epithelial and mesenchymal population (Additional file 4: Table S3). Within the immune compartment, we detected relatively fewer cells with variable representation of cell types between patients, so these populations were grouped as either myeloid or lymphoid. Interestingly, the zone of the prostate tissue did not have a clear effect on the transcriptome (Additional file 1: Fig. S8B, C), as previously reported [5].

Fig. 4

Reference plots for human prostate scRNA-seq data. A Aggregated composite UMAP plot for samples of benign human prostate and adjacent benign prostate. B–E Plots of individual datasets. B UMAP plot corresponding to primarily LumAcinar cells taken from the peripheral zone of 1 patient [4]. C Plot containing primarily basal and PrU cells from 1 patient with PCa [5]. D Dataset containing primarily basal, PrU, and LumDuctal cells from 1 patient without PCa [11]. E Dataset containing mixture of prostate and seminal vesicle, originating from 2 organ donor patients with no history of prostate disease [13]. F Dot plot of select top differentially expressed genes (among genes that are expressed in more than 60% of the population and have the highest mean expression in that population) for the epithelial and stromal cell populations from the reference aggregated normal human prostate. The lung club cell marker SCGB1A1 and hillock cell marker KRT13 are highlighted, indicating that these do not clearly correspond to single, distinct prostate cell types. G Heatmap comparing the total gene expression profiles of the cell types in the normal human prostate dataset [13] to those of the aggregated normal mouse prostate, using Wasserstein distance as a metric. Darker color indicates greater transcriptomic similarity. Tables listing the most similar mouse and human epithelial populations based on gene expression, generated by overlaying the mouse cell type signatures onto the human populations (H) and vice versa (I)

Since the nomenclature of human prostate epithelial populations differs between publications, we compared our previous nomenclature [4] to an alternative system that uses “Club” and “Hillock” lung terminology [11], using the Tabula Sapiens as a source of normal tissue (Fig. 4E, Additional file 1: Fig. S8D, E). Notably, we found that most of the “Club” cells corresponded to LumDuctal and PrU cells (Additional file 1: Fig. S8B, C), as they expressed common luminal genes and more specific markers like RARRES1, but did not consistently express the defining marker SCGB1A1 (Fig. 4F). Similarly, most “Hillock” cells corresponded to PrU cells (Additional file 1: Fig. S8B, C), as they expressed common luminal and basal genes as well as more specific markers such as KRT7, PSCA, RARRES1, LYPD3, and AQP3; moreover, expression of the Club- and Hillock-defining markers were not specific (Fig. 4). The remaining luminal cells corresponded to LumAcinar cells (Additional file 1: Fig. S8D, E), expressing common luminal cytokeratins as well as more specific markers including KLK3, MSMB, FOLH1, and TGM4 (Fig. 4F). These transcriptional similarities were separately confirmed by plotting the expression of each of these genes on the prostate single-nuclei RNA-seq data from the GTEx project portal [12]. Based on these analyses, we find that our descriptive nomenclature of human prostate epithelial populations correlates with lung terminology but appear to align more accurately with distinct cell types in the prostate.

To perform an updated cross-species comparison of cell type identities [4], we calculated the Wasserstein distance between gene expression profiles for each population in the aggregated mouse and human datasets in transcriptomic latent space (“Methods”) (Fig. 4G). While human and mouse basal cells have notably different profiles, human basal and PrU populations most closely resemble mouse PrU, human LumDuctal most closely resembles mouse LumP, and human LumAcinar most closely resembles mouse LumL followed by LumD (Additional file 1: Fig. S8D). To test the robustness of this analysis, we individually removed the top 20 differentially expressed genes from the mouse LumL expression profile and repeated the comparisons, which revealed that the greater similarity of human LumAcinar to LumL was primarily dependent on differential expression of a single gene, Msmb; otherwise, the transcriptomes of the different LumDist populations had similar marker overlap with human LumAcinar. Consequently, we suggest that human LumAcinar cells, which are distributed throughout different zones of the human prostate, may correspond more generally to mouse LumDist populations of all lobes. We additionally plotted the signatures of each population on the aggregated data of the other species to see how the differentiating genes versus the whole transcriptome compare across species (Fig. 4H, I). Together, these results suggest a clear correlation across species.

Distinguishing human prostate cancer progression by AR signaling levels

To examine alterations of the human prostate due to disease, we combined the normal prostate scRNA-seq datasets with those from patients with benign prostate hyperplasia (BPH) [7] and treatment-naïve prostate cancer [6, 14,15,16,