IntroductionThe ovarian architecture exhibits a dual-compartment organization comprising the parenchyma and stroma1,2. Each follicle, which serves as the fundamental functional unit of the ovarian parenchyma, consists of an oocyte enveloped by concentric somatic cell layers (granulosa and theca cells) that coordinate follicular expansion through paracrine/autocrine regulatory circuits3. These follicular structures reside within a stromal microenvironment and display pronounced spatial zonation, featuring dynamic cellular subtypes and region-specific extracellular matrix (ECM) remodeling4. Although the ovarian stroma dynamically participates in cellular crosstalk during follicular development and relies on structural and functional complexity, its functional zonation and regulatory mechanisms remain poorly resolved. Defining stromal specialization is critical for understanding reproductive competence and pathologies such as premature ovarian failure.The ovarian stroma constitutes a dynamic multicellular niche that integrates immune populations, vascular networks, neural projections, and molecularly heterogeneous stromal subtypes within an ECM-rich matrix5,6,7. Functioning as a follicular interface, this compartment orchestrates folliculogenesis through biomechanical constraints and paracrine signaling8. Experimental models have demonstrated that three-dimensional culture systems supplemented with defined ECM components9,10 and soluble factors (gonadotropins, FGF, and IGF) partially recapitulate oocyte maturation, yet persistent functional disparities persist between in vitro-generated and in vivo-developed follicles11. This functional divergence may reveal that perifollicular stromal cells (SCs) actively contribute to follicular development. Single-cell RNA sequencing (scRNA-seq) analyses further reveal spatially partitioned stromal subtypes with discrete developmental trajectories12,13, although the current models fail to explain how localized signaling coordinates follicular maturation, partly due to the limited integration of transcriptomic and architectural data.Through integrated single-cell transcriptomics and spatial profiling, we generated a multimodal developmental atlas spanning five postnatal stages (postnatal days (P) 3, 5, and 7 and weeks (W) 3 and 8) in murine ovaries14. Our analysis mapped four molecularly distinct stromal subtypes exhibiting discrete spatial zonation and functional specialization. Mechanistically, we identified a perifollicular SC subtype that orchestrates theca cell–granulosa cell (GC) coordination through midkine signaling (MK signaling) to sustain follicular progression. These findings establish perifollicular SCs as previously unrecognized functional components of the follicular developmental unit, which contributes to the definition of the functional ovarian niche and provides a framework for targeting stromal–GC interactions in individuals with fertility disorders.ResultsA single-cell map and spatiotemporal dynamics of gene expression during follicle development in the postnatal ovaries of miceWe established a multimodal cellular atlas by integrating single-cell transcriptomic profiles from five developmental stages (P3, P5, P7, W3, and W8) with spatial transcriptomic profiles from three developmental stages (P7, W3, and W8) to elucidate the spatiotemporal dynamic changes in postnatal folliculogenesis (Fig. 1a). This design displays follicular developmental trajectories from primordial to antral stages across the postnatal stage to sexual maturity while systematically mapping niche-specific cellular interactions and enabling the systematic mapping of cellular states across development and their microenvironmental dependencies.Fig. 1: Single-cell map of follicle development in the postnatal ovaries of mice.Full size imagea Temporal and spatial dynamic transcriptome sequencing strategy for mouse ovaries. Ovaries were isolated and disaggregated into single-cell suspensions; cells were barcoded and used for library construction and sequencing using the 10× Genomics platform, and the data produced after sequencing were analyzed with dedicated software (P3, P5, P7, W3, and W8). b UMAP plot of ovarian cells colored by 19 clusters. c UMAP plot of 10 ovarian cell populations. d UMAP plot of cellular markers of different cell populations shown as a function of density. e Percentages of the 10 ovarian cell populations at P3, P5, P7, W3, and W8. f, g UMAP plot and SpatialFeaturePlot of ovarian spots colored by cluster at 10× magnification. h Weights/strength of interactions of the 10 ovarian cell populations during the time course from P3 to P5, P7, W3, and W8.Our integrated single-cell atlas encompassed 58,319 ovarian cells and resolved 19 transcriptionally distinct clusters that were annotated based on canonical marker expression (Fig. 1b). A cluster-specific marker analysis (FindAllMarkers) further classified cells into nine established lineages and one undefined somatic population (Fig. 1c; Supplementary Table S1), with GCs and SCs dominating the cellular proportions (Fig. 1d, e). Spatial transcriptomic analyses (10× Visium) of samples collected at P7, W3, and W8 mapped scRNA-seq-identified lineages onto ovarian niches. The uniform manifold approximation and projection (UMAP) plot revealed nine spatiotemporal clusters (Fig. 1f), with spatial localization patterns revealing a conserved topological organization: endothelial cells, SCs, and GCs formed concentric layers around the oocytes (Fig. 1g; Supplementary Fig. S1a). This resource provides a systematic spatiotemporal annotation of ovarian cellular ecosystems.The quantitative assessment showed stage-specific variations in interaction frequency and intensity across developmental time points (P3, P5, P7, W3, and W8), with peak activity observed at P3, P7, and W8 (Supplementary Fig. S1b), according to the result of the CellChat analysis. Dynamic intercellular communication patterns during folliculogenesis revealed that SCs predominantly exhibit signaling activity (Fig. 1h). An analysis of different signaling roles indicated that SCs possess distinct functional and structural characteristics compared to other ovarian cell populations (Supplementary Fig. S1c), which suggested that SCs are the most important cell types involved in the transmission of biological signals in the ovary. These findings suggest that SCs function as functional orchestrators of follicular niche dynamics and are essential for coordinating developmental signaling with other cell types in the microenvironment of follicles.Identification and genetic dynamics of ovarian SC subtypesLeveraging CellChat-inferred signaling networks, we investigated ovarian SC heterogeneity across five developmental time points. Single-cell sequencing resolved eight transcriptomically distinct stromal subtypes (Fig. 2a, b), with spatially restricted clusters (e.g., Clusters 1, 2, and 6) exhibiting niche-specific distributions (Supplementary Fig. S2a). Cluster-defining markers (Star, Cyp17a1, Aldh1a2, Gstm2, Ptn, Col1a1, Enpep, and Tmem100) were systematically identified (Fig. 2c; Supplementary Fig. S2c and Table S2), and spatial mapping with the 10× Visium platform combined with orthogonal validation by in situ hybridization (ISH), immunohistochemistry (IHC), and immunofluorescence (IF) staining, confirmed the compartmentalized expression of Ptn mRNA, Enpep mRNA, STAR protein, and ALDH1A2 protein (Fig. 2d). Functional zonation emerged: Enpep mRNA marked perifollicular stromal niches correlated with follicular growth; STAR protein localized to both theca-adjacent and distal stromal domains; and ALDH1A2 protein dominated the cortical stroma and epithelial interfaces. These findings establish ovarian SCs as molecularly and spatially stratified subtypes that are classifiable through integrated single-cell and spatial omics analyses.Fig. 2: Genetic dynamics of ovarian SC subpopulations.Full size imagea UMAP plot of ovarian SCs colored by time at P3, P5, P7, W3, and W8. b UMAP plot of ovarian SCs colored by the 8 clusters. c Heatmap of the top 5 marker genes expressed in the SC clusters. d The locations of the cellular markers Ptn, Enpep, STAR, and ALDH1A2 in sections based on the results of the 10× Visium analysis and ISH, IHC, IF staining. e CytoTRACE analysis showing the reconstruction of the cellular differentiation trajectories of the 8 clusters. f Single-cell trajectories of the 5 SC states identified through a pseudotime analysis as a function of the developmental timeline. g Single-cell trajectories of SCs along the development time points (P3, P5, P7, W3, and W8). h Single-cell trajectories of SCs along the 8 clusters. i Expression of Ptn, Cyp17a1, Enpep, and Aldh1a2 along single-cell pseudotime trajectories. The red circles highlight cells with high expression of the genes on the pseudotime axis.We integrated the CytoTRACE-predicted differentiation potential with marker expression to determine the SC fate during folliculogenesis and identified Cluster 6 as putative less differentiated SCs (Fig. 2e). Pseudotemporal ordering via Monocle2 revealed five differentiation states and two bifurcation points (Fig. 2f, g). Early postnatal stages (P3) predominantly occupied root states, with stromal diversification initiating at P7 — a critical window coinciding with secondary follicle expansion (Fig. 2h). The trajectory analysis positioned Cluster 6 and a subtype of cells in Cluster 4 at differentiation origins. The stage-resolved expression patterns of Ptn, Cyp17a1, Enpep, and Aldh1a2 (Fig. 2i) further supported the partitioning of SCs into four functionally distinct subtypes defined by differentiation trajectories and the expression of effector gene modules.SC classification and functional analysis-based pseudotime trajectoriesIntegrating multimodal data, we defined four ovarian stromal subtypes — structural SCs (extracellular matrix organization), perifollicular SCs (proliferative niches), less differentiated SCs (differentiation potential), and steroidogenic SCs (hormone synthesis) (Fig. 3a), and UMAP plots of the results obtained at different postnatal stages (P3–W8) are presented in Fig. 3a and Supplementary Fig. S2b. Dynamic shifts in abundance revealed stage-dependent subtype predominance — the abundances of structural SCs and less differentiated SCs decreased postnatally, whereas the abundances of prefollicular SCs and steroidogenic SCs increased after P7 (Fig. 3b). Multimodal intersection analysis (MIA) validated the spatial context of the stromal classification (Fig. 3c). The Gene Ontology (GO) enrichment analysis of subtype-specific markers further resolved the functional divergence (Fig. 3d): Structural SCs were enriched for ECM organization and regulation of cell growth (Ptn15 and Col1a216,17,18); Perifollicular SCs for mitotic nuclear division and oxidative phosphorylation (Teme10019 and Enpep20); Steroidogenic SCs for steroid hormone biosynthetic process and aging (Star and Cyp17a121); Less differentiated SCs for pattern specification process and microvillus organization (Gstm222 and Aldh1b123).Fig. 3: Identification of ovarian SC subtypes.Full size imagea UMAP plot of 4 specific SC subtypes in mouse ovaries. b Percentages of the 4 SC subtypes at P3, P5, P7, W3 and W8. c MIA of the specific ovarian SC subtypes mapping to ovarian tissue. d Heatmap of the top 5 marker genes expressed in the 4 SC subtypes and GO analysis of enriched biological processes. e Spatial localization of 4 SC subtypes on slices analyzed using 10× Visium and an RCTD analysis. f UMAP plot of 6 ovarian cell populations in the human ovary scRNA-seq dataset (GSE255690); colors indicate the different populations. g UMAP plot of SCs colored by the 5 clusters identified in the human ovary scRNA-seq data. h–k UMAP plots of the cellular markers (COL1A2, CYP17A1, ALDH1A2, and ENPEP) of 4 SC subtypes in human ovaries shown as a function of density. Red circles highlight cells with high expression of the indicated genes in different SC clusters in human ovaries. l GO analysis of differentially expressed genes in perifollicular SCs in human ovaries.We mapped transcription factor (TF) regulon activity via SCENIC to elucidate the transcriptional regulatory mechanism underlying the functions of SC subtypes and identified dynamic TF networks orchestrating folliculogenesis. Leveraging 51 regulon activities (5481 target genes), we identified stage-resolved regulon activation patterns, with key developmental regulators (Maf, Wt1, Maff, Bclaf1, and Stat1) showing subtype- and stage-specific activity (Supplementary Fig. S2d). Robust Cell Type Decomposition (RCTD)-based spatial mapping demonstrated subtype-specific niche localization, with functional enrichment patterns mirroring their topological distributions (Fig. 3e). These spatially and functionally resolved stromal subtypes redefine ovarian niche organization, linking cellular identity to the regulation of folliculogenesis.We validated the evolutionary conservation of the ovarian SC classification by analyzing cynomolgus monkey (GSE130664)24 and human (GSE255690)2 ovarian datasets. Monkey SCs were resolved into four conserved clusters (Supplementary Fig. S2e, f), with the COL1A2+ and TMEM100+ subtypes indicating spatially distinct niches in the dimensionality-reduced clusters (Supplementary Fig. S2g, h). In human ovaries, integrated transcriptomics classified SCs into five clusters (Fig. 3f, g). Cross-species alignment using mouse stromal markers (COL1A2, CYP17A1, ALDH1A2, and ENPEP; Fig. 3h–k) revealed a conserved stromal zonation: Human ENPEP+ SCs (Cluster 1) exhibited functional congruence with murine perifollicular SCs and were enriched in the regulation of kinase activity (GO:0043549) and reactive oxygen species metabolic processes (GO:0072593) according to the results of the GO enrichment analysis, as shown in Fig. 3l. These results establish the biological similarity of stromal niche organization across different mammals, linking molecular signals to the regulation of folliculogenesis.Specific regulon networks of ovarian SCs during follicle developmentOur analysis revealed pervasive SC communication across all five developmental stages, with SCs engaging in frequent signaling crosstalk with neighboring cell subtypes. Stage-resolved ligand–receptor mapping revealed seven dominant stromal-derived pathways (PTN, MK15, noncanonical WNT25, PROS26, GAS27, ANGPTL28, and TWEAK29) that orchestrate folliculogenesis (Fig. 4a). MK signaling persisted throughout folliculogenesis, while PTN activity was restricted to pre-pubertal stages (P3–W3). A comparative pathway analysis revealed the following stage-specific regulatory modules (Supplementary Fig. S3a): P3 (nWNT, PROS, APELIN, GRN, and CALCR), P5 (IL-1, ANNEXIN, CCL, BMP, and FGF), P7 (TWEAK, APELIN, ANGPT, CALCR, and CXCL), W3 (NPR2, EGF, HH, BMP, and AMH), W8 (LIFR, FGF, KIT, SPP1, and MIF).Fig. 4: Stromal-derived MDK co-culture promotes follicle development.Full size imagea The outgoing signaling patterns of the 10 ovarian cell populations at P3, P5, P7, W3, and W8. b MK signaling pathway network among the 10 ovarian cell populations. c The contribution of each L−R pair to the MK signaling pathway. d Representative images showing the location and RNA expression levels of Mdk, Ncl and Lrp1 by ISH in ovarian sections from W8 mice. Arrowheads indicate positive signals in ovaries. e The intensity of signaling pathways related to follicular development among perifollicular SCs, GCs, and oocytes in the follicular microenvironment. f The weights/strength of interactions of the 6 ovarian cell populations in the human ovary. g The contribution of each L−R pair to the MK signaling pathway in the human ovary. h Violin plot showing the expression profiles of ligands and receptors in the MK pathway in six ovarian cell populations. i Representative images showing the location and RNA expression levels of MDK detected using ISH in a human ovarian section. Arrowheads indicate positive signals in ovaries.Notably, stromal-derived MK signaling spanned all developmental windows, suggesting its role as a constitutive niche maintenance signal. The spatiotemporal regulation of the PTN/MK (early) to TWEAK/ANGPT (mid) to NPR2/HH (late) pathways delineates a stage-specific stromal signaling code. These dynamic secretory programs position SCs as central signaling hubs that coordinate follicular developmental transitions through temporally constrained paracrine modules.MK signaling exhibited sustained activity throughout folliculogenesis (Fig. 4b), with SCs as primary midkine (MDK) protein sources and GCs as dominant targets. The analysis of receptor–ligand pairs identified Ncl, Lrp1, Sdc4, and Sdc1 as core components of the MK pathway (Fig. 4c). Spatial mapping of Mdk mRNA via ISH revealed that stromal-enriched expression peaked at P5, with perifollicular localization diminishing around the antral follicles and corpora lutea (Supplementary Fig. S3b), consistent with stage-specific roles in follicular maturation. In mature follicles (W8), Ncl and Lrp1 localized predominantly to GCs (Fig. 4d). Single-cell ligand–receptor mapping confirmed that Mdk is involved in the strongest SC–GC signaling axis, which is complemented by bidirectional MIF-mediated GC–SC crosstalk (Fig. 4e). The cross-species analysis revealed conserved MK pathway activity in human ovaries, although with expanded cellular sources — both GCs and SCs emitted MK signals (Fig. 4f–h). ISH using a human MDK probe showed developmental stage-dependent differences in localization between GCs and SCs, potentially reflecting species-specific niche adaptations. These results establish MK signaling as an evolutionarily conserved SC‒GC communication axis that is dynamically regulated across follicular transitions. While human inter-individual variability necessitates further validation, our multi-species framework positions MK signaling as a central regulator of follicular niche homeostasis.Spatiotemporal dynamics of GCs cooperate with SC subtypes to regulate folliculogenesisOur findings established SC-derived MK signaling as a critical regulator of follicular development through GC modulation. Given the spatiotemporal dynamics of SC subtypes, we asked whether GCs exhibit complementary heterogeneity to coordinate niche functions30. UMAP-based clustering resolved four GC subtypes across developmental stages: differentiated cumulus/mural GCs (Cyp19a1 and Slc38a3), less-differentiated GCs (Wt1 and Wnt6), proliferative GCs (Tpx2 and Mki67), and steroidogenic GCs (Cyp11a1 and Lhcgr) (Fig. 5a, b; Supplementary Table S3). Spatial transcriptomics confirmed the histological zonation of these subtypes, with less differentiated GCs localized to primordial follicle niches and steroidogenic GCs enriched in antral follicles (Fig. 5c). Proportional shifts revealed developmental priming: the abundance of less differentiated GCs progressively decreased from P3 to W3, whereas the abundance of steroidogenic GCs increased after P7 (Fig. 5d). Pseudotemporal ordering revealed progressive differentiation from less differentiated (trajectory root) to steroidogenic states (Fig. 5e, f), mirroring follicular maturation. These spatiotemporally resolved GC states synergize with stromal dynamics to coordinate folliculogenesis.Fig. 5: Genetic dynamic signatures of the GC lineage during follicle development.Full size imagea UMAP plot of GCs colored by analysis at P3, P5, P7, W3, and W8. b Dot plot of special cellular markers of GCs in mice. c SpatialFeaturePlot of ovarian spots colored by the expression of special cellular markers of GCs at 10× magnification. d UMAP plot of 4 GC subtypes and percentages of the 4 GC subtypes at P3, P5, P7, W3, and W8. e Single-cell trajectories of the GC subtypes as a function of the developmental timeline and cell state. f Single-cell trajectories of SC subtypes along with the 4 GC cell types. g The incoming signaling patterns of the 10 ovarian cell types at P3, P5, P7, W3, and W8. h Comparison of the cellular roles between GCs and 4 SC subtypes in the MK signaling pathway at P3, P5, P7, W3, and W8. The red dashed rectangles highlight the roles played by GCs and perifollicular SCs in the MK signaling pathway at P3, P5, P7, W3, and W8.We elucidated the cellular origin of MK signaling targeting GCs by prioritizing ligands and receptors across SC subtypes. MK emerged as the dominant incoming signaling pathway in GCs throughout folliculogenesis (Fig. 5g), with perifollicular SCs identified as the primary MDK protein source through a spatiotemporal signaling network analysis (Fig. 5h). These results establish perifollicular SC-derived MK signaling as a dynamically regulated paracrine axis that is essential for follicular niche maturation, where stage-specific MDK secretion might coordinate GC differentiation and proliferation while maintaining stromal ecosystem connectivity.Stromal-derived MDK co-culture promotes follicle developmentFunctional validation of perifollicular stromal-derived MDK through exogenous supplementation (1 mg/L) in alginate-encapsulated secondary follicles revealed stage-specific developmental regulation (Fig. 6a)31. MDK enhanced follicular growth kinetics on Days 6 and 8 (Day 6: control treatment, n = 113, 295.3 + 59.15 μm vs MDK treatment, n = 109, 315.4 + 49.67 μm, P = 0.0048; Day 8: control treatment, n = 113, 325.8 + 61.29 μm vs MDK treatment, n = 109, 347.1 + 58.41 μm, P = 0.0086), accelerating the increase in diameter despite the equivalent terminal size threshold (400 μm) (Fig. 6b). As shown in Fig. 6c, d, treated follicles displayed the activation of GCs through increased levels of proliferative indices (Ki67+, control treatment: n = 12 vs MDK treatment: n = 12, P = 0.0354), reduced levels of apoptotic markers (Cleaved C3, control treatment: n = 16 vs MDK treatment: n = 17, P = 0.0482), and increased gonadotropin responsiveness (FSHR+, control treatment: n = 18 vs MDK treatment: n = 18, P = 0.0464). This evidence suggests that MDK may become one of the critical paracrine effectors that coordinate follicular niche maturation. A spatial analysis of ENPEP+ stromal layers revealed MDK-dependent proliferative amplification by Day 4 (Supplementary Fig. S3c, d; control treatment: n = 46 vs MDK treatment: n = 71, P = 0.0418), demonstrating autonomous expansion mechanisms independent of follicular volume scaling. We also examined the levels of hormone response (Fshr and Lhcgr) and steroid hormone synthesis (Cyp11a1, Cyp19a1, and Hsd17b1) in GCs cultured in vitro. The results showed that the responsiveness of GCs to FSH and LH was increased by MDK supplementation, but the efficiency of steroid hormone synthesis was not significantly increased (Supplementary Fig. S3e; control treatment: n = 9 vs MDK treatment: n = 9), which may also be related to the heterogeneity of follicle culture in vitro. However, these peri-follicular stromal subtypes maintain structural homeostasis through mitotic synchronization, suggesting an intrinsic self-organizing program that dynamically calibrates the size of the stromal compartment to meet the stage-specific developmental demands. The MDK-mediated expansion of ENPEP+ cells indicates either autocrine regulatory circuits or follicle-derived signal potentiation underlying stromal niche plasticity.Fig. 6: Stromal-derived MDK co-culture promotes follicle development.Full size imagea Analysis of secondary follicles after culture in the alginate system in vitro with/without the addition of 1 mg/L MDK. b Analysis of the diameter of the follicles in 8-day cultures (**P