A Comprehensive Multiomics Signature of Doxorubicin‐Induced Cellular Senescence in the Postmenopausal Human Ovary

De‐identified human ovarian tissue was obtained from the Northwestern University Reproductive Tissue Library (NU‐RTL) under Institutional Review Board approved protocols (STU00215770, STU00215938). Ovaries were obtained from females aged 50–70 years old (Table S1; Average 61.3 ± 6.2 years) undergoing bilateral salpingo‐oophorectomies and/or total laparoscopic hysterectomies (Table S1). Females with BRCA mutations, diagnoses of endometriosis, ovarian neoplasia, complex ovarian cysts, adnexal masses, or a history of breast cancer, radiotherapy, and chemotherapy were excluded. Upon collection, the tissue was divided into cross‐sections (3–5 mm thick) that were generated perpendicular to the long axis of the ovary (Figure 1b.i) (Devrukhkar et al. 2023). In the absence of significant gross pathology as assessed by a certified gynecologic pathologist, up to two ovarian cross‐sections were designated for research (Figure 1b.ii) and transported to the laboratory on ice in ORIGIO Handling IVF medium (Cooper Surgical Inc., Trumbull, CT, USA).

Ovarian tissue was kept in the ORIGIO Handling IVF medium for all the processing steps which were performed at room temperature (Anvari et al. 2023). The gross ovarian tissue slice was cut to remove the inner medulla, and a Stadie‐Riggs tissue slicer (Thomas Scientific, Chads Fort Township, PA, USA) was then used to generate 500 μm‐thick slices of tissue (Figure 1b.iii–vi). As confirmed by histology (Figure 1b.vii), the first slice was cortex‐enriched tissue, whereas the last slice was medulla‐enriched tissue. The cortex and medulla slices were processed separately and cut into 1 × 1 × 0.5 mm pieces with a scalpel (Figure 1b.viii). Ovarian tissue pieces were then transferred into cell culture inserts (Millipore Sigma, Burlington, MA, USA) with three pieces/insert for either the cortex (Figure 1b.ix with corresponding histology) or medulla (Figure 1b.x with corresponding histology). The inserts were then placed in 24‐well culture plates with each well containing 400 μL of pre‐equilibrated growth media. Growth media was made with MEM Alpha + GlutaMAX (Thermo Fisher Scientific, Waltham, MA, USA) supplemented with 20 mIU/mL recombinant Follicle Stimulating Hormone (Gonal F RFF Redi‐ject, Rockland, MA, USA), 1 mg/mL fetuin (Sigma‐Aldrich, St. Louis, MO, USA), 1 μg/mL Insulin‐Transferrin‐Sodium Selenite (Thermo Fisher Scientific, Waltham, MA, USA), and 3 mg/mL human serum albumin (Cooper Surgical Inc., Trumbull, CT, USA). This media was prepared with or without doxorubicin hydrochloride (Tocris Bioscience, Minneapolis, MN, USA) at concentrations of 0, 0.1, and 1 μg/mL depending on the experiment.

Explants were cultured at 37°C in a humidified atmosphere of 5% CO₂ according to the static culture paradigm (Figure 1a). Cortex and medulla explants from the same individual were used in parallel for both doxorubicin and control treatment in each experiment. In this manner, we minimized the impact of confounding variables due to participant heterogeneity. For doxorubicin treatment, explants were cultured in the compound for 24 h followed by culture in a doxorubicin‐free medium for up to 2, 5, or 9 additional days for short‐term and long‐term culture studies, respectively. Explants cultured without doxorubicin throughout the culture served as controls. For short‐term cultures, complete media changes were performed on Days 1 and 3. For long‐term cultures, complete media changes were performed on Day 1 at the end of doxorubicin treatment followed by half‐media changes every other day. Conditioned media was saved throughout culture for downstream analysis. Following the culture period, the explants were either fixed or snap frozen and used for downstream analyses as described. For the SASP analysis, additional steps were taken to ensure minimal background from protein supplements in the conditioned media (Anvari and Duncan 2023). On Day 10, a complete media change was performed with two washes with basal medium without protein supplements (MEM Alpha + GlutaMAX) in the original culture plate. The explants were separated from the inserts by removing the mesh bottom with sterile forceps, and then they were washed in basal medium. Explants were then transferred onto new inserts in plates with fresh basal medium and cultured for an additional 24 h, at which point the conditioned media was snap frozen for proteomic analysis.

Explants were fixed in Modified Davidsons fixative (mDF) (Electron Microscopy Sciences, Hatfield, PA, USA) at room temperature for 2 h and then overnight at 4°C. After overnight fixation, the explants were washed in and transferred to 70% ethanol and stored at 4°C until further processing. The explants were then dehydrated in an automated tissue processor (Leica Biosystems, Buffalo Grove, IL, USA), embedded in paraffin, and sectioned (5 μm thickness) with a microtome (Leica Biosystems, Buffalo Grove, IL, USA). To assess tissue architecture, standard hematoxylin and eosin (H&E) staining was performed using a Leica Autostainer XL (Leica Biosystems, Buffalo Grove, IL, USA). Tissue sections were cleared with xylene (Mercedes Scientific, Lakewood Ranch, FL, USA) in three 5‐min incubations and mounted with Cytoseal XYL (Epredia Thermo Fisher Scientific, Waltham, MA, USA).

