Modeling of senescence-related chemoresistance in ovarian cancer using data analysis and patient-derived organoids

Background Ovarian cancer (OC) is a malignant tumor associated with poor prognosis owing to its susceptibility to chemoresistance. Cellular senescence, an irreversible biological state, is intricately linked to chemoresistance in cancer treatment. We developed a senescence-related gene signature for prognostic prediction and evaluated personalized treatment in patients with OC. Methods We acquired the clinical and RNA-seq data of OC patients from The Cancer Genome Atlas and identified a senescence-related prognostic gene set through differential and cox regression analysis in distinct chemotherapy response groups. A prognostic senescence-related signature was developed and validated by OC patient-derived-organoids (PDOs). We leveraged gene set enrichment analysis (GSEA) and ESTIMATE to unravel the potential functions and immune landscape of the model. Moreover, we explored the correlation between risk scores and potential chemotherapeutic agents. After confirming the congruence between organoids and tumor tissues through immunohistochemistry, we measured the IC50 of cisplatin in PDOs using the ATP activity assay, categorized by resistance and sensitivity to the drug. We also investigated the expression patterns of model genes across different groups. Results We got 2740 differentially expressed genes between two chemotherapy response groups including 43 senescence-related genes. Model prognostic genes were yielded through univariate cox analysis, and multifactorial cox analysis. Our work culminated in a senescence-related prognostic model based on the expression of SGK1 and VEGFA. Simultaneously, we successfully constructed and propagated three OC PDOs for drug screening. PCR and WB from PDOs affirmed consistent expression trends as those of our model genes derived from comprehensive data analysis. Specifically, SGK1 exhibited heightened expression in cisplatin-resistant OC organoids, while VEGFA manifested elevated expression in the sensitive group (P<0.05). Intriguingly, GSEA results unveiled the enrichment of model genes in the PPAR signaling pathway, pivotal regulator in chemoresistance and tumorigenesis. This revelation prompted the identification of potential beneficial drugs for patients with a high-risk score, including gemcitabine, dabrafenib, epirubicin, oxaliplatin, olaparib, teniposide, ribociclib, topotecan, venetoclax. Conclusion Through the formulation of a senescence-related signature comprising SGK1 and VEGFA, we established a promising tool for prognosticating chemotherapy reactions, predicting outcomes, and steering therapeutic strategies. Patients with high VEGFA and low SGK1 expression levels exhibit heightened sensitivity to chemotherapy.


Introduction
Ovarian cancer (OC) is the primary contributor to gynecologic carcinoma worldwide.The first-line OC treatments, outlined by the NCCN guidelines, encompass comprehensive debulking surgery and platinum-based chemotherapy.In advanced stages, adjunctive antiangiogenic agents are recommended (1).While chemotherapy remains pivotal in OC treatment, about 70%-80% of patients experience relapse after treatment eventually culminating in chemotherapy resistance (2,3).
Cellular senescence is an irreversible biological state in which cells lose their ability to proliferate and transition from the cell cycle into a relatively stable state, an indispensable mechanism for tumor suppression (4).Chemotherapy-induced senescence (CIS), represents a subtype of cellular senescence, triggered by platinumbased chemotherapy, such as cisplatin (5).Eluding CIS might serve as a plausible explanation for drug resistance (6).It has been observed that CIS cells can break free from their arrested state, re-enter the cell cycle, and exhibit significantly elevated tumor initiation potential (7).In OC, both spontaneous and druginduced senescent cells contribute to cancer progression (8).Leveraging death receptor 5 (DR5)-selective agonists to augment treatment-induced apoptosis in senescent cancer cells may impede tumor progression (9).Consequently, prognosticating patients susceptible to senescence-induced chemotherapy resistance could hold the key to improving OC outcomes.
Organoids represent 3D multicellular structures akin to organs, enabling more precise replication of microenvironments for disease modeling, drug development, regenerative medicine, toxicology research, and personalized medicine.In 2009, Sato's team achieved successful cultivation of the first organoid in the small intestine that can undergo long-term proliferation and passage, marking the beginning of contemporary organoid research (10).In 2017, Soragni et al. successfully constructed OC organoids and revealed the potential therapeutic effects of ReACp53 through highthroughput screening (11).Phan, meanwhile, suggested the utility of high-throughput drug screening based on organoid technology in OC treatment (12).
In our study, we performed extensive data analysis to ascertain senescence-related differential expression genes in patients with OC, distinguishing between platinum-resistant and platinumsensitive groups.Subsequently, we formulated a predictive signature and validated its efficacy by assessing the expression of senescence genes in OC patients-derived organoids (OC PDOs).Further categorizing clinical samples into resistant and sensitive groups allowed us to scrutinize the cisplatin sensitivity and verify the prognostic potential of the signature, thereby paving the way for more personalized treatments.

