Single-Cell RNA-Sequencing Atlas Reveals the Tumor Microenvironment of Metastatic High-Grade Serous Ovarian Carcinoma

Abstract

Ovarian cancer is the most common and lethal gynecological tumor in women worldwide. High-grade serous ovarian carcinoma (HGSOC) is one of the histological subtypes of epithelial ovarian cancer, accounting for 70%. It often occurs at later stages associated with a more fatal prognosis than endometrioid carcinomas (EC), another subtype of epithelial ovarian cancer. However, the molecular mechanism and biology underlying the metastatic HGSOC (HG_M) immunophenotype remain poorly elusive. Here, we performed single-cell RNA sequencing analyses of primary HGSOC (HG_P) samples, metastatic HGSOC (HG_M) samples, and endometrioid carcinomas (EC) samples. We found that ERBB2 and HOXB-AS3 genes were more amplified in metastasis tumors than in primary tumors. Notably, high-grade serous ovarian cancer metastases are accompanied by dysregulation of multiple pathways. Malignant cells with features of epithelial-mesenchymal transition (EMT) affiliated with poor overall survival were identified. In addition, cancer-associated fibroblasts with EMT-program were enriched in HG_M, participating in angiogenesis and immune regulation, such as IL6/STAT3 pathway activity. Compared with ECs, HGSOCs exhibited higher T cell infiltration. PRDM1 regulators may be involved in T cell exhaustion in ovarian cancer. The CX3CR1_macro subpopulation may play a role in promoting tumor progression in ovarian cancer with high expression of BAG3, IL1B, and VEGFA. The new targets we discovered in this study will be useful in the future, providing guidance on the treatment of ovarian cancer.

Keywords: T cells; high-grade serous ovarian carcinoma; myeloid cells; scRNA-seq; tumor microenvironment.

Figure 1

Overview of TME in primary HGSOCs, metastatic HGSOCs, and ECs. (A) Workflow of the samples collected and the data analysis strategy. (B) Cell populations identified. The UMAP projection of 55802 single cells from HG_P (n=5), HG_M (N=3), EC_P (n=2), HG_Nor (n=1) samples shows the 10 main clusters with annotation. Each dot corresponds to a single cell, colored according to cell type. (C) Canonical cell markers are used to identify the clusters. (D) Barplots of the cell type and cancer subtypes for all 11 tumors. (E) The cancer subtypes proportion for each pathological group.

Figure 2

Copy number profiles, intertumoral heterogeneity, and EMT signature subpopulations are identified. (A) The chromosomal landscape of copy number for 13,634 epithelial cells of seven primary tumors; amplification (red) and deletions (blue). (B) The chromosomal landscape of copy number for 2849 epithelial cells of metastatic tumors and primary tumors of the HG3 patient (L_HG3_P means the primary tumor from the left ovary in the HG3 patient; R_HG3_P means the primary tumor from the right ovary in HG3 patient). (C) The UMAP projection of 17,551 epithelial cells from 10 tumors of six patients (indicated by labels and colors) reveals tumor-specific clusters. (D) Differentially expressed genes of the top 10 genes (rows) that are differentially expressed in each cluster (columns). (E) Differentially expressed genes between Scissor⁺ cells and all other cells in HGSOCs, each point represents a gene. Red: significant genes; Black: NS genes. avg_logFC: log 2 fold-change of the average expression between the two groups. ((log-FC > 0.25, FDR <0.05) (F) Enrichment of significant genes related to Reactome and Hallmark pathways. (G) Kaplan-Meier plot shows that high expression of EMT signature has shorter overall survival in ovarian cancer. The high and low patients are split by the mean expression of the EMT-related gene set.

Figure 3

Trajectory reconstruction during metastatic HGSOCs. (A) Monocle2 infers the development of epithelial cells along with pseudo-time (from patients HG3 and HG4 respectively, L_HG3_P means the primary tumor from the left ovary). Pseudo-time legend from dark to bright indicates cancer progression from the early to late stage. (B) Genswitches deduces the genes switch between cell states (left: L_HG3; right: HG4). (C, D) The heatmap displays the dynamic gene expression profiles during metastasis of ovarian cancer (from patients HG3 and HG4 respectively). The color key from blue to red indicates relative expression levels from low to light. The top annotated GO and KEGG terms in each cluster are shown. (E, F) Top 100 differentially expressed transcription factor genes (TFs; left) and the expression of specific TFs are on view along with the pseudo-time curve in (right). (G) Overexpression of proliferation and metastasis-related genes predicts poor prognosis in HGSOCs.

Figure 4

Diversity of fibroblasts in HG_M. (A) The UMAP projection of 12,236 fibroblast cells of 11 samples from six patients (indicated by labels and colors). (B) Proportion and cell number of each fibroblast subtype in 11 samples. (C) Heatmap of marker genes expression. (D) Heatmap of functional gene sets. (E) GSVA analysis of differential pathways is scored per cell among five fibroblast subsets. (F) Active regulons in each fibroblast subsets.

Figure 5

Subpopulations of tumor-infiltrating lymphocytes (TILs) in HG_M. (A) The UMAP projection of 7967 TILs of 11 samples from six patients (indicated by labels and colors). (B) Proportion and cell number of each subtype in 11 samples. (C) Dot plot (left) and UMAP-plot (right) display canonical cell markers. (D) Hierarchical clustering heatmap groups the tumors between HG_P, HG_M, and EC_P. (E) Reconstruction trajectory of CD8⁺ T cells inferred by Monocle2 (color by subtypes, expression of signature genes, and pseudotime). (F) Heatmap of the functional gene sets in TILs. (G) Active Regulons in each TILs. (H) Overexpression of the PRDM1 gene predicts a worse prognosis in ovarian cancer. (I) Cumulative distribution of cytotoxic CD8⁺ T cells between HG_P, HG_M, and EC_P. The cytotoxic score is calculated based on the average expression of cytotoxic markers. P-value was calculated by a two-sided unpaired Kruskal-Wallis rank-sum test.

Figure 6

Subpopulations of myeloid cells in HG_M. (A) The UMAP projection of 7265 myeloid cells of 11 samples from six patients (indicated by labels and colors). (B) Proportion and cell number of each myeloid subtype in 11 samples. (C) RNA velocity of each myeloid subtype. (D) The dot plot displays canonical cell markers. (E) Trajectory reconstruction of monocyte evolved into macrophages. (F) Dynamics gene expression profile during monocyte-to-macrophage terminal differentiation. (G, H) Biological processes enrichment analysis of module 2 and module 4. (I) Heatmap of significant genes in each subtype. (J) Activate regulons in each myeloid subtype.

Figure 7

The intricate intercellular interplay in HG_P, HG_M, and EC_P. (A) Circos plot shows the intercellular interactions in HG_P, HG_M, and EC_P. Each line represents an interaction where one end represents a ligand that is expressed in one cell type and the other end represents a receptor that is expressed in another cell type. The thickness of each line corresponds to the number of distinct interacting pairs. (B, C) Dot plot shows the means of the average expression levels and the possibility of occurrence in selective interaction pairs.