IHC on ovarian explants was performed with the following antibodies: cleaved caspase‐3 (CC3), Ki67, p21^CIP1/WAF1, p16^INK4A, Superoxide Dismutase 2 (SOD2), nonmuscle myosin IIA (MYH9), Lumican, and Periostin (refer to Table S2 for all information including source of antibodies, dilutions, and concentrations) according to a previously established protocol by our laboratory (Machlin et al. 2021). Optimization of all antibodies was performed on native postmenopausal ovarian tissue containing both ovarian cortex and medulla along with positive controls (High Grade Serous Ovarian Carcinoma tissue sections) and nonimmune controls (Rabbit and Mouse IgG control antibodies, Vector Laboratories Inc., Burlingame, CA, USA) at the same concentrations as the corresponding primary antibody. In brief, slides were cleared using CitriSolv (Decon Labs, King of Prussia, PA, USA) and rehydrated in decreasing concentrations of ethanol. Antigen retrieval was performed by heat induced epitope retrieval (HIER) using Reveal Decloaker 10× at pH 6 (Biocare Medical, Pacheco, CA, USA) for all antibodies except Periostin (HIER at pH 9.0) and Lumican (no antigen retrieval). Slides were incubated in primary antibodies diluted in 3% Bovine serum albumin (BSA, Sigma‐Aldrich, St. Louis, MO, USA) in Tris‐buffered saline (TBS) at 4°C overnight, followed by incubation in secondary antibody (biotinylated goat anti‐rabbit, 1:200, or biotinylated goat anti‐mouse, 1:200, Vector Laboratories, Burlingame, CA, USA) for 1 h at room temperature. Antibody detection was performed using a 3′,3′‐diaminobenzidine (DAB) Peroxidase Substrate kit (Vector Laboratories, Burlingame, CA) that resulted in a brown precipitate. As soon as a brown staining was visible in the experimental section, the DAB reaction was stopped by transferring slides to distilled water. The slides were then counterstained with hematoxylin (Mercedes Scientific, Lakewood Ranch, FL, USA), cleared with CitriSolv, and mounted with Cytoseal XYL.

Ovarian explant samples were imaged on a RebelScope Imaging System (Discover ECHO Inc., San Diego, CA, USA) using a 40× objective with 200% optical zoom. Image analysis was performed using FIJI/ImageJ (ImageJ2 Version 2.14.0/1.54f, Madison, WI, USA) (Rueden et al. 2017). Images were processed with background subtraction followed by color deconvolution to split the images into Hematoxylin and DAB‐only images. These images were converted to binary images and then the following steps were performed: Watershed → Threshold (to highlight all nuclei in red) → Analyze particles (Size in pixels = 2500‐Infinity). For markers with nuclear localization (CC3, p21^CIP1/WAF1, and p16^INK4A) the percentage of positive cells was calculated by dividing the number of DAB‐positive cells by the total number of nuclei, and the data from three images per explant were averaged. For extracellular matrix (ECM) proteins or those with cytosolic localization (SOD2, MYH9, Lumican, Periostin), the percentage of DAB‐positive area per total area of each explant was reported.

Native ovarian tissue samples were scanned in brightfield with a 20× Plan Apo objective using the NanoZoomer Digital Pathology whole slide scanning system (HT‐9600) (Hamamatsu City, Japan) at the University of Washington Histology and Imaging Core. The digital image analysis (DIA) platform Visiopharm Integrator System (VIS); (Ver. 2023.01.1.13563) (Visiopharm, Hørsholm, Denmark) was used to analyze the IHC. Positive staining was detected by binary thresholding. The percent positive staining was calculated by determining the area of the positive stain label relative to the whole tissue section area.

Glucose levels in conditioned media were measured to assess glucose consumption by tissues as a measure of tissue viability (Elson et al. 2015; Prill et al. 2014). Conditioned media was thawed and vortexed, and glucose levels were measured using a GlucCell glucometer (CESCO Bioengineering, Marina Del Rey, CA, USA) according to the manufacturer's instructions. Baseline measurements on control and doxorubicin‐containing media were performed on media samples prior to culture. Glucose levels were then plotted as a mean value from all wells per condition for each participant.