Collection of data and tumor tissues
We initiated our study by sourcing clinical information and RNA-seq profiles of patients with OC from TCGA database (https://www.cancer.gov/).We retrieved senescence-related genes from the GeneCards database (https://www.genecards.org).Singlecell RNA-seq data were acquired from GSE154600 (https:// www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE154600).Tumor tissues were procured from patients with OC who underwent primary resection without prior chemotherapy, confirmed through pathology at the Department of Gynecology, Clinical Oncology School of Fujian Medical University, Fujian Cancer Hospital.The study secured approval from the ethical committee of Fujian Cancer Hospital (K2022-052-01).Informed consent was obtained from all participants, apprised of the research objectives.

DEGs of senescence in platinumresistant and sensitive groups
The analysis encompassed 2740 differentially-expressed genes (DEGs) derived from the examination of raw count RNA-seq data in 197 patients with platinum sensitivity and 90 patients with platinum resistance to OC from TCGA using R package "DEseq2".Employing Hallmark enrichment analysis, we scrutinized mechanisms or pathways with potential relevance to chemoresistance.After standardization of data by vst function in R package "DEseq2", 18 prognostic senescence-related DEGs were extracted through univariate Cox analysis, and protein-protein interaction (PPI) was evaluated using the STRING database (https://version-11-5.string-db.org/ ).

Construction of the risk model
Through multivariate Cox regression analysis, we identified three independent prognostic-related genes.The normality of the expression data was verified using Kolmogorov-Smirnov test and visualized by QQ plots.To test the significance of differences between the two groups, we performed Student's t-test for normally distributed data, and Wilcoxon rank sum test for non-normally distributed data.Subsequently, we constructed a prognosis model following the formula: risk score=SGK1*coef (SGK1)-VEGFA*coef (VEGFA).Patients were categorized into high-and low-risk groups based on the median score value, and the prognosis for these groups was scrutinized.We utilized the "survival" R package to execute Kaplan-Meier (K-M) analysis, probing for the survival differences between OC patients in the low-and high-risk groups.Based on the somatic mutation data from TCGA, we conducted gene mutation through "maftools" package.Single-cell RNA-seq data was integrated by "Seurat" and "SingleR" R packages.We normalized data with Log Normalize method, then utilized t-distributed stochastic neighbor embedding (t-SNE) via the "RunTSNE" function to cluster and visualize cell populations.Finally, we explored the expression of these two genes across various cell types.

Immune infiltration landscape and chemotherapy sensitivity analysis
To estimate the proportion of 22 immune cell types between the low-and high-risk groups, we employed the CIBERSORT algorithm (13).The R package "ESTIMATE" was employed to evaluate stromal, immune, and tumor purity scores of OC, based on the proportion of immune and stromal cells.We forecasted the chemotherapy response of commonly used drugs for OC patients through the Genomics of Drug Sensitivity in Cancer (GDSC; https://www.cancerrxgene.org/)database.The half-maximal inhibitory concentration (IC 50 ) was determined using the "pRRophetic" package to assess chemotherapy response (14).

Human ovarian cancer organoid establishment, culture, and verification
Specimens, sized 1-2 cm 3 , were minced into 1-mm 3 fragments and incubated in a tissue digestion solution (K601003; Tumor Tissue Digestion Solution; bioGenous ™ ) at 37°C for 20-40 minutes.Digestion duration was determined through microscopic observations of cells dissociating into 10-to 20-µm small clusters.The suspension was then filtered through a 100-µm nylon cell strainer and centrifuged at 300 g for 5 minutes.After lysing red blood cells using the Red Blood Cell Lysis Solution (E238010; bioGenous ™ ) and centrifuging again to wash the cell pellet with DMEM (Gibco), a portion of the suspension was mixed with organoid cryopreservation medium (E238023; bioGenous ™ ).The mixture was then placed in a programmed cooling box at -80°C for 24 hours and subsequently stored in liquid nitrogen.Approximately 2000-5000 cells per well were mixed with 35µL basement membrane extract (BME; bioGenous ™ ) and seeded into pre-warmed 24-well plates.Following BME solidification, each well was incubated with 500 µL human OC organoid medium (K2168-OC; bioGenous ™ ) at 37°C.The medium was replenished every 2-3 days.Capture images every 3 days with a microscope (TS2; Nikon) and a camera (SYA-C20).Organoids were passaged approximately every 14 days by dissociation with an organoid dissociation solution (E238001; bioGenous ™ ) for 10 minutes at 37°C.The success rate for establishing OC PDOs was 50%.After organoids were propagated to the first generation, immunohistochemistry was performed to validate whether their origin was consistent with tumor tissue.The process is outlined in Figure 1.

Statistical analysis
Statistical analyses were carried out using Perl scripts (v5.30.0),R software (v.4.1.0),and GraphPad Prism (v8.0.2).Kaplan-Meier (K-M) analysis was employed to assess survival differences.Differential functions were analyzed using the Wilcoxon ranksum test between two groups.All experiments were repeated in triplicate.The data were expressed as mean ± standard.P < 0.05 was considered as statistical significance.The flowchart for constructing and passaging ovarian cancer organoids.

Validation of senescence-related genes in the platinum-sensitive and resistant OC groups
Our investigation harnessed clinical data from the TCGA database, yielding 197 cases in the platinum-sensitive group and 90 cases in the drug-resistant group.Clinical pathological parameters were shown in Supplementary Table 1.Subsequently, RNA-seq data for these cases were acquired and subjected to analysis, resulting in the identification of 2740 DEGs (|log2FC|>0 and P<0.05), including 1373 upregulated and 1367 downregulated genes, which were then visualized through heat-maps and volcano plots (Figures 2A, C).We further extracted 43 aging-related DEGs by intersecting 279 aging-associated genes with the 2740 DEGs (Figure 2B).Employing univariate Cox analysis, we identified 18 genes exhibiting significant correlations with prognosis, including 9 genes with negative associations and 9 with positive correlations (Figure 2D).These 18 prognosis-linked genes were used to construct PPI networks, from which 12 pivotal genes formed a central network marked by interactions (Figure 2E).Functional enrichment analysis indicated their enrichment in pathways associated with hypoxia, inflammation, and glycolysis processes closely intertwined with cellular aging (Figure 2F).

Construction of a senescence-related gene signature in platinum-resistant and -sensitive OC groups
Through multivariate Cox regression analysis on these 12 core genes, we obtained 3 genes significantly link to prognosis, depicted in Figure 3A.Distinct survival ability of OC patients with different expression level of IFNG, SGK1 and VEGFA were discerned (P<0.05)(Supplementary Figure 1A) (Figure 3B).The normality of the expression data was verified using Kolmogorov-Smirnov test and visualized by QQ plots.The expression of SGK1 and VEGFA was verified to be normally distributed (Figure 3B), whereas IFNG was not (Supplementary Figure 1B).Student's t-test was used to exam the differences in expression levels between the platinumresistant and sensitive groups for SGK1 and VEGFA, and IFNG was using the Wilcoxon rank sum test (Supplementary Figure 1C).Notably, SGK1 and VEGFA demonstrated significantly different expression (P<0.05), while IFNG presented no significant expression differences between two groups (Figure 3C).Therefore, we included SGK1 and VEGFA in our risk model construction following the formula, risk score=SGK1*coef (SGK1)-VEGFA*coef (VEGFA) (coef values were in Supplementary Table 2).Segregation of patients based on risk scores yielded the high-and low-risk groups, with those in the high-risk category facing a worse prognosis (P<0.05, Figure 3D).Gene mutation waterfall plots for both groups are showcased in Figure 3E.In the GSE154600 dataset, there are 9 cell types marked mainly including cancer, immune, and stromal cells, and a t-SNE plot visualized cell subtypes in patients with OC (Figure 4A), with subtype annotations provided in Figure 4B.Both SGK1 and VEGFA exhibited increased expression in macrophages, suggesting a correlation between these risk genes and macrophage dysregulation in OC (Figures 4C, D).

GSEA and immune microenvironment status
Based on KEGG database, genes with higher expression levels in resistant group were involved in the peroxisome proliferatoractivated receptor (PPAR) signaling pathway, Cytokine Receptor Interaction pathway, and Neuroactive Ligand Receptor Interaction pathway (P<0.05), the same trend was found in Interferon Alpha/ Gamma Response pathway based on hallmark gene sets.While elevated involvement of oxidative phosphorylation pathway was found in genes highly expressed in sensitive group (P<0.01;Figures 5A, B).We carried out the CIBERSORT algorithm to analyze the proportion of the 22 immunocyte.Discrepancies in the microenvironment between high-and low-risk groups primarily centered around CD8 + T cell, monocytes, and macrophages (Figure 6A).Estimate analysis depicted correlations between immune, stromal and tumor purity scores with high-and lowrisk scores (Figure 6B). Figure 6C illustrated the heatmap of differences in infiltrating immune cells between the two groups.

Prediction of potential drugs
We further investigated the correlations between the IC 50 of potential target drugs and risk scores.The smaller IC 50 presented, the more sensitive patients were to drugs.Patients with high risk appeared to benefit from drugs such as gemcitabine, dabrafenib, epirubicin, oxaliplatin, olaparib, teniposide, ribociclib, topotecan, and venetoclax (Figure 7), while others seemed more beneficial for patients in the low-risk group.

Cultivation and verification in human ovarian cancer organoids
Our study successfully obtained and cultured three OC samples into organoids, visually documented through growth progression as indicated by the red arrows in Figure 8A.The development of organoids from P0 to P2 showed in Figure 8A (Scale bar = 100 µm).Furthermore, patients' clinical data are summarized in Table 1.Immunohistochemical testing of p53, WT-1, and MUC16 in the organoids, the routine marker in diagnosing OC, exhibited striking concordance with tissue expression, (Figure 8B, Scale bar for tissue = 70 µm, scale bar for organoid = 40 µm), affirming the homology between the organoids and the corresponding tumor tissue.We exposed organoids from the P1 generation to cisplatin diluted with medium across varying concentrations and captured their growth state after 72 hours.Figure 9A showed the representative growth status changes under 12.5 µM cisplatin medium in drugresistant and sensitive PDOs (Scale bar = 100 µm).Notably, the sensitive group's organoids exhibited shrinkage, rupture, and loss of original structure following cisplatin addition.By contrast, organoids in the resistant group displayed only mild shrinkage.The dose-response curve reflecting drug inhibition rates on the organoids is illustrated in Figure 9B.Organoids were categorized into cisplatin-resistant and -sensitive groups based on IC 50 and IC 90 .Specifically, PDO1121 was recognized as cisplatin-resistant sample, whereas PDO0213, and PDO0315 were designated as the cisplatin-sensitive ones.In these groups, VEGFA exhibited low expression in cisplatin-resistant organoids, while SGK1 was prominently expressed (Figure 9C). Figure 9D illustrated the protein expression and the gray scale analysis of SGK1 and VEGFA in OC PDOs.SGK1 were higher expressed in PDO1121 and VEGFA were higher in PDO0213 and PDO0315.

Discussion
Cellular senescence, instigated by DNA damage response, unfolds distinct characteristics: an irreversible growth halt; augmented lysosomal activity (16), amplified damage response signaling pathways; the emergence of macromolecular damage (17), and a senescence-associated secretory phenotype (SASP) (18,19).Senescence operates as a dual-edged sword.In normal tissue, it can impair repair and regeneration, accentuating aging,as well as maintain normal tissue homeostasis by immune-mediated clearance (20).In tumor tissue, it serves as a robust anti-tumor mechanism, thwarting the proliferation of cancerous cells and the genetic transfer of damaged cells on the one hand (21), but also promotes cellular reprogramming into a stem-like state, resulting in drug-resistant and invasive clones on the other hand (22).Contemporary studies spotlight chemotherapy and radiation as catalysts for treatment-induced senescence (TIS) within tumor cells (23).And CIS could interact with tumor micro environment (TME).In breast cancer it is confirmed CIS could change TME and increase the aggressiveness through the CXCL11 signaling pathway (24).Nicolas Adele et al. also revealed that the senescence of inflammatory cancer-associated fibroblasts (iCAFs) induced by chemoradiotherapy was tightly correlated with poor prognosis in rectal cancer (25).Meanwhile, cellular senescence might induce EMT and drug resistance through the PI3K/AKT pathway in colorectal cancer (26).Additionally, in lung cancer, TIS-induced cisplatin resistance may be pivotal to autophagy and hypoxia (27).
The study of OC manifests TIS within stromal cells and exhibits SASP through cancer-associated fibroblasts (CAFs), which underlies resistance of PARP inhibitors (28).Hence, the study of cellular senescence-related therapeutic targets may help to overcome tumor drug resistance.Taking death receptor 5 (DR5) as an example, the activation of DR5 may impede tumor progression, not only through the senescent cells but also nonsenescent cells adjacent by bystander effect (29).Furthermore, a novel nanoplatform has been studied targeting tumor necrosis Different proportion of 22 immune cells between high-and low-risk groups analyzed by CIBERSORT algorithm (A) and heatmap (C).(B) Estimate analysis of immune, stromal, and tumor purity score in high-and low-risk group.
factor-related apoptosis-inducing ligand (TRAIL) to further strengthen DR5-induced apoptosis, which provides a promising approach to clinically overcome tumor resistance (30).
In our investigation, a novel senescence-related gene signature sourced from patients with platinum-resistant and platinumsensitive OC, featuring genes such as SGK1 and VEGFA, was constructed with potential prognostic predictive value.The biological functions of these model genes are listed in Table 2. Patients were classified into low-and high-risk cohorts based on the median scores.Notably, patients in the low-risk group displayed superior prognoses, while those at higher risk had a significantly worse prognosis.SGK1 exhibited elevated expression, while VEGFA presented a low expression trend in the platinum-resistant group with a poor prognosis.We speculate that this may be related to poor angiogenesis in the platinum-resistant OC patients and the chemotherapeutic agents are unable to reach the tumor site.The low expression of VEGFA may also be strongly associated with the acquisition of resistance to bevacizumab leading to poorer outcomes in patients with platinum-resistant (39).Through the construction of OC PDOs, the credibility of this predictive signature Correlation between risk score and IC 50 of potential chemotherapeutics predicted by the "pRRophetic" package.was validated.Nonetheless, further inquiry remains imperative to elucidate the precise effect of these senescence-related genes on the prognosis of patients with OC.

B A
The GSEA unveiled significant enrichment of high expressed genes in chemotherapy-resistant cohorts in the PPAR signaling pathway.PPAR encompasses three isoforms: PPARa, PPARb/d, and PPARg-a nuclear receptor family pivotal in regulating and transcribing target genes, energy metabolism, cellular dynamics, inflammation, and carcinogenesis (42,43).The interplay between the PPAR pathway and chemoresistance has garnered attention across diverse cancers, including diffuse large B-cell lymphoma (DLBCL) (44), breast cancer (45), hepatocellular carcinoma (HCC) (46), non-small cell lung cancer (NSCLC) (47), and OC (48).In addition, activation of the PPAR gamma signaling pathway can balance the release of inflammatory and anti-inflammatory cytokines, orchestrating a pre-malignant microenvironment that promotes cell senescence, and alleviating tumor burdens (49).In endometrial carcinoma, inhibition of PARR through Bcl-2/Caspase3 pathways could hinder apoptosis, fostering carcinogenesis (50).In our investigation, the activation of the PPAR pathway was correlated with poor prognosis in OC patients.
To further validate the prognostic potential of our risk model, we constructed the OC PDOs.Organoids are intricate 3D multicellular structures that distinguish themselves from patientderived xenograft (PDX), which entail drawbacks such as high costs and extended durations, as well as 2D cell lines that suffer from limitations such as the absence of cell morphology, original tumor tissue structure, and intercellular interactions.Tumor organoids possess distinctive attributes: they adeptly mimic microenvironments, foster in vitro miniature tumor formations, and closely mirror the differentiation and expression characteristics of the original tumor tissue.This renders them highly suitable for investigating tissue responses to drugs, injuries, or mutations, conducting drug candidate screenings, and enabling precision clinical treatments (51).Phan et al. previously indicated the potential application of high-throughput drug screening based on organoid technology in OC treatment (12).In our study, we developed three human OC organoids, predicted patient responses to cisplatin treatment through ATP viability assay, and corroborated the expression of model genes (SGK1 and VEGFA) within cisplatin-resistant and cisplatin-sensitive groups.Notably, our findings aligned consistently with the analysis conducted on the TCGA database.
In this endeavor, a senescence-linked risk model was devised, underscoring the connection between model gene and PDO drug responses to cisplatin treatment.A broader sample set of clinical trials is necessary to adequately assess the clinical significance of our features.However, few patients have primary resistance in the clinical setting, while platinum-resistant patients tend to relapse six months after initial treatment, with fewer opportunities for surgery and difficult access to specimens.Due to the difficulty of obtaining enough platinum-resistant tissue to culture organoids, there are limitations to our organoid validation cohort.We plan to conduct further fundamental experiments to explore the complex mechanisms of senescence-related genes in OC's chemoresistance.

Conclusion
Through a holistic analysis and rigorous in vitro validation, our study conceived a senescence-related gene signature (SGK1 and VEGFA) in OC.Our model demonstrated the potential to forecast chemotherapy outcomes and prognosis, guiding therapeutic interventions.

Funding
The author(s) declare financial support was received for the research, authorship, and/or publication of this article.This work was sponsored by the Joint Funds for the Innovation of Science and Technology, Fujian Province (Grant number 2021Y9209).

2
FIGURE 2 Identification of 12 vital differentially expressed senescence genes in ovarian cancer.Heatmap (A) and volcano plot (C) of DEGs between the platinum-resistant and sensitive OC patients based on TCGA.(B) Venn diagram of ageing-genes and DEGs, (D) Univariate Cox analysis of the agingrelated DEGs expression.(E) PPI network illustrated the relationship among the 12 DEGs.(F) Hallmark enrichment analysis of the aging-related DEGs.(*P < 0.05).

3 5 GSEA 4
FIGURE 3 Construction and verification of the prognostic index.(A) Forest plot of aging-related DEGs via multivariate Cox regression analysis.(B) The survival ability of OC patients with high and low expression level of VEGFA (P=0.0168) and SGK1 (P=0.0137).Q-Q plot of VEGFA and SGK1 to test the normal distribution of data.(C) The expression of VEGFA and SGK1 in resistant and sensitive groups.(D) Differences in overall survival between high-and low-risk groups (P=0.0198).(E) Mutation waterfall maps show the gene mutation differences in high-and low-risk groups.(* P<0.05).
FIGURE 9 (A) The different state of organoids in two groups after culturing 72 hours with 12.5 mM cisplatin.The red arrows indicate organoid changes.(B) The inhibition ratio plot of organoids under different concentrations cisplatin.Data presented as mean ± SD. (C) Relative expression level of SGK1 and VEGFA in organoids.Data presented as the mean ± SEM, one-way ANOVA was used for statistical calculation, n=3 independent experiments.(D) Western blot and gray values of SGK1 and VEGFA expression in PDOs.(Scale bar = 100 µm, ns: No significance **P < 0.01.***P < 0.001, ****P < 0.0001).

TABLE 1
Clinical information of patients for culturing organoids.
(TC: Paclitaxel+ Carboplatin, * As of the time of paper writing, the patient has received 3 cycles of chemotherapy without CA125 normal).