Altered tumor signature and T-cell profile after chemotherapy reveal new therapeutic opportunities in high-grade serous ovarian carcinoma

Abstract

Chemotherapy combined with debulking surgery is the standard treatment protocol for high-grade serous ovarian carcinoma (HGSOC). Nonetheless, a significant number of patients encounter relapse due to the development of chemotherapy resistance. To better understand and address this resistance, we conducted a comprehensive study investigating the transcriptional alterations at the single-cell resolution in tissue samples from patients with HGSOC, using single-cell RNA sequencing and T-cell receptor sequencing techniques. Our analyses unveiled notable changes in the tumor signatures after chemotherapy, including those associated with epithelial–mesenchymal transition and cell cycle arrest. Within the immune compartment, we observed alterations in the T-cell profiles, characterized by naïve or pre-exhausted populations following chemotherapy. This phenotypic change was further supported by the examination of adjoining T-cell receptor clonotypes in paired longitudinal samples. These findings underscore the profound impact of chemotherapy on reshaping the tumor landscape and the immune microenvironment. This knowledge may provide clues for the development of future therapeutic strategies to combat treatment resistance in HGSOC.

Abbreviations
• BCR
• B-cell receptor
• CNVs
• copy number variations
• DEGs
• differentially expressed genes
• EMT
• epithelial–mesenchymal transition
• HGSOC
• high-grade serous ovarian carcinoma
• NACT
• neoadjuvant chemotherapy
• PCA
• principal component analysis
• scRNA-seq
• single-cell RNA sequencing
• TCR
• T-cell receptor
• UMAP
• Uniform Manifold Approximation and Projection

INTRODUCTION

High-grade serous ovarian cancer (HGSOC) is one of the most lethal cancers due to challenges in early diagnosis, leading to advanced-stage disease in most patients.^1,2 Currently, neoadjuvant chemotherapy (NACT) is universally accepted as the standard treatment for patients with HGSOC. However, approximately 70% of patients experience relapses with drug resistance within 5 years.³ As second-line treatment, anti-angiogenic agents and poly(ADP-ribose) polymerase inhibitors are utilized, while immunotherapy has not demonstrated the same degree of effectiveness witnessed in non-small cell lung cancer or prostate, kidney, and urothelial cancer.^4–6 Therefore, understanding the repercussions of post-chemotherapy on both tumor cells and the surrounding microenvironment is essential to discover novel strategies for treatment-resistant HGSOC.

Chemotherapy is known to elicit an anticancer response by directly eliminating tumor cells and immune activation. Although some studies have attempted to measure the effects of chemotherapy on the immune microenvironment,^7–9 there is a lack of comprehensive cellular-level studies capable of fully elucidating the complex tumor microenvironment.

In this study, we characterized post-chemotherapeutic tumors and their surrounding microenvironments, focusing on transcriptional alterations at single-cell resolution. Single-cell RNA sequencing (scRNA-seq) combined with T-cell receptor sequencing (TCR-seq) for six treatment-naïve and four post-treatment HGSOCs revealed the activation of epithelial–mesenchymal transition (EMT) and repair programs in the tumor and stroma. Additionally, we identified alterations in immune phenotype toward less suppressive states. These alterations provide new therapeutic targets and windows for potentiating the efficacy of cytotoxic agents and immunotherapies in patients with chemotherapy-resistant HGSOC.

MATERIALS AND METHODS Patient information and sample preparation process

In this study, we included HGSCs of ovary, fallopian tube, and peritoneum, because HGSOC now refers to the high-grade serous carcinoma (HGSC) of pelvis arising from the ovary, fallopian tube, or peritoneum, which are considered to be originated from distal fallopian tube epithelium¹⁰ (Table S1). This study was approved by CHA University, CHA Bundang Medical Center (IRB No. 2019-08-039). All patients provided signed informed consent for specimen collection. A total of ten fresh tumor tissue samples were obtained from nine patients, with a paired sample from one patient. Briefly, tumor tissues were minced and dissociated using the enzyme mix sourced from the Tumor Dissociation Kit, human (130-095-929, Miltenyi Biotec) with the h_tumor_02 program twice on a gentle MACS dissociator/C tube (130-093-237, Miltenyi Biotec). At the end of each program, the samples were incubated in a prewarmed 37°C incubator for 30 min. The reaction was halted with a medium containing penicillin–streptomycin and fetal bovine serum (FBS). The sample was then filtered through a 70-μm strainer filter (352350, Life Technologies) and processed with Ficoll-Paque PLUS (17-1440-02, GE Healthcare) to isolate mononuclear cells. The resultant cell pellets were washed, resuspended in a cell-freezing medium, and stored in liquid nitrogen.