SA‐β‐Gal activity was evaluated in frozen tissue sections. Ovarian tissue explants were first embedded in Tissue Tek optimum cutting temperature (OCT) compound (Sakura Finetek, VWR, Torrance, CA, USA) and frozen over a mixture of 2‐methylbutane (Fisher Scientific, Hampton, NH, USA) and dry ice. The frozen blocks were stored at −80°C until further use. Cryosections of 10 μm thickness were obtained through the Pathology Core Facility (Northwestern University). SA‐β‐Gal activity was assessed using a Senescence Detection Kit (Biovision, Milpitas, CA, USA) according to the manufacturer's instructions. In brief, the sections were allowed to equilibrate at room temperature for 1 min, fixed with fixative solution for 5 min, washed twice in 1× phosphate‐buffered saline (PBS) (Fisher Scientific, Hampton, NH, USA), and incubated overnight (15 h) in the SA‐β‐Gal‐staining solution at 37°C. The slides were then washed in PBS, counterstained with Nuclear Fast Red (Vector Laboratories, Newark, CA, USA) for 5 min, and mounted in aqueous mounting medium (CC/mount, Sigma‐Aldrich, St. Louis, MO, USA). The stained tissues were imaged under the RebelScope Imaging System using a 40× objective with 200% optical zoom. The presence of blue staining was considered indicative of SA‐β‐Gal positivity and was analyzed similar to IHC explants.

Human ovary explants were transported on dry ice from Northwestern University to the Buck Institute for Research on Aging and transferred to −80°C upon arrival. Explants were transferred into a prechilled 1.5 mL Eppendorf tube and immediately dissociated into single nuclei suspension using the 10× Genomics nuclei isolation kit (PN: 1000494) with a 15‐min lysis and two resuspension washes. A cordless motor pestle (VWR, Radnor, PA, USA) was used to dissociate the tissue into a single nuclei suspension. Final nuclei concentration was determined using the countess II automated cell counter (Thermo Fisher, Waltham, MA, USA) with propidium iodide (Invitrogen, Waltham, MA, USA). Single‐cell libraries were then prepared using the Chromium Next Gem Automated Single 5′ Library and Gel Bead Kit v.2 (PN100290) on a Chromium Connect robot (PN1000171), following the manufacturer's instructions. The cDNA and final gene expression libraries were quantified using a tape station (Agilent Technologies Inc., Santa Clara, CA, USA) and submitted for sequencing.

Data processing was performed using 10× Genomics Cell Ranger v6.1.2 pipeline. The "cellranger count" was used to perform transcriptome alignment, filtering, and UMI counting from the FASTQ (raw data) files. Alignment was done against the human genome GRCh38‐2020‐A. Cell numbers after processing were: 10‐day explants: doxorubicin‐treated cortex 7207 cells, control cortex 7980 cells, doxorubicin‐treated medulla 3774 cells, and control medulla 5324 cells; 6‐day explants: doxorubicin‐treated cortex 24,702 cells, control cortex 11,332 cells, doxorubicin‐treated medulla 22,852 cells, and control medulla 12,808 cells.

Downstream analyses were performed in R (version 4.2.0), primarily using the Seurat R package (version 4.1.1) (Hao et al. 2021; Satija et al. 2015) and custom analysis scripts. First, we executed a quality control step that removed the cells containing > 10% mitochondrial RNA and ≤ 250 genes/features. The doublet cells were identified and removed from the downstream analysis by using the DoubletFinder R package (version 2.0.3) (McGinnis et al. 2019) with parameters PCs = 1:30, pN = 0.25, and nExp = 7.5%. A total of 22,055 cells from 10‐day explant cultures and 71,694 cells from 6‐day explant cultures remained for subsequent analysis. Raw RNA counts were first normalized and stabilized with the SCTransform v2 function (SCT), then followed by the RPCA integration workflow for joint analysis of single‐cell datasets. In doing so, the top 3000 highly variable genes/features among the datasets were used to run SCT; and then 3000 highly variable genes/features and the 30 top principal components (PCs) with k.anchor = 5 were used to find "anchors" for integration. The clustering step was executed by using the 30 top PCs summarizing the RNA expression of each cell with a resolution parameter of 0.8.

To identify putative cell types, the singleR R package (version 2.0.0) (Aran et al. 2019) was used with the reference dataset of human primary cell atlas data (HPCA). Cell type annotation results from singleR were further refined by checking the manually curated marker gene list for the main cell types present in the human ovary (Data S1) (Fan et al. 2019; Lengyel et al. 2022; Wagner et al. 2020). The correlation and enrichment analyses of marker gene expression were conducted to assist in determining the cell types. Differential expression gene (DEG) analyses were done by functions in Seurat PrepSCTFindMarkers then FindAllMarkers/FindMarkers functions with the MAST algorithm (Finak et al. 2015). For the overall analyses and each cell type, the comparisons of DEGs of doxorubicin against control in either cortex or medulla ovarian tissue (Data S2 and S3) were calculated by using the FindMarker function with the parameter min.pct = 0.1 and logFC = 0.2. Rank‐Rank Hypergeometric Overlap (RRHO) analysis (Cahill et al. 2018; Plaisier et al. 2010) was performed by using the RRHO2 R package (version 1.0) to compare the differential expression patterns between doxorubicin and control of cortex versus medulla tissues. The ranks of the genes in the two gene lists were determined by calculating −log10(adj.pvalue)*Log2FC.

Following differential expression, ingenuity pathway analysis (IPA, Qiagen) was used to discover changes in enriched pathways in each comparison. DEGs with adjusted p values < 0.05 and |Log2FC| > 0.2 were incorporated into the IPA canonical pathway analysis. As the ovary explants were derived from three individuals across three distinct collections, where the donor effect was inherently mixed with the batch effect, we accounted for the donor effect in our DEG analysis to mitigate variability introduced by individual differences. Refer to the source code ov_snSeq_combined.R for the statistical analysis (Available on the GitHub webpage: https://github.com/FEI38750/Ovary_snRNAseq↗).

The degree of cellular senescence was quantified utilizing the AUCell R package (Aibar et al. 2017) which scores cell activity based on gene expression profiles. We applied multiple gene sets associated with cellular senescence, sourced from various peer‐reviewed databases and publications, as well as SASP factors from this study ("Buck ovary SASP Cortex and Medulla") (Data S4) (Ansere et al. 2021; Basisty et al. 2020; de Magalhães et al. 2024; de Magalhães et al. 2009; Gao et al. 2023; Kiss et al. 2020; Landry et al. 2022; Lengyel et al. 2022; Milacic et al. 2023; Saul et al. 2022; Shen et al. 2019; Thomas et al. 2022), to ensure comprehensive coverage and robustness of the senescence scoring. The two‐sided Student's t‐test was used to compare the difference in senescence scores between doxorubicin and control groups. Statistical significance was established at p value < 0.05.

To assess SASP factors in conditioned media by proteomics, 400 μL of conditioned media from human ovarian cortex and medulla explant tissue cultures in both doxorubicin treatment and control conditions was concentrated to ~30 μL with 0.5 mL 3 kDa filters (Millipore Sigma, Burlington, MA, USA). Aliquots of concentrated secretome (15 μL) for each sample were reduced using 20 mM dithiothreitol in 50 mM triethylammonium bicarbonate buffer (TEAB) at 50°C for 10 min, cooled to room temperature (RT) and held at RT for 10 min, and alkylated using 40 mM iodoacetamide in 50 mM TEAB at RT in the dark for 30 min. Samples were acidified with 12% phosphoric acid to obtain a final concentration of 1.2% phosphoric acid. S‐Trap buffer consisting of 90% methanol in 100 mM TEAB at pH ~7.1 was added, and samples were loaded onto the S‐Trap micro spin columns. The entire sample volume was spun through the S‐Trap micro spin columns at 4000 × g and RT, binding the proteins to the micro spin columns. Subsequently, S‐Trap micro spin columns were washed twice with S‐Trap buffer at 4000 × g at RT and placed into clean elution tubes. Samples were incubated for 1 h at 47°C with 2 μg of sequencing grade trypsin (Promega, San Luis Obispo, CA) dissolved in 50 mM TEAB. Afterwards, trypsin solution was added again at the same amount, and proteins were digested overnight at 37°C.

Proteolytic peptides were sequentially eluted from micro S‐Trap spin columns with 50 mM TEAB, 0.5% formic acid (FA) in water, and 50% acetonitrile (ACN) in 0.5% FA. After centrifugal evaporation, samples were resuspended in 0.2% FA in water and desalted with Oasis 10 mg Sorbent Cartridges (Waters, Milford, MA). The desalted protein lysates were then subjected to centrifugal evaporation and resuspended in 30 μL of 0.2% FA in water. Finally, indexed Retention Time standard peptides (iRT, Biognosys, Schlieren, Switzerland) (Escher et al. 2012) were spiked into the samples according to the manufacturer's instructions.

LC–MS/MS analyses were performed on a Dionex UltiMate 3000 system online connected to an Orbitrap Eclipse Tribrid mass spectrometer (both Thermo Fisher Scientific, San Jose, CA). The solvent system consisted of 2% ACN, 0.1% FA in water (solvent A) and 98% ACN, 0.1% FA in water (solvent B). Proteolytic peptides (2 μL of 1:40‐diluted sample) were loaded onto an Acclaim PepMap 100 C₁₈ trap column (75 μm × 20 mm, 3 μm particle size; Thermo Fisher Scientific) for 10 min at 5 μL/min with 100% solvent A. Peptides were eluted on an Acclaim PepMap 100 C₁₈ analytical column (75 μm × 50 cm, 3 μm particle size; Thermo Fisher Scientific) at 300 nL/min using the following gradient of solvent B: 2% for 10 min, linear from 2% to 20% in 125 min, linear from 20% to 32% in 40 min, up to 80% in 1 min, 80% for 9 min, and down to 2% in 1 min. The column was equilibrated with 2% of solvent B for 29 min, with a total gradient length of 215 min.

All samples were acquired in data‐independent acquisition (DIA) mode (Bruderer et al. 2017; Collins et al. 2017; Gillet et al. 2012). Full MS spectra were collected at a resolution of 120,000 (AGC target: 3e6 ions, maximum injection time: 60 ms, 350–1650 m/z), and MS2 spectra at a resolution of 30,000 (AGC target: 3e6 ions, maximum injection time: Auto, NCE: 27, fixed first mass 200 m/z). The isolation scheme consisted of 26 variable windows covering the 350–1650 m/z range with an overlap of 1 m/z (Data S5) (Bruderer et al. 2017).

DIA data were processed in Spectronaut (version 17.6.230428.55965) using directDIA. Data was searched against a human database containing all UniProt‐SwissProt entries extracted on 06/30/2023 (20,423 entries). Trypsin/P was set as the digestion enzyme and two missed cleavages were allowed. Cysteine carbamidomethylation was set as a fixed modification while methionine oxidation and protein N‐terminus acetylation were set as dynamic modifications. Data extraction parameters were set as dynamic and nonlinear iRT calibration with precision iRT was selected. Identification was performed using 1% precursor and protein q value. iRT profiling was selected. Quantification was based on the peak areas of extracted ion chromatograms (XICs) of the 3–6 best fragment ions per precursor ion, and q value sparse data filtering was applied. Interference correction was selected, and no normalization was applied. Differential protein abundance analysis was performed using an unpaired t‐test, and p values were corrected for multiple testing, using the Storey method (Burger 2018; Storey 2002). Protein groups were required to have at least two unique peptides. Protein groups with q value < 0.05 and absolute Log₂(foldchange) > 0.58 were considered significantly altered (Data S6).

Partial least square‐discriminant analysis (PLS‐DA) of the proteomics data was performed using the package mixOmics in R (version 4.0.2) (Rohart et al. 2017).

All graphs were generated using GraphPad Prism Software Version 9.4.1. Values were represented as mean ± SD. For immunohistochemistry analysis, normality was assessed by the Shapiro–Wilk test, and statistical significance was determined using unpaired t‐tests and one‐way ANOVA, with p values < 0.05 considered statistically significant. For transcriptomics data, the two‐sided Student's t‐test was used to compare the difference in senescence scores between doxorubicin and control groups, and statistical significance was established at p value < 0.05. For proteomics data, differential protein abundance analysis was performed using an unpaired t‐test, and p values were corrected for multiple testing using the Storey method (Burger 2018; Storey 2002). Protein groups with q‐value < 0.05 and absolute Log₂(foldchange) > 0.58 were considered significantly altered.

Postmenopausal human ovarian tissue contains numerous cell types, including fibroblasts, endothelial, epithelial, smooth muscle, and immune cells (Lengyel et al. 2022). Therefore, we developed an ovarian tissue explant culture model to interrogate doxorubicin‐induced senescence in an intact tissue that maintained this cellular heterogeneity. The postmenopausal ovary exhibits distinct architectural regions consisting of a cellular dense outer cortex and a vascular rich inner medulla (Figure 1b.vii). Therefore, we analyzed cortical and medullary tissues separately. We first determined an optimal dose that induced senescence but did not impact explant viability in culture. Ki67 and cleaved caspase‐3 (CC3) are biomarkers of cell proliferation and apoptosis, respectively, which are commonly used to evaluate cell/tissue viability, and the baseline expression of both is low in the postmenopausal ovary (Figure S1a,b). Since apoptosis is one of the main pathways involved in doxorubicin‐induced cell death in the ovary (Spears et al. 2019), we used it as a primary marker of viability in our explant model. First, we evaluated tissue morphology by H&E staining and apoptosis levels, quantifying CC3 expression in explants exposed to different concentrations of doxorubicin (0, 0.1, and 1 μg/mL) for 24 h followed by 2 days of culture in doxorubicin‐free control medium (Figure S2). Doxorubicin treatment did not impact the tissue morphology of ovarian cortical (Figure S2a) or medullary explants (Figure S2b) nor show any histological evidence of tissue necrosis (nuclear condensation, eosinophilic cytoplasm, pyknotic nuclei) when compared to control tissues (Laronda et al. 2014; Otala et al. 2002). The explant cultures also did not exhibit appreciable cell death on Day 1 or Day 3 of the culture relative to controls, suggesting our treatment paradigm did not affect explant viability or cause ovarian damage (Figure S2c–f). Ultimately, we selected a concentration of 0.1 μg/mL doxorubicin for subsequent long‐term cultures based on the fact that this dose has been used in other cell types to induce cellular senescence and the viability we established in our explant culture model (Alimirah et al. 2020; Demaria et al. 2017; Kitada et al. 2019).

For long‐term cultures, explants were cultured with or without doxorubicin (0 and 0.1 μg/mL) for 24 h, followed by an additional 5 and 9 days of culture in a doxorubicin‐free medium (Figure 1a). Doxorubicin did not affect gross tissue morphology or show histological evidence of overt tissue necrosis in both cortical and medullary explants as compared to Day 0 and control tissues at both 6 (Figure 2a) and 10 days (Figure 2b) of culture. Additionally, the explants displayed ovarian surface remodeling, a characteristic of wound healing and healthy stroma as seen by the smooth epithelial edge at the later stages of culture as compared to the rough edges at Day 0 (Figure 2a,b) (Jackson et al. 2009; Laronda et al. 2014). The level of apoptosis in cultured explants was also not significantly changed for both 6‐ and 10‐day cultures, suggesting that the cultured explants were viable (Figure 2c,e). We further confirmed tissue viability by assessing metabolic activity via measurement of glucose levels in conditioned media, which is commonly used as a noninvasive method to assess viability in long‐term cultures (Elson et al. 2015; Prill et al. 2014). Glucose levels in our cultures decreased throughout culture, indicating that the tissues were consuming glucose (Figure 2d,f). Taken together, we confirmed ovarian explant viability through a combination of histological, immunohistochemical, and metabolic assays.

To gain deeper and unbiased molecular insight into an ovarian senescence signature induced by doxorubicin treatment, we performed single nuclei RNA sequencing (snRNA‐seq) to assess the transcriptomic signature of cultured ovarian explants (Figure 1a). A total of 71,694 cells was profiled from 6‐day explants (cortex control = 11,332, cortex doxorubicin‐treated = 24,702, medulla control = 12,808, medulla doxorubicin‐treated = 22,852), and a total of 22,055 cells were profiled from 10‐day explants (cortex control = 7980, cortex doxorubicin‐treated = 7207, medulla control = 5324, medulla doxorubicin‐treated = 3744) as shown in the UMAP (Figure 3a). We utilized 11 senescence gene sets compiled from peer‐reviewed databases and publications to calculate a cellular senescence score based on gene expression fold change profiles across 6‐ and 10‐day cortex and medulla explants (Henceforth "senescence score" refers to the log2 fold change in senescence score) (Figure 3b). Based on this analysis, the 10‐day explants had consistently higher cellular senescence scores relative to the 6‐day explants in both the cortex and medulla (Figure 3b). Interestingly, the cortex had higher senescence scores than the medulla at Day 6 and Day 10, which were driven by distinct gene sets (cortex: 4, 11, 2, 1, 3, 6; medulla: 9, 5, 7, 10, 8) (Figure 3c and Figure S4a,b).

Given the time‐dependent increase in cellular senescence response, further analyses were restricted to the 10‐day explants. At a global level, doxorubicin treatment relative to the control resulted in 693 downregulated and 279 upregulated differentially expressed genes (DEGs) in the cortex, and 72 downregulated and 200 upregulated DEGs in the medulla (Figure 3d,e). Of the upregulated DEGs, 27 were shared between the cortex and medulla (10% of cortex DEGs and 14% of medulla DEGs) (Figure 3f, Data S2 and S3). Of the downregulated DEGs, 9 were shared between the cortex and medulla (1% of cortex DEGs and 13% of medulla DEGs) (Figure 3g, Data S2 and S3). We identified the top 20 upregulated and downregulated DEGs in the cortex (Figure 3h) and medulla (Figure 3i), showing a distinct senescence regional signature. However, we also identified the key upregulated and downregulated DEGs common between the cortex and medulla (Figure 3j). Importantly, cyclin‐dependent kinase inhibitor 1A (CDKN1A), better known as p21, a common marker of senescence, was equally upregulated in both regions (Figure 3j), consistent with the trend in IHC staining observed in Figure S3b,d. Based on ingenuity pathway analysis (IPA), doxorubicin altered genes in the cortex that are enriched in pathways associated with inflammation, fibrosis, oxidative stress, and immune responses, in particular IL‐6 and HMGB1 signaling, which are common hallmarks of senescence (Figure S5a). However, there was a reversal in activated pathways between the two regions, such as the unfolded protein response and the immunogenic cell death signaling pathways (Figure S5a,b). These were among the highest activated pathways in the cortex but the most inactivated pathways in the medulla (Figure S5a,b). Collectively, these results implicate a compartment‐specific response to senescence induction in the human ovary.

A major strength of our model is the ability to induce cellular senescence within the context of an intact tissue with cellular heterogeneity. To further understand the senescence response at a cellular level, we analyzed cellular composition in Day 10 cultured explants using the Uniform Manifold Approximation and Projection (UMAP) algorithm to plot the eight identified transcriptionally distinct cell clusters based on cell type‐specific markers (Figure 4a and Data S1). Processing the cortex and medulla separately in our paradigm allowed us to observe proportional cell type differences between these compartments. Except for epithelial cells in the cortex and immune cells in the doxorubicin‐treated explants, the proportion of each cell type was similar in each region and experimental condition (Figure 4b). Two distinct stromal cell populations were observed: Stroma 1 and Stroma 2 (Figure 4b). These populations appear to be biologically meaningful given that they were observed in all individuals (Figure S6a). Based on analysis of DEGs and gene ontology, the top genes for the Stroma 1 population were involved in lipid metabolism and cell‐matrix adhesion, suggesting stromal cells with some remnant steroidogenic activity, whereas the Stroma 2 population appears to be matrix fibroblasts (Figure S6b,c and Data S7). However, more studies are needed to fully understand the differences between these populations and their biological relevance. As expected, given the presence of the ovarian surface epithelium, epithelial cells were greatly enriched in the control cortex relative to the medulla. This cell type, however, appeared to be sensitive to doxorubicin and showed a decrease following treatment (Figure 4b). Additionally, immune cells were relatively higher in doxorubicin‐treated cortex and medulla explants. However, these trends in cell composition changes were not significant due to participant heterogeneity (Figure S6d).

A senescence score was calculated to determine the enrichment of senescence genes across the eight distinct cell clusters in the cortex and medulla, permitting us to localize potential cell types driving the senescence signature we identified in the cortex and medulla (Figure 3). Our results revealed that epithelial and stromal cells (Stroma 1) in the cortex (Figure 4c and Figure S4c,d) and stromal cells (Stroma 1) in the medulla (Figure 4d and Figure S4e) exhibited the highest senescence scores when compared to the overall tissue score.

To further characterize the cell types that may be driving the senescence signatures, we examined the common DEGs shared between the cortex and medulla and the epithelial and stromal cells. Although cortical epithelial cells had the highest senescence score in the cortex, only 15 upregulated genes (Figure S7a) and 10 downregulated genes (Figure S7b) were common with the overall cortex DEGs. Nonetheless, 3 DEGs out of the top 20 upregulated cortex genes are shared with the epithelial cells, PLAT, CFB, and SOD2 (Figure 3h and Figure S7c). Whereas 5 DEGs out of the top 20 downregulated cortex genes are shared with the epithelial cells, A2093722.1, A2024230.1, THSD4, EYA2, and CCSER1 (Figure 3h and Figure S7d). However, cortex stromal cells (Stroma 1) contributed more shared DEGs with the overall cortex DEGs. Of the upregulated genes, 88 were shared (32% of overall cortex DEGs and 53% of stromal DEGs) (Figure 4e), whereas 12/20 of the top 20 upregulated cortex DEGs were shared, including TIMP1, SOD2, and SERPINE 1 (Figures 3h and 4i). There were 54 shared downregulated genes between cortex and cortex stroma (8% of overall cortex DEGs and 43% of stromal DEGs) (Figure 4f,i). In the medulla, there was a greater overlap of DEGs with medulla stromal cells (Stroma 1). There were 51 shared upregulated genes (26% of overall medulla DEGs and 65% of stromal DEGs) (Figure 4g), whereby 9/20 of the top 20 upregulated medulla DEGs were shared, including SLPI, PDE10A, and SPON1 (Figures 3i and 4j). Between the medulla and medulla stroma, 20 downregulated genes were shared (28% of overall medulla DEGs and 49% of stromal DEGs) (Figure 4h,j). These results suggest that epithelial and stromal cells (Stroma 1) in the cortex and stromal cells (Stroma 1) in the medulla might be driving the senescence signature in these regions.

IPA provided an overview of enriched pathways in the stromal cells of the cortex and medulla revealing contrasting activated pathways (Figures ). The top significantly activated pathway in the cortex stroma was "Integrin signaling," which regulates cell–cell and cell‐extracellular matrix adhesion, including cellular proliferation/migration and activation/release of cytokines (Figure ). The top significantly activated pathway in the medulla stroma was "EIF2 signaling," which regulates mRNA translation and controls proteome expression of downstream pathways such as the Integrated Stress Response (ISR) and inflammatory production of cytokines (Figure ). These results unveil the cell type‐specific heterogeneity in response to doxorubicin‐induced senescence in human ovaries. S8a,b S8a S8b

Senescent cells are characterized by the production and release of a prominent secretome known as the senescence‐associated secretory phenotype (SASP) (Di Micco et al. 2021; Hernandez‐Segura et al. 2018; Herranz and Gil 2018; Kumari and Jat 2021; Wiley and Campisi 2021). We were able to profile the SASP in our model by assessing the conditioned media using proteomics (Figure 5a, Table S1). On Day 10 of culture, tissues were thoroughly washed with serum‐free basal media (alpha‐MEM) and transferred to inserts in a clean plate with pre‐equilibrated serum‐free basal media (Figure 5a). The conditioned media were then collected after 24 h for proteomic analysis (6 replicate media wells/condition), except "medulla control" that had media from four replicate wells (400 μL of conditioned media/well). Proteins from the conditioned media were concentrated and prepared for proteomic analysis as described using directDIA to identify SASP factors (Figure 5a).

Supervised clustering using partial least squares‐discriminant analysis (PLS‐DA) performed with all 314 quantified protein groups (with ≥ 2 unique peptides) revealed treatment‐based grouping where the SASP from control explants ("Control") clustered separately from SASP from doxorubicin‐treated explants ("Doxo") for both cortex and medulla (Figure 5b,e). In the cortex, 164 protein groups (135 upregulated and 29 downregulated) and in the medulla, 217 protein groups (184 upregulated and 33 downregulated) were significantly altered with doxorubicin exposure as compared to controls, including Basement membrane‐specific heparan sulfate proteoglycan core protein (HSPG2), Myosin‐9 (MYH9), and Periostin (POSTN) (q value < 0.05 and absolute log₂(foldchange) > 0.58) as highlighted in the volcano plots (Figure 5c,f). These included factors such as Decorin (DCN), Lumican (LUM), Clusterin (CLU) in the cortex, HSPG2, mitochondrial Stress‐70 protein (HSPA9) in the medulla, and POSTN, MYH9, Serpin H1 (SERPINH1), and Tenascin (TNC) in both cortex and medulla (Figure 5c,f). An over‐representation analysis performed using ConsensusPathDB revealed Gene Ontology (GO) biological processes (BP) upregulated in the cortex including the regulation of coagulation and response to wounding, and in the medulla including protein‐DNA complexes and chromatin assembly (q‐value < 0.05 and term level ≥ 4) (Figure 5d,g). Between the 135 proteins upregulated in the cortex and 184 proteins upregulated in the medulla, 120 proteins overlapped across both compartments with doxorubicin exposure (Figure 5h). The top five upregulated proteins were serpin H1, tenascin, periostin, apolipoprotein A1 (APOA1), and superoxide dismutase 2 (SOD2), several of which are known to be involved in cellular senescence and aging. For example, Periostin and Tenascin‐C were revealed to be crucial players in ovarian aging (Dipali et al. 2023). Additionally, Serpin H1 and Periostin appear to be biomarker candidates for chronic inflammation‐associated cancers (Bons et al. 2023; González‐González and Alonso 2018). Several of these proteins are part of the established core SASP, including Periostin, Tenascin, SerpinH1, HSPG2, and COL12A1 (Basisty et al. 2020). However, our study also revealed potential novel SASP factors specific to the postmenopausal human ovary, including hemopexin, complement C5, small leucine‐rich proteins (SLRPs), nidogens, and fibulin‐1 (Data S6).

To obtain an integrated doxorubicin‐induced ovarian senescence senotype in the human postmenopausal ovary, we identified the upregulated transcripts and secreted proteins that overlapped between the tissue transcriptome and conditioned media proteome for both cortex and medulla explants (doxorubicin vs. control). A total of 26 unique markers overlapped between the transcriptome and the SASP, and several of these genes/proteins have been implicated in cellular senescence, aging, and/or ovarian cancer (Table 1, Figure 6a). The senotype included regulated transcripts and secreted proteins involved in extracellular matrix (ECM) organization and remodeling (LUM, LRP1, LAMA4), reactive oxygen species (SOD2, PRDX1, P4HB), metabolic pathway regulators (NAMPT, NT5E, PKM, ALDOA), cytoskeletal and actin‐binding proteins (MYH9, ACTB, GSN), the protein disulfide isomerase (PDI) family of ER proteins (PDIA3, PDIA4, PDIA6), heat shock proteins (HSP90AB1, HSP90B1, HSPA5, HSPA8, HSPB1), and others (TGM2, UBC, VCP, ANXA2, CFI) (Table 1). Unbiased IPA revealed enrichment of genes in pathways including "unfolded protein response," "protein ubiquitination," and "cellular response to heat stress," and this is particularly interesting given the known involvement of these pathways in cellular senescence (Abbadie and Pluquet 2020; Hamazaki and Murata 2024) (Figure 6b).

To determine the physiologic relevance of these findings, we mapped a subset of these proteins back onto native postmenopausal ovarian tissue (Figure 6c and Figure S9). Among these were lumican (LUM), superoxide dismutase (SOD2), and nonmuscle myosin IIA, all of which showed the highest overlap between the transcriptome and the proteome. Additionally, we mapped the SASP factor Periostin (PSTN) which was among the top five upregulated proteins in the cortex and medulla SASP along with SOD2 (Figure 5). All these markers were expressed in a subset of cells throughout the ovarian cortex and medulla (Figure 6c and Figure S9). Lumican was most widely expressed (~10%–25% of tissue area), whereas periostin showed the least expression (~0.2%–2% of total tissue area) (Figure 6c; Figures S1e,h, S9a,d). SOD2 was expressed in both the cortex and medulla (3%–5% of total tissue area) with particularly strong expression in the ovarian surface epithelium (Figure 6c; Figures S1f and S9b). MYH9 was highly expressed in the vessel walls across both ovarian compartments (4%–7% of total tissue area) (Figure 6c; Figures S1g and S9c). It is important to note that not all cells that stain positive for individual markers are likely to be senescent, but it may instead be regions of overlapping expression that define senescent cells (e.g., white box in cortex and black box in medulla, Figure 6c). Future studies to map and quantify these factors across a broader aging series of human ovaries will provide novel insight into how the senotype changes across age, will aid in identifying senescence biomarkers specific to the ovary, and delineate the potential role of senescent cells in ovarian aging.