scRNA-seq

After post-thaw washing, single-cell suspensions were subjected to 5′ scRNA-seq using Single Cell 5′ Library and Gel Bead Kit V1, Single Cell A Chip Kit, and i7 Multiplex Kit (10× Genomic) according to the manufacturer's protocol. The library for T92 cells was constructed using Chromium Next GEM Single Cell 5′ Library & Gel Bead Kit v1.1, and Chromium Next GEM Chip G Single Cell Kit (10× Genomics). Each target for cell recovery was set at 5000–7000 cells. TCR-enriched and B-cell receptor (BCR)-enriched libraries were prepared using the Single Cell V(D)J Enrichment Kit, Human T-cell, and Single Cell V(D)J Enrichment Kit, Human B Cell (10× Genomics) according to the manufacturer's protocol.

Subsequently, the libraries were sequenced using Illumina NovaSeq 6000 (T04, T34, T39, T66, and T70), HiSeq X (T107 and T110), or HiSeq 2500 (T17, T18, T92, TCR, and BCR libraries). Sequencing reads were mapped to the GRCh38 human reference genome using the Cell Ranger software (version 5.0.0). The sequencing information for each sample is summarized in Table S2. Notably, scRNA-seq data generated from the nine samples were previously analyzed,¹¹ while one paired sample, T92, was newly included in this study.

scRNA-seq data processing

Cell Ranger outputs were refined for each sample using R package SoupX to eliminate ambient RNA.¹² Except for the batch correction function, all data processing steps were performed using the Seurat R package (Version 4).¹³ Low-quality cells were filtered out based on the counts per cell, number of genes per cell, and mitochondrial gene ratio per cell (Figure S1A). The filtered cells were then normalized using default parameter values, and gene expression was scaled to z-scores. During the scaling process, regression was applied to cell cycle score and percentage of mitochondrial RNA. Principal component analysis (PCA) was performed using the RunPCA function. For batch correction and clustering, BBKNN algorithm was employed using the RunBBKNN function from the bbknnR package¹⁴ along with the FindClusters function. The principal component (PC)/resolution options used for tumor cell and T-cell clustering were 10/1.0 and 13/1.5, respectively.

Gene set enrichment analysis for tumor cells

Hallmark pathway gene sets were obtained from the Molecular Signatures Database (MSigDB),¹⁵ and cell-level enrichment scoring was performed using the VISION R package.¹⁶ Sample level gene scoring was conducted using the gene set variation analysis (GSVA) R package with average gene expression values.¹⁷ The heatmap showed the sample-level GSVA results, and after performing hierarchical clustering, the pathways were divided into four distinct groups.

Differentially expressed genes (DEGs) analysis

We used the FindMarkers function of the Seurat package to identify DEGs between the pre- and post-chemotherapy groups. DEGs were filtered based on adjusted p-value <0.05 and abs(log₂foldchange) > 1. We then performed GSEA using the fgsea R package.¹⁸ The leading-edge gene pathways with adjusted p-values <0.05 were highlighted on the volcano plots.

Copy number variation (CNV) analysis in scRNA-seq

The CNV in each cell was estimated using the inferCNV R package (https://github.com/broadinstitute/inferCNV). Within this package, we adopted 10× default parameter values and estimated the CNV of epithelial cells using myeloid cells as a reference.

Processing of public ovarian cancer datasets

To compare sampling sites, we obtained scRNA-seq data (count matrix and metadata) from multiple tissues (ovary, omentum, and peritoneum) of a single patient with HGSOC, as provided by Olbrecht et al. (http://blueprint.lambrechtslab.org/).¹⁹ We classified T cells using pre-annotated metadata and then compared the scores of activated (CXCL13, CD200, ICOS), resting (CCR7, IL7R, GPR183, LMNA), and regulatory (IL2RA, TNFRSF9, FOXP3, CTLA4) T cells across different sites. We computed the scores using the AddModuleScore function of Seurat and compared them on a per-site basis.

To validate T-cell phenotype, we obtained ovarian cancer scRNA-seq data (collected from Gene Expression Omnibus [GEO], GSE165897) from 22 samples.²⁰ Cells pre-annotated as “T cells” were separated, and clustering was conducted using the same process as previously described in Section 2.3, with the pc/resolution option set at 21/1.6. Two patients were excluded from the analysis because either their before or after samples contained fewer than 100 cells each.

Trajectory analysis of T-cell clusters

We performed trajectory analysis using Monocle3 (version 1.3.1).²¹ To directly compare pseudo-time and trajectory with the pre-annotated clusters, we utilized Uniform Manifold Approximation and Projection (UMAP) constructed using Seurat as the reduction method and ordered the cells based on the pseudo-time. A naïve T-cell cluster was selected as the root cell.

T-cell receptor repertoire analysis

We harnessed the Cellranger V(D)J pipeline to generate TCR clonotype and complementarity-determining region 3 (CDR3) sequence information from TCR-seq data using the Cellranger V(D)J pipeline. Using the R package scRepertoire,²² we generated clonotype information (alpha and beta chain sequences) for each cell barcode based on the filtered contig outputs. We used the “gene + nt” option in the cloneCall function to assign T-cell clonotypes.

RESULTS Outline of single-cell RNA-seq data on patients with HGSOC

Figure 1A illustrates the streamlined research workflow. Patient specimens, obtained through surgical resection of ovaries, fallopian tubes, omentum, and para-aortic lymph nodes, all exhibited HGSC histologically (Table S1). Among the ten samples collected from nine patients, six originated from patients who had not undergone chemotherapy before surgery (T92, T66, T17, T18, T39 and T34). Out of the four cases obtained post chemotherapy, two underwent NACT (T107 and T110), with one retrieved from a recurring tumor 10 months after debulking operation and additional post-op chemotherapy (T107). The other two cases involved recurring tumors after post-op chemotherapy (T70 and T04). Notably, T92 and T110 were collected from the same patient as pre- and post-chemotherapy samples. The chemotherapy protocol included Paclitaxel and Carboplatin, with three cycles for NACT and six cycles for post-op chemotherapy.

FIGURE 1. Outline of single-cell RNA-seq data on patients with high-grade serous ovarian cancer (HGSOC). (A) Schematic illustration of the experimental workflow. Six pre-treatment and four post-treatment samples were analyzed from nine patients (created with biorender.com). (B) Uniform Manifold Approximation and Projection (UMAP) plots of 29,794 cells showing 6 cell type lineages, 24 clusters, and 10 samples. (C) Average expression levels and fraction of canonical marker gene expression in 24 clusters. (D) Proportion bar plot of 10 samples for 6 cell types, grouped by treatment state.

After enzymatic dissociation of tumor specimens, we applied 5′ scRNA-seq and paired V(D)J TCR-seq²³ to characterize both the cellular landscape and T-cell clonotypes (Figure 1A and Table S2). In the analysis, ambient RNAs were removed before integrating the multipatient data. Following quality control, we acquired 29,039 cells (Figure S1A) and performed clustering analysis with batch correction. The resulting 25 clusters were categorized into six cell types based on the canonical marker gene expression (Figure 1B,C). These cell types encompassed epithelial cells (EPCAM, CD24, KRT8, KRT18, KRT19), myeloid cells (LYZ, CD163, CD68, FCGR3A), fibroblasts (DCN, BGN, THY1, COL1A1, COL1A2), T/NK cells (CD3D, CD3E, TRAC, NCAM1, KLRD1), B/plasma cells (CD79A, BANK1, MZB1), and endothelial cells (PECAM1, CLDN5, FLT1, RAMP2). The UMAP space illustrated the distribution of these identified cell types (Figure S1B,C). Upon analyzing the cellular composition within each sample, the cellular landscape showed heterogeneity in epithelial cell (potential tumor cell population) content as well as in stromal and immune cell components (Figure 1D). Afterward, we compared pre- and post-chemotherapy samples within the segregated cell types.

Chemotherapy-induced alterations in the tumor signature

At first, we performed a sub-clustering analysis on the epithelial cell type and identified patient-specific clusters without batch correction (Figure 2A). These epithelial sub-clusters exhibited distinct chromosomal gene expression patterns (Figure S2A), indicative of copy number aberrations, marking them as tumor cell clusters. Notably, those epithelial tumor cells from the same patient (T92 and T110) occupied an identical cluster in the UMAP and shared similar chromosomal gene expression patterns, demonstrating minimal clonal evolution.

FIGURE 2. Chemotherapy-induced alterations in the tumor signature. (A) Uniform Manifold Approximation and Projection (UMAP) plots of 9891 epithelial cells colored by sample (left) and treatment (right). (B) Heatmap plot of the gene set variation analysis (GSEA) showing the average pathway score in each sample. (C) UMAP plots of Vision pathway cell-level scoring result. (D, E) Volcano plots (left) showing cell-level differentially expressed genes (DEGs) between post- and pre-treatment samples (D) and paired sample (T110 vs. T92, E). Red and blue dots indicate up- and downregulated genes, respectively. Bar plots (right) show normalized enrichment scores (NES) from GSEA for the DEGs. Bold and black colored letters of hallmark pathways indicate overlapped pathways in the pair and total both (NES > 1.5).

To elucidate the distinct tumor cell signatures after chemotherapy, we compared gene expression profiles in pre- and post-treatment samples at both pseudo-bulk and individual cell levels. This dual approach was tailored to account for patient-specific variations, allowing the identification of chemotherapy-induced alterations.

In the pseudo-bulk comparison, we conducted a GSVA¹⁷ on the hallmark gene sets curated from MSigDB using sample-level counts generated by the AverageExpression function in Seurat (Figure 2B). T107 was excluded from the analysis due to a low tumor cell count (45 cells compared with >500 cells in other samples). Group_1 gene set, enriched in pre-treatment samples, represented E2F_TARGETS and G2M_CHECKPOINT, suggestive of robust proliferation and cell cycle progression. Group_2 reflected diverse metabolic pathways without a clear inclination toward a particular sample group. Group_3 showed enrichment in post-treatment samples (ES > 2) and signified EPITHELIAL_MESENCHYMAL_TRANSITION and ANGIOGENESIS. Group_4 encompassed immune pathways such as INTERFERON_ALPHA_ or GAMMA_RESPONSE, along with IL2_STAT5_SIGNALING but demonstrated no sample group inclination. Representative pathway scoring within the feature plot projected the pseudo-bulk assessment into single-cell levels for the Group_1 and Group_3 gene sets (Figure 2C).

Next, we conducted a DEG analysis between tumor cells from pre- and post-treatment samples for cell-level comparisons (Figure 2D, left and Table S3). Subsequently, a gene set enrichment analysis was performed on the DEGs (Figure 2D, right). The same procedure was repeated for the longitudinal samples from the singular patient (Figure 2E). The normalized enrichment scores (NES) of the two comparisons showed a significant positive correlation (R = 0.71, p = 9.7e−09) (Figure S2B). In the post-treatment group, pathway enrichments in TNFA_SIGNALING_VIA_NFKB, COAGULATION, and EPITHELIAL_MESENCHYMAL_TRANSITION were observed. In the pre-treatment group, E2F_TARGETS, OXIDATIVE_PHOSPHORYLATION, and G2M_CHECKPOINT were among the most enriched pathways (Figure 2D,E, right). Notably, the CellCycleScoring function in Seurat indicated that the majority of post-chemotherapy samples remained in the G1 phase of the cell cycle (Figure S2C). Overall, both pseudo-bulk and cell-level analyses yielded consistent results, suggesting that the transcriptional signatures for TNFA_SIGNALING_VIA_NFKB, COAGULATION, and EPITHELIAL_MESENCHYMAL_TRANSITION are associated with chemotherapy.

Chemotherapy-induced remodeling of the tumor immune microenvironment

We proceeded to assess the impact of chemotherapy on the tumor microenvironment, by comparing the pre- and post-treatment groups within each non-epithelial cell type identified in Figure 1B. Cell-level comparisons were carried out through DEG analysis within fibroblasts, endothelial cells, myeloid cells, B/plasma cells, and T/NK cells (Figures S3–S6A, upper left and Table S4). Further investigation involved subtype clustering analysis within each cell type (Figures S3–S6A, upper right, B). In this phase, all subtypes were annotated, and their proportions were compared between the pre- and post-treatment groups (Figures S3–S6C,D).

Firstly, fibroblasts were categorized into cancer associated fibroblasts (CAFs, FB0_MMP11 and FB3_STAR), adipogenic fibroblasts (FB1_CFD), myofibroblasts (FB2_MYH11), and mesothelial cells (FB4_CALB2) (Figure S3B).¹⁹ In the post-chemotherapy group, fibroblasts demonstrated upregulation of adipogenic genes such as C7, CFD, MGP, and ADIRF (Figure S3A,D).^24,25 Secondly, endothelial cells were clustered into those of post capillary venules (EC0_ACKR1), tip cells (EC1_ANGPT2), capillary (EC2_FABP4), pericytes (EC3_ACTA2), and lymphatic vessels (EC4_PROX1) (Figure S4B).^19,26 Interestingly, post-chemotherapy endothelial cells exhibited an upregulation of genes related to fatty acid uptake, including FABP4 and CD36 (Figure S4A).²⁷ As half of post-chemotherapy samples were collected from the omentum, the adipogenic gene expression is likely influenced by omental adipocytes in addition to the chemotherapy.²⁸ Thirdly, myeloid cells were clustered into macrophages (MC0_APOE, MC1_STMN1), monocytes (MC2_FCN1), and conventional or plasmacytoid dendritic cells (MC3_LAMP3 or MC4_IRF4) (Figure S5B).^29,30 Among the myeloid immune cell types, macrophages were the predominant subpopulation, marked by an increase in immediate-early gene expression (JUN, FOS) after chemotherapy (Figure S5A).³¹ Fourthly, B/plasma cell types were clustered into B cells, plasma cells, and cycling plasma cells. Among them, plasma cell populations were slightly underrepresented in the post-treatment samples (Figure S6). Lastly, T/NK cells were clustered into 14 clusters (Figure S7A, upper right) representing 11 distinct subpopulations (Figure 3A). Genes associated with functional differentiation showed disparate expression, such as NKG7, GZMB, ISG15, and LAG3 being over-represented in the pre-treatment group and IL7R and CCR7 in the post-treatment group, respectively (Figure S7A, upper left and Table S5).

FIGURE 3. Alterations in T-cell phenotype after chemotherapy. (A) Uniform Manifold Approximation and Projection (UMAP) plot of 7102 cells showing the T- and NK-cell subsets. (B) Average expression levels and fraction of canonical marker gene expression in assigned 10 T- and NK-cell subsets. (C, D) Bar plots for the proportion of CD4 (C) and CD8 (D) T-cell subsets in each patient grouped by chemotherapy treatment. Box plots indicate the proportion of T-cell subsets in the treatment groups (Exclude T39 and T107). Differences between the pre- and post-chemotherapy samples were assessed by t-tests; *p [less than] 0.05, **p [less than] 0.01, ***p [less than] 0.001. (E) UMAP feature plots showing interferon-stimulated gene (ISG) score (IFIT3, ISG15 and MX1) and specific gene (GZMK, GZMB, LAG3, PRF1, and GNLY) expression. (F) Heatmap showing the gene expression changes in CD8 T-cell subsets along the pseudo-time inferred from Monocle3.

Collectively, these alterations in gene expression and subpopulation proportions across diverse cellular components underscore the systemic remodeling effect of chemotherapy on the tumor microenvironment.

Alterations in T-cell phenotype after chemotherapy

In the subsequent analysis, we focused on T/NK subpopulations (Figure 3A). The T/NK cell type was largely categorized into CD4+ T cells (CD4, CD40LG), CD8+ T cells (CD8A, CD8B), proliferating cells (MKI67, TOP2A), and natural killer (NK) cells (KLRD1, NCAM1, FGFBP2). CD4+ and CD8+ T cells were further subcategorized using canonical markers for naïve/memory (SELL, CCR7, IL7R, GPR183, LMNA), interferon-stimulated gene (ISG)-positive (IFIT3, ISG15), follicular helper (CXCL13, CD200), regulatory (FOXP3, IL2RA), effector memory (GZMK, KLRG1, GZMH), and exhausted T cells (PDCD1, HAVCR2, LAG3) (Figure 3B and Figure S7B).³² Naïve CD4+ T cells were highly enriched in the T107 sample from the metastatic lymph nodes, indicating a tissue-specific phenotype (Figure 3C).³³ Consequently, T107 was excluded from further statistical analysis, along with T39 sample, which harbored a limited number of cells (31 cells).

The thorough UMAP visualization highlighted a clear separation of CD4+ T, CD8+ T, and NK cell populations along the y-axis (Figure 3A and Figure S7). Along the x-axis, marker expression associated with functional diversity was segregated, which coincided with cluster enrichment for the pre- or post-treatment groups even after batch correction (Figure S7A, lower right). Upon conducting a comparative analysis between the two groups, we observed a significant enrichment of CD4 memory T cells (Tm, p = 0.0002), but a depletion of CD4 follicular helper (Tfh) and CD4 regulatory (Treg) T-cell subpopulations (p = 0.018 and 0.03, respectively) in the post-treatment group (Figure 3C). Among the CD8+ T cells, we observed an enrichment of CD8 effector memory T cells (Tem, p = 0.046) but a depletion of CD8+ T cells expressing ISG or exhaustion marker genes (Tex) in the post-treatment group (p = 0.0054 for combined clusters) (Figure 3D). This enrichment pattern suggests that chemotherapy might induce selective depletion of late-differentiation-stage T cells, particularly with a regulatory and/or suppressive phenotype. We ruled out tissue-specific effects in the T-cell enrichment analysis by demonstrating similar T-cell phenotype scores in the ovarian, omental, and peritoneal tissues using a public scRNA-seq dataset (Olbrecht et al., 2019)¹⁹ (Figure S7C).

To gain insights into the relationship between T-cell clusters, we performed a trajectory analysis of CD4+ and CD8+ T-cell lineages. First, CD4+ T-cell differentiation initiated from the naïve/naïve-like T-cell (Tn) subpopulation, branching to Tm and ISG-positive T cells, and bifurcated to Tregs or Tfh (Figure S7D, upper). Second, CD8+ naïve/memory (Tn/m) subsets transitioned through Tem, which is CD8 and ISG positive, and finally to Tex states (Figure S7D, lower). These differentiation patterns were consistent with those observed in previous studies.³² Furthermore, effector CD8 T cells in ovarian cancer tissues could be classified into two major types: effector memory cells expressing GZMK and pre-dysfunctional T cells expressing GZMB (Figure 3E).³⁴ Interestingly, we found a subpopulation expressing both GZMK and GZMB, which was also ISG positive (CD8 T ISG). This subpopulation also showed strong expression of the LAG3 gene with a simultaneous expression of the cytotoxic marker genes PRF1 and GNLY (Figure 3E). The pseudo-time heatmap revealed that ISG expression peaks at an intermediate point between the dominant expression states of GZMK and GZMB (Figure 3F). When the ISG-positive subpopulation was combined with CD8 exhausted T cells, we observed a significant decrease in the post-treatment group compared with the pre-treatment group (p = 0.0054) (Figure 3D, lower right).

The altered T-cell profiles following chemotherapy were further validated in publicly available scRNA-seq data from 18 paired samples (treatment-naïve and post-NACT) originating from nine patients.²⁰ We utilized consistent clustering criteria from our dataset for the analysis of T-cell subpopulations, resulting in the identification of 21 clusters and annotation of 13 distinct cell subtypes (Figure S8A,B and Table S6). When comparing the overall T-cell subpopulation between the treatment-naïve and post-NACT groups (Figure S8C), we observed significant enrichment of CD4 memory T cell (CD4 Tm, p = 0.048) and CD8 central memory T cell (CD8 Tcm, p = 0.036), depletion of CD4 follicular helper T cell (CD4 Tfh, p = 0.057), CD8 exhausted T cell (CD8 Tex, p = 0.0039), and proliferating T cell (Proliferating, p = 0.011). Taken together, these findings correspond to outcomes observed within our dataset and suggest that the chemotherapy influences T-cell profiles within the tumor microenvironment in a consistent direction.

Tracking T-cell phenotype with the TCR clonotype

Next, we performed a TCR repertoire analysis to integrate T-cell clonotypes with gene expression phenotypes. Within the resting or memory phenotype T-cell clusters (CD4 Tn, CD8 Tn/m, and CD4 Tm), a significant portion of cells exhibited either a unique clonotype or low levels of clonal expansion (Figure 4A,B). In contrast, the activated subsets (Proliferative, CD4 Tfh, CD4 Treg, CD8 Tem, CD8 ISG, and CD8 Tex) showed higher degrees of clonal expansion, with CD8 T cells exhibiting higher overall expansion levels than in CD4 T cells. Notably, CD8 Tex cells displayed the highest clonal expansion and expressed residential memory T-cell markers ZNF683 and ITGAE, along with proliferative cells (Figure 4C). These findings align with previous studies suggesting that resident memory T cells, recognizing tumor-specific antigens, undergo extensive expansion and exhaustion.³⁵

FIGURE 4. Tracking T-cell phenotype with the T-cell receptor (TCR) clonotype. (A) Uniform Manifold Approximation and Projection (UMAP) plot showing the T-cell clonotype classification split by treatment, colored by clonal TCR frequency ranges. (B) Proportion bar plot showing clonotype distributions for each T-cell subset. (C) UMAP plots and violin plots showing residential memory T-cell marker gene expression (ITGAE and ZNF683). (D) Box plots showing clonotype proportions of T cells, grouped by treatment states. (E) Proportion bar plot showing clonality distribution for each sample. Differences between the pre- and post-chemotherapy samples were assessed by t-tests; *p [less than] 0.05, **p [less than] 0.01, ***p [less than] 0.001. (F) Venn diagram (left) of unique and shared T-cell clones in the paired sample (T92 and T110) and dot plot of two shared clonotypes showing subtypes and frequencies (right).

When comparing pre- and post-treatment groups, we observed prominent clonal expansion in the pre-treatment samples (Figure 4D,E). In contrast, unique T-cell clones dominated the post-treatment tissues (p = 0.031). Interestingly, previous reports have indicated a decrease in highly expanded T-cell clones in PBMC samples from a patient with HGSOC after chemotherapy,³⁶ suggesting a systemic depletion of dividing T cells.

For a more refined analysis, we examined T92 and T110 samples collected longitudinally from the same patient. In this pair, the effects of chemotherapy were evident in the depletion of highly expanded clones (Figure 4E). The two samples shared 22 T-cell clones, many of which were CD8 Tem that managed to withstand the impacts of chemotherapy (Figure 4F). Conversely, CD8 Tex/ISG clones were largely depleted after chemotherapy.

In summary, these results suggest that chemotherapy resets the immune system toward less suppressive states, offering a promising opportunity to rebuild the immune microenvironment.

DISCUSSION

In this study, we compared treatment-naïve and post-chemotherapy ovarian cancers using scRNA-seq to better understand the effects of such treatments on HGSOC. The identified cellular and molecular alterations offer insights that could guide the development of treatment strategies for patients resistant to chemotherapy.

In the previous studies, the application of scRNA-seq in ovarian cancer has revealed tumor cell heterogeneity as well as the immune and stromal landscapes of primary and metastatic tumor sites.^34,37–39 In these high-resolution studies, HGSOC molecular subtypes from previous tissue-level findings were refined and attributed to both tumor cell characteristics and the surrounding microenvironment.^37,39 In these studies, cancer-associated fibroblasts were the major population conferring high EMT molecular characteristics. Furthermore, immune cell profiling unveiled the quantitative dynamics and activation states of tumor-infiltrating T cells and macrophages.^32,34 Similarly, macrophages with anti-inflammatory M2 features and T cells exhibiting regulatory or exhausted phenotypes were abundant in HGSOC tumor tissues. The immunosuppressive microenvironment suggested that patients with HGSOC may benefit from immunotherapies designed to reprogram the suppressive immune milieu.

Our study recapitulated the HGSOC landscape, that is, heterogeneity in tumor signatures as well as suppressive/exhausted immune phenotype, and further demonstrated chemotherapy-induced alterations. First, many EMT signature genes (CLU, CRYAB, S100A1, CCN2, IGFBP7, and GADD45B) were upregulated in post-treatment tumor cells, which is often associated with resistance to chemotherapy (CLU, CRYAB, CCN2, IGFBP, and GADD45B).^40–45 In a study elucidating the role of IGFBP7 in conferring chemotherapy resistance in T-cell acute lymphoblastic leukemia, an IGF1-R inhibitor restored chemotherapy sensitivity.⁴³ Thus, expression or functional modulation of the identified EMT genes could serve as a promising avenue for intervention to overcome chemotherapy resistance. Second, genes and subpopulations associated with lipid uptake and metabolism were upregulated in the post-treatment fibroblasts and endothelial cells. Lipid metabolism in tumor cells and stroma has been linked to tumor progression and drug resistance.^46,47 It is also noteworthy that endothelial expression of FABP4, PRCP, and KLF4, which are upregulated after chemotherapy, has been linked to vascular injury repair.^48,49 Thus these genes may represent potential targets for anti-angiogenic therapies. Third, T cells at the late differentiation stage were reduced following chemotherapy, and naive-like or early-activated T-cell-populated post-chemotherapy tumors. These findings are consistent with those of previous studies demonstrating chemotherapy-induced immune activation.^7,34,50–52 The observed phenotypic change of T cells may indicate a reset within the exhausted or suppressive immune microenvironment,⁵² opening up an opportunity for therapeutic intervention to activate tumor-specific immune responses.

The limited sample size posed challenges to our comparative analysis. To address this limitation, we utilized publicly available scRNA-seq data,²⁰ which focused on the stress-associated gene expression program in tumor after chemotherapy. From this data, we extracted T cells and validated the reduction of suppressive T-cell populations following chemotherapy. Nevertheless, we emphasize the need for a more extensive sample set, particularly including longitudinal specimens from chemotherapy treatments, to strengthen our conclusions and identify more robust therapeutic targets.

Another key limitation in our study is the absence of clear mechanistic insights into the chemotherapy-induced changes observed. Specifically, it remains uncertain whether chemotherapy drugs directly trigger the upregulation of EMT-associated genes in tumors, or if this is a secondary effect mediated by cellular damage and repair processes. Although we hypothesize that the observed reduction in suppressive T-cell populations results from their preferential elimination post chemotherapy, this theory requires empirical validation. Given the constraints in obtaining human patient samples, these hypotheses could be rigorously tested in experimental mouse models over an extensive timeframe.

Collectively, we dissected chemotherapy-induced alterations in cellular resolution and identified candidate genes, pathways, and subpopulations that potentially confer chemotherapy resistance. The targets were derived from diverse cellular components, including the tumor, stroma, vasculature, and T cells, suggesting different therapeutic combinations. However, further research is needed to fully understand the effects and mechanisms of action of each target and how they can be leveraged to improve cancer treatment outcomes.

AUTHOR CONTRIBUTIONS

Huiram Kang: Formal analysis; visualization; writing – original draft. Sohyun Hwang: Writing – review and editing. Haeyoun Kang: Writing – review and editing. Areum Jo: Resources. Ji Min Lee: Resources. Jung Kyoon Choi: Conceptualization. Hee Jung An: Conceptualization; funding acquisition. Hae-Ock Lee: Conceptualization; funding acquisition; project administration; writing – original draft; writing – review and editing.

ACKNOWLEDGMENTS

We acknowledge the Basic Medical Science Facilitation Program through the Catholic Medical Center of the Catholic University of Korea funded by the Catholic Education Foundation and the KREONET/GLORIAD service provided by the Korea Institute of Science and Technology Information (KISTI).

FUNDING INFORMATION

This work was supported by National Research Foundation of Korea (NRF) grants funded by the Ministry of Science and ICT (RS-2023-00220840 and NRF-2019M3A9B6064691) and the Korean Health Technology R&D Project funded by the Ministry of Health & Welfare (HI16C1559).

CONFLICT OF INTEREST STATEMENT

The authors declare no conflict of interest.

DATA AVAILABILITY STATEMENT

Raw and processed single-cell RNA sequencing data were deposited in the Gene Expression Omnibus (GEO) under the accession numbers GSE192898¹¹ and GSE241221. The data sources and handling of the publicly available ovarian cancer dataset used in this study are described in Section 2.

ETHICS STATEMENT

Approval of the research protocol by an Institutional Reviewer Board: N/A.

Informed Consent: N/A.

Registry and the Registration No. of the study/trial: N/A.

Animal Studies: N/A.