Single-cell analysis revealed that IL4I1 promoted ovarian cancer progression
Journal of Translational Medicine volume 19, Article number: 454 (2021)
Ovarian cancer was one of the leading causes of female deaths. Patients with OC were essentially incurable and portends a poor prognosis, presumably because of profound genetic heterogeneity limiting reproducible prognostic classifications.
We comprehensively analyzed an ovarian cancer single-cell RNA sequencing dataset, GSE118828, and identified nine major cell types. Relationship between the clusters was explored with CellPhoneDB. A malignant epithelial cluster was confirmed using pseudotime analysis, CNV and GSVA. Furthermore, we constructed the prediction model (i.e., RiskScore) consisted of 10 prognosis-specific genes from 2397 malignant epithelial genes using the LASSO Cox regression algorithm based on public datasets. Then, the prognostic value of Riskscore was assessed with Kaplan–Meier survival analysis and time-dependent ROC curves. At last, a series of in-vitro assays were conducted to explore the roles of IL4I1, an important gene in Riskscore, in OC progression.
We found that macrophages possessed the most interaction pairs with other clusters, and M2-like TAMs were the dominant type of macrophages. C0 was identified as the malignant epithelial cluster. Patients with a lower RiskScore had a greater OS (log-rank P < 0.01). In training set, the AUC of RiskScore was 0.666, 0.743 and 0.809 in 1-year, 3-year and 5-year survival, respectively. This was also validated in another two cohorts. Moreover, downregulation of IL4I1 inhibited OC cells proliferation, migration and invasion.
Our work provide novel insights into our understanding of the heterogeneity among OCs, and would help elucidate the biology of OC and provide clinical guidance in prognosis for OC patients.
Ovarian cancer (OC) was one of the most fatal and aggressive tumors of the female reproductive system and had emerged at an increased incidence in recent years. Generally acknowledged treatment for OC is surgery followed by platinum and taxane-based combination chemotherapy. However, nearly 25 percent of OC patients were found to relapse within six months after combination therapy . Most patients finally died from metastasis, due to a lack of other treatments aimed at improving the prognosis of OC patients. Thus, it is necessary to recognize OC-associated risks and make accurate prognostic of OC.
Effort made in understanding prognosis [2, 3] and response to platinum-based chemotherapy [4,5,6] in OC was highly focused on profiling gene expression and genetic aberrations. A prior report pointed out that OC patients studied by TCGA, revealed that critical gene mutations drive the pathogenesis of OC, which include the TP53 driver mutation (95%), and other major target genes including CCNE1, MYC, TERT, and NF1. Known risk factors that can accelerate ovarian carcinoma progression include BRCA1/BRCA2 mutations, family history, pregnancy, and other factors [8, 9]. Although carcinogenic- and metastatic-specific mutations were confirmed to accelerate carcinogenesis, dysregulated signal transduction or genetic variation of tumor cells were also critical for cancer progression.
Exploring the mechanisms of OC progression using bulk transcriptomics are confounded by a variety of factors, including copy number variation and infiltration by non-cancerous cells. Understanding the relationship of cancer cells and the tumor microenvironment (TME) is dependent on the accurate identification and characterization of individual cell states. In addition, intra-tumorigenic heterogeneity represents a key mechanism for both OS and progression of cancer [10, 11]. Moreover, the extensive intra-tumor heterogeneity that exists between OC cells, makes accurate identification of genetic diversity based on bulk mRNA sequencing highly controversial. Actionable diagnostic markers and identified therapeutic targets were based on bulk profiling technologies, with a complete disregard for intra-tumoral heterogeneity, which was not suitable for all patients. Recent advances in single-cell sequencing provide powerful tools to explore genetic and functional heterogeneity, which should assist in resolving the problem. In addition, scRNA-seq studies have revealed new insights into intra-tumor heterogeneity and have identified distinct sub-populations, which have provided novel mechanisms in our understanding of both carcinogenesis and in revealing strategies for treatment [12,13,14,15]. However, few studies have explored OC at the single cell level.
One recent scRNA-seq study of OC performed by Shih et al.,  investigated intra- and inter-tumorigenic heterogeneity at the cellular resolution with a large number of samples. That study not only identified several cell clusters, but also found that specific cell types were correlated with a well-known cancer subtype. Their single-cell assessment of patient samples enhanced our cognition of OC, and provided critical information that was helpful in advancing our understanding of OC progression. Regretfully, the clinical application of markers across OC progression at the single cell level was not explored, and these markers might display greater precision in the context of personalized anti-cancer therapy with consideration for intra-tumoral heterogeneity.
To overcome this limitation, in present study we performed further bioinformatics analysis using the data from the Shih et al.,  study, with the aim of identifying several diverse clusters. The design of this study and workflow was summarized in Fig. 1. Significantly, we identified the dominant M2-like TAMs in OC and explored the critical pathways across tumor progression in OC. In addition, we developed a RiskScore that was associated with robust prognosis in OCs based on markers of OC progression. Interleukin 4 Induced 1 (IL4I1) was an important gene in Riskscore. Notably, a variety of bioinformatic methods and experimental assays were conducted, revealing that IL4I1 accelerated cell proliferation, migration and invasion. Our work would help elucidate the biology of OC based on single-cell RNA-sequencing, which might provide clinical guidance in prognosis for OC patients.
Ovarian cancer and other cancer datasets
Single-cell RNA-seq for ovarian cancers (GSE118828) was downloaded from GEO database (https://www.ncbi.nlm.nih.gov/); bulk RNA-seq data and corresponding clinic-pathological data of multiple cancer patients in TCGA were obtained from UCSC Xena (https://xenabrowser.net/datapages/). Data retrieved from multiple GEO databases was used for integrated analysis using the Combat with the sva package . All public data used in this study is described in Additional file 1: Table S1.
Single-cell RNA-seq data preprocessing
The barcode matrix was processed with Seurat v3 , a toolkit for single-cell RNA-seq data analysis. All functions were run with default parameters, unless otherwise specified. Low quality cells (< 300 genes/cell, and < 3 cells/gene) were excluded. The UMI count data was normalized by log-transformation. The top 2000 highly variable genes (HGVs) were selected to aggregate samples into a merged dataset. Next, the merged cells-by-genes matrix was scaled by dividing the centered expression by the standard deviation. Batch effects among patients were eliminated using the “RunHarmony” function with the harmony package . The top 20 principal components, along with HGVs were used in this process.
Subsequently, the main cell clusters were identified using the “FindClusters” function of Seurat and visualized using the t-distributed stochastic neighbor embedding (tSNE) function. DEGs were appraised using the “FindMarkers” or “FindAllMarkers” function with the default parameter. For sub-clustering analysis, we applied the same procedure of finding variable genes, dimensionality reduction, and clustering.
Cluster annotation was based on classical markers. We characterized the identities of cell types based on known markers: EPCAM and KRT18  (epithelium), ACTA2  (mesenchyme), CD163 and CD68  (macrophage), CCL5  and IL7R  (T cell), VWF and CDH5  (endothelium), CD79A and MS4A1  (B cell).
The chromosomal copy number variation (CNV) estimation
Initial CNVs for each region were estimated by the infer-CNV package . The CNV of total cell types were calculated by expression levels from single-cell sequencing data for each cell with a cut-off 0.1. The CNV score of each cell was calculated as the mean of the CNV region.
Single cell trajectories were performed to explore cell-state transitions using the Monocle2 package . Differentially expressed genes over the Pseudo-time among cluster cell transitions were calculated by the “differentialGeneTest” function (q value < 0.1). The “DDRTree” was applied to reduce dimensions and visualization functions, and the “plot_cell_trajectory” was used to plot the minimum spanning tree on cells.
Recognition of malignant transcription factors (TFs)
In order to identify malignant TFs, we extracted a list of all identified TFs from Animal TFDB 2.0 (http://bioinfo.life.hust.edu.cn/AnimalTFDB2/about.shtml). We compared the TFs list with 2397 DEGs from 2 different epithelial cells (C0 and C6), and identified the malignant TFs.
Gene set functional analysis
Gene set functional analysis was done by the clusterProfiler package , gsva package , and DAVID (https://david.ncifcrf.gov/home.jsp). The h.all.v7.0.symbols.gmt and c2.cp.reactome.v7.0.symbols.gmt were down-loaded from the Molecular Signatures Database (http://www.broad.mit.edu/gsea/msigdb/). To investigate IL4I1 mediated biological parameters in OC, gene set enrichment analysis (GSEA) was conducted by clusterProfiler package. False discovery rate (FDR) < 0.05 and P < 0.05 were utilized as the enriched terms. In addition, 308 OC patients in the TCGA dataset were divided into high-expression group and low-expression group according to the median value of IL4I1.
Construction a single-cell transcriptome network
To explore the relationship between the clusters, a Python-based computational analysis tool for single-cell RNA-seq data analysis—CellPhoneDB  was used to construct the single-cell transcriptome network. Namely, ligand-receptor pairs for each cluster with other clusters were generated.
Construction a RiskScore with malignant marker genes
To build a novel RiskScore associated with survival in OC, firstly, malignant marker genes associated with survival were selected at P < 0.05 with the univariate Cox proportional hazards regression model, 117 of 2397 genes were identified with statistical significance in the GEO OC meta-dataset1 (an integrated OC cohort: GSE14764, GSE23554, and GSE26712 with GPL96). Subsequently, the 117 genes were narrowed down using the penalized logistical least absolute shrinkage and selector operation (LASSO) algorithm. The GEO OC meta-dataset1 as the training cohort and TCGA OC dataset used as the testing cohort. Using OS as the predictor variable, this procedure was repeated 10,000 times to construct the RiskScore, which generated with a linear combination of expression values and LASSO coefficients of signature genes using the following formula:
Kaplan Meier survival analysis and time-dependent ROC curves were used to assess the performance of gene signatures. Patients were divided into a high- and low-RiskScore group based on the median value of RiskScore. The Kaplan–Meier survival curves of the RiskScore were generated using the log-rank test. The time-dependent receiver operating characteristic (ROC) curve was generated with the timeROC package. The prognostic value of the RiskScore was externally validated with GEO OC meta-dataset2 (an integrated OC cohort: GSE18520, GSE26193, GSE30161, GSE63885, GSE54388, and GSE9891 with GPL570), and other TCGA cancers.
Cell culture and siRNA transfection
All OC cell lines (SKOV3, A2780 and CAOV8) were obtained from ATCC. SKOV3 and A2780 were cultured in Roswell Park Memorial Institute (RPMI)-1640 medium supplemented with 10% fetal bovine serum (FBS) and 100 U/ml penicillin/streptomycin. CAOV8 was cultured in high-glucose Dulbecco’s modifed Eagle’s medium (DMEM) with 10% FBS and 100 U/ml penicillin/streptomycin. All cell lines were incubated at 37 ℃ with 5% CO2.
All cell lines are transfected with Lipofectamine™ RNAmax according to the manufacturers’ instructions. IL4I1-target specific small interfering RNA (siRNA) is synthesized by JTSBIO Co., Ltd. (Wuhan, China). The sequences of IFI6-target-specifc-siRNA (siIL4I1) were as follows: siRNA1, 5ʹ-GCAUGCAGGAUCCUGACUATTUAGUCAGGAUCCUGCAUGCTT-3ʹ; siRNA2, 5ʹ-GCGAUGAAGAAGUUUGAAATTUUUCAAACUUCUUCAUCGCTT-3ʹ, siRNA3, 5ʹ-GCUUCUUCUAUCUCAGCUUTTAAGCUGAGAUAGAAGAAGCTT-3ʹ; and the sequence of control is 5ʹ-UUCUCCGAACGUGUCACGUTTACGUGACACGUUCGGAGAATT-3ʹ. The transfected cells were cultured in a fresh medium for 24 h for the subsequent assays.
Western blot analysis
Total protein was extracted using RIPA buffer (Thermo Fisher Scientific, Waltham, MA, USA) and detected utilizing BCA assay (Thermo Fisher Scientific, Waltham, MA, USA). 30 μg of protein per sample was separated by SDS-PAGE, then transferred onto PVDF membrane (Gene Molecular Biotech, Inc., Shanghai, China). After obturated with 5% milk for 2 h at room temperature, the membrane was incubated overnight at 4 °C with primary antibodies as follow: GAPDH (1:1000, CST), MMP2 (1:1000, CST), MMP9 (1:1000, CST). Afterwards, the membrane was incubated with HRP-conjugated rabbit IgG secondary antibodies (1:7500, CST) for 1 h at room temperature, the expression level was measured with an ECL kit (Roche Diagnostics, Basel, Switzerland) by Western blot imaging system.
Cell proliferation assays
1 × 103 cells were incubated in 96-well plates. Then, 10 uL Cell Counting Kit-8 (CCK-8; Dojindo, Rockville, MD, USA) solution was added to each well and incubated for 2 h for evaluating cell proliferation. The absorbance of each well was measured at OD450 with a Tecan Infinite M1000 PRO (Tecan, Switzerland) from days 1 to 4.
Cell wound healing assay
The transfected SKOV3 were cultured in 6-well plates with 1 × 106 cells/well for 24 h. Then, the wounds were produced using a 100 μL pipette tip. The images were photographed with a microscope after 0 and 96 h. The scratch area was measured to evaluate cell migration ability by image J software.
The transwell migration and invasion assays were conducted with Corning Transwell Inserts (8.0 μm). For the transwell migration assay, 1.5 × 104 transfected cells suspended in 50 μL serum-free medium were placed in the upper chamber and 600 μL medium (10% FBS) was filled in the lower compartment. The cells were incubated at 37 ℃ for 24 h. The successfully translocated cells were fixed with 4% paraformaldehyde (PFA) and stained with 0.1% crystal violet, and counted in four randomly chosen fields (200×) under a microscope.
For the transwell invasion assay, 1.2 × 105 cells were seeded on transwell coated with 50 μL Matrigel (dilution of 1:4 with 0.2% BSA). It is worth noting that Matrigel was used to coat membranes for 12 h at 37 °C prior to invasion assays. The culture condition was the same as the transwell migration assay. The cells on the lower surface were fixed, stained, and photographed microscopically after 48 h.
All the statistical analyses were performed using R 3.6.1 and Graphpad prism 8.0 software. The assay was repeated at least three times and the data were presented as mean ± standard deviation (SD). Two-tailed Student’s t-test was used to assess the differences between two groups. P < 0.05 was considered statistically significant.
A single-cell atlas of ovarian cancers and control samples
Eighteen samples obtained from GSE118828 were studied, namely a total of 17 neoplasms including 1 benign cancer, 1 peritoneal cancer, 3 LGSOC, 12 HGSOC, and 1 control ovarian. After removing low quality cells, a total of 3066 cells were finally acquired. In Shih’s  paper, they only studied 14 samples and captured 2911 cells. A reason for disparate cell count between our study and Shih’s might be the inequable samples we analyzed. These cells were classified into 9 main cell lineages namely C0-C8 (Fig. 2A). The clinical information of each cell population was illustrated in Fig. 2B–F. The corresponding proportions for each cluster with different clinical characteristics were discrepant (Fig. 2G–K). Based on the expression of well-known markers (Fig. 2L), we found that the atlas mainly comprised two different types of epithelial cell (i.e., C0 and C6), two mesenchyme clusters (i.e., C1 and C4), a T cell cluster (C2), a macrophage cluster (C5), an endothelium cluster (C7), and a B cell cluster (C8). C3 was an XIST highly expressed cluster. The gene ontology (GO) analysis for each cluster also confirmed the cluster annotation (Fig. 2M).
To explore the relationship between these clusters in OC, we used CellPhoneDB to calculate potential ligand-receptor pairs of each cluster and construct the single-cell transcriptome network. Network was visible using Cytoscape (Fig. 2N). Notably, macrophages possessed the most interaction pairs with other clusters, revealing the dominant role of macrophages in OC (Fig. 2O). As shown in Fig. 2O, TNF, TNFRSF10B, LGALS9, CX3CR1, VEGFA, and LAMP1 secreted by C5 interact with receptors expressed on epithelial cells, mesenchymal cells, endothelial cells, and other immune cells. These ligand-receptor pairs may be related to immune, angiogenesis, and CAF proliferation.
Macrophages exhibit M2 polarization in OC
For discussing the heterogeneity among macrophages, 236 macrophages were reclustered into 2 subclusters basis on tSNE analysis, including 158 cells in S0 and 78 cells in S1 (Fig. 3A, B). Top 10 genes of each subcluster were shown Fig. 3C. Pseudotime analysis suggested that a differential process existed from S1 to S0 (Fig. 3D). In addition, we explored the stem score with GSVA based on multiple terms, it showed that higher stem score enriched in S1 (Fig. 3E). Thus, we inferred that S0 might originate from S1. Based on the trajectory graph, it seemed that S0 had different branches of differentiation, meaning that cells in S0 had different differential fate, and the heterogeneity existed in S0. Thus, we reclustered S0 cells to 4 cluster (Fig. 3F), mac_0 was consisted of most S0 cells (Fig. 3G). Subsequently, we explored their characteristics, top 5 markers of each cell types were shown in Fig. 3H. Obviously, the monolike marker, FCN1, was highly expressed in mac_1; the DC markers, CD1C/E, were highly expressed in mac_2; markers associated with proliferation related genes and functions were highly enriched in mac_3, such as TOP2A, G2M_CHECKPOINT, and etc. (Fig. 3H–J). To deeply expounded their clinical features in OC, we calculated the GSVA score using top 5 genes of each celltype. As a results, OC patients with higher mac_0 gsvascore demonstrated worse survival (Fig. 3K, L). According to violinplot, the M2-like TAMs marker, TGFBR2, was highly expressed in mac_0. Thus, we inferred that mac_0 represented an M2-like TAMs cluster, which was the dominant type of macrophages in OC. Furthermore, we explored the gsvascore of four celltypes in different clinical phenotypes (Fig. 3M). As a results, mac_3 gsvascore was highest in proliferation_subtype, which was accorded with our previous finds that mac_3 was associated with proliferation.
Malignant epithelial cells were distinguished in OCs
As we were aware, ovarian cancer mainly originated from the epithelium. In the present study, the presence of two different types of epithelial cells (C0 and C6) encouraged us to investigate their malignant status. The chromosomal copy number variation (CNV) score of each cell help to identify the malignant clusters. First, we calculated large-scale CNV in each cell type based on averaged expression patterns across intervals of the genome. For all cell groups and different subtypes, we found that C0 exhibited remarkably higher CNV levels than C6 and other types of cells (Fig. 4A). The CNV score of each subtype and each sample was shown in Additional file 1: Fig. S1.
Then, we analyzed the gene expression profiles along the trajectory of C0 and C6 with pseudotime analysis and demonstrated a differentiated process from C6 to C0 (Fig. 4B). Heatmap analysis showed dynamically changed genes from C6 to C0 and identified 3 clusters (Fig. 4C). Functional enrichment analysis illustrated that multiple oncogenesis pathways were activated during OC progression, such as the cAMP signaling pathway, the T cell receptor signaling pathway, pathways in cancer development and progression, and the HIF-1 signaling pathway (Fig. 4C). We also explored the expression of multiple TFs that associated with the tumorigenesis of OC, such as the oncogene MYC (Additional file 1: Fig. S2).
We compared discrepant genes between C0 and C6, identifying 1264 upregulated and 1133 down-regulated genes in C0 (Fig. 4D, Additional file 1: Table S2). The peculiar marker used in the clinical identification of the malignant OCs, PAX8, was highly expressed in C0 (Fig. 4D, E). In addition, GSVA demonstrated that C0 was significantly related to carcinogenic terms, such as MYC_TARGETS, G2M_CHECKPOINT, EPITHELIAL_MESENCHYMAL_TRANSITION, etc. (Fig. 4F). Thus, C0 was the malignant epithelial cluster. In addition, we explored the ligand-receptor pairs interactions between C0 and other cell types. As a results, C0 had the most relationship with C5 (Fig. 4G), certifying the close contact between tumor cells and macrophages in OC carcinogenesis.
Distinct subgroups in the epithelial cluster
To identify the heterogeneity of epithelial cells in OC, we reclustered the two former epithelial clusters. Four distinct subgroups in the malignant epithelial cluster were identified on basis of tSNE analysis (Fig. 5A). We noticed that specific markers of S0-S2 were related to immunity; CAV1 was related to carcinogenesis and was specifically shown in S3 (Fig. 5B). Pseudotime graphs demonstrated that S1 represented the original cells that could differentiate into S0 and S2 (Fig. 5C). Functional terms confirmed that S0-S2 were related to immunity, and S3 was related to carcinogenesis, including the Hippo signaling pathway and the PI3K-AKT signaling pathway (Fig. 5D, E).
Another epithelial cluster, C6, was reclustered to two subtypes (Fig. 5F). S0 was found to highly express ITGB4 and DDX3X, and was associated with the functions of RNA splicing and cell cycle checkpoint (Fig. 5G, H); and S1 specifically expressed BIRC3 and TM4SF1, which was related with protein targeting functions to the ER and neutrophil-mediated immunity (Fig. 5G, I). The above observations indicated that epithelial cells were closely associated with host immunity.
Construction a robust prognostic model based on malignant epithelial markers in OC
To explore the clinical application of gene expression patterns in malignant epithelial cells, we used univariate Cox proportional hazards regression and LASSO algorithms to narrow down 2397 malignant genes. The GEO OC meta-dataset1 (an integrated OC cohort: GSE14764, GSE23554, and GSE26712 with GPL96) was used as the training set and the TCGA OC was used as the testing set. Ultimately, we selected 10 prognosis-specific genes to construct the RiskScore. The formula of the RiskScore was Y = [ZNF440 × (− 0.237) + USP53 × (− 0.455) + TSPAN12 × 0.469 + RARS × (− 0.153) + NFX1 × (− 0.246) + LRRC6 × (− 0.369) + IL4I1 × (− 0.258) + GPC1 × (0.269) + CD59 × (− 0.386) + ARID5B × 0.334]. In addition, the prognostic value of the RiskScore was externally validated with GEO meta-dataset2 (an integrated OC cohort: GSE18520, GSE26193, GSE30161, GSE63885, GSE54388, and GSE9891 with GPL570).
RiskScore analysis for 10 specific biomarkers in OC patients was shown in Fig. 6A–C. Significantly, RiskScore was negative association with survival in OC patients (log-rank P < 0.01; Fig. 6D–F). RiskScore was reliable and robust to predict the survival of OCs based on time-dependent ROC curves (Fig. 6G–I). As a result, the area under curve (AUC) in GEO OC meta-dataset1 validation set was 0.666, 0.743 and 0.809 in 1-year, 3-year and 5-year survival, respectively. In TCGA testing set, the AUC was 0.629, 0.644 and 0.623 for 1-year, 3-year and 5-year survival, respectively. Meanwhile, we found that RiskScore also had a high predictive accuracy of survival in GEO OC meta-dataset2 validation set. Furthermore, we explored the prognostic capability of the RiskScore of multiple cancers in TCGA, and the forest plot demonstrated the diverse OS across multiple cancers (Additional file 1: Fig. S3). This indicated the good potential of RiskScore in survival monitoring.
IL4I1 accelerated OC cells proliferation, invasion and migration
IL4I1, an important gene in Riskscore model, was reported to play important roles in immunoregulation and tumor progression [29,30,31,32]. In present study, IL4I1 was higher in multiple cancers than paired normal samples in GTEx (Fig. 7A). Noteworthy was the observation that the representative protein expression level of IL4I1 was positive in OCs based on the Human Protein Atlas database (HPA) (Fig. 7B). In truth, IL4I1 upregulation was associated with poor OS in OC (Fig. 7C). GSEA shown that high-IL4I1 group were mainly associated with G2M Checkpoint and EMT (Fig. 7D). As validation, downregulation of IL4I1 notably inhibited the proliferation of SKOV3, A2780 and CAOV8 (Fig. 7E).
Given the effect of IL4I1 on OC cells metastasis, the expression of MMP2 and MMP9 were lower in siIL4I1 group (Fig. 7F). Furthermore, siIL4I1 significantly inhibited wound closure in SKOV3 compared to the control cells based on wound healing assay (Fig. 7G). In the transwell migration and invasion assay, migrated (Fig. 7H) and invasive (Fig. 7I) cells transfected with siIL4I1 were remarkably decreased, as compared to the control groups. Thus, IL4I1 was carcinogenesis in OC.
Ovarian cancer was characterized as having a poor rate of survival and being deficient in an effective treatment because of the inherent intra-tumoral heterogeneity. Studies have attempted to identify distinct sub-populations and to explore mechanisms in disease carcinogenesis and strategies for the treatment of multiple cancers with scRNA-seq [12,13,14,15]. However, the gene expression profiles of different gene clusters in OC are unclear. It is highly desirable to explore the heterogeneity of OC and to explore the underlying mechanisms that are crucial in improving OC prognosis.
According to Shih’s  study, they captured 2911 cells from 14 samples. However, we downloaded 18 samples from the GEO datasets and ultimately obtained 3066 cells which was different from them. In this study, we conducted a comprehensive single-cell expression atlas of the 3066 cells and identified nine diverse cell types of OC, including epithelial cluster, mesenchyme cluster, macrophage cluster, T cell cluster, endothelial cluster, and B cell cluster. Despite our cluster number was disparate with Shih’s 16 clusters, the cell types in our study were similar with that described by Shih et al., demonstrating the reliability of our methods. In our study, we undertook a more thorough analysis using their OC single-cell datasets. We identified M2-like TAMs were the dominant type of macrophages in OC; in addition, we explored the potential mechanism of tumorigenesis and clinical application of malignant epithelial markers in OC.
In this study, macrophages exhibited the most ligand-receptor pairs with other clusters, revealing its’ important role in OC. Macrophages were reclustered to two subtypes that revealed the heterogeneity of macrophages. Studies revealed that TAMs could promote tumorigenesis through TGFβ signaling in tumor cells [33,34,35]. In present study, we observed some carcinogenesis genes such as TGFBR2 were highly expressed in S0. GSVA demonstrated that more carcinogenesis hallmark terms existed in S0, such as EMT, angiogenesis, and PI3K_AKT_MTOR_SIGNALING, etc. In addition, pseudotime analysis illustrated a differentiated trajectory from S1 to S0. Thus, we inferred that S0 was M2-like TAMs, which accounted for a greater proportion of macrophages in OC. Crosstalk between high EMT, angiogenesis and TAMs might exist in OC, influencing OC progression. TAM-derived exosomes enriched miRNA, lncRNA, and specific proteins that were contribute to tumor cell dissemination in gastric cancer . A high density of CD206 + TAM was significantly associated with worse survival in colon cancer . CAF might have crosstalk with TME and contribute to cancer biology by inducing the EMT process . Tumour angiogenesis could be indirectly regulated by promoting M2 polarization of macrophages in liver cancer .
Epithelial ovarian cancer (EOC) accounts for the majority OC cases and advanced EOC eventually develops into a recurrent platinum-resistant disease. In the present study, we identified two epithelial cell types with different gene expression patterns, namely C0 and C6, which demonstrated heterogeneity in epithelial cells. C0 predominantly comprised HGSOC cells and C6 mainly comprised fallopian cells. Pseudotime analysis demonstrated a differentiated trajectory from C6 to C0; furthermore, CNV-based analysis demonstrated that C0 was a relatively malignant epithelial cluster. GSVA identified some carcinogenic terms that were enriched in C0, and the classical clinical marker for identification of malignant OCs, PAX8, was only found in C0. The above observations collectively demonstrated that C0 belonged to malignant epithelial cells and displayed wide divergence in diverse epithelial cells because of their differential origins. A similar study exploring an OC differentiated trajectory with single-cell sequencing has not previously been reported. Furthermore, our study found that the malignant epithelial cells possessed the most interaction pairs with macrophages, and macrophages play important roles in tumorigenesis and had crosstalk with tumor cells. Macrophages M2 polarization was related to malignancy events in cancers [40, 41]. Meaningfully, some immune-associated subtypes were identified for two epithelial clusters, which considered the heterogeneity of epithelial cells in OC and supported the important roles of immunity in OC disease progression.
To explore the clinical application of markers in malignant epithelial cells, we construct a Riskscore consisted of ten genes, which presented encouraging prognostic value in predicting survival in OCs. During the past years, several signatures have been identified for prognostic prediction based on bulk mRNA transcription dataset [42, 43]. However, they ignored the heterogeneity of tumors. In this study, RiskScore generated with carcinogenic genes in single-cell levels which was more credible. Interleukin-4 Induced gene 1 (IL4I1), an important gene in Riskscore model, was reported to play important roles in immunoregulation and tumor progression [29,30,31,32].
IL4I1 was strongly detected in the tumor bed of most human tumor types  and was identified as a prognostic biomarker [45, 46]. Our results were similar with them. In the current research, prognostic IL4I1 was higher in ovarian tumors compared to the normal ovary. Dysregulated signal transduction of tumor cells was critical for cancer progression, and carcinogenic- and metastatic-specific genes that were confirmed to accelerate carcinogenesis. In our study, we linked IL4I1 to tumor intrinsic malignant properties. As a result, we disclosed that IL4I1 promotes OC cell proliferation, migration and invasion. IL4I1 was a secreted L-phenylalanine oxidase expressed by antigen-presenting cells. IL4I1 increased the threshold of T-cell activation, inhibiting T-cell proliferation or differentiation [29, 44, 47], thus influencing immune microenvironment. Imbalance of tumor cell and immune cell might promote cancer progression. By enhancing systemic Trp-catabolism, IL4I1 contributed to a systemic tumor-promoting environment that allowed tumor cells to migrate and protected them from immune destruction . Thus, IL4I1 might promote tumor progression through influencing tumor cell motility and adaptive immunity, and regulating the priming of tumor-specific immune cell, such as T cells [29, 31, 32, 47], B cell [48, 49] and macrophage .
Our findings provided a new perspective for understanding the progression of OC. We identified the dominant M2-like TAMs in OC and recognized some novel markers of tumorigenesis. In addition, we developed a RiskScore that was associated with robust prognosis in OCs based on markers of OC progression. Furthermore, we demonstrated that IL4I1 was an oncogene and promoted OC progression. This approach was potentially helpful for personalized anti-cancer strategies in the setting of OC.
Availability of data and materials
All the datasets used in the present research are summarized in Additional file 1: Table S1.
The area under curve
Chromosomal copy number variation
Dulbecco’s modifed Eagle’s medium
Epithelial ovarian cancer
False discovery rate
Gene set enrichment analysis
Gene set variation analysis
The Human Protein Atlas database
Interleukin 4 Induced 1
The penalized logistical least absolute shrinkage and selector operation
Receiver operating characteristic
Roswell Park Memorial Institute
Fetal bovine serum
The t-distributed stochastic neighbor embedding
Miller DS, Blessing JA, Krasner CN, Mannel RS, Hanjani P, Pearl ML, Waggoner SE, Boardman CH. Phase II evaluation of pemetrexed in the treatment of recurrent or persistent platinum-resistant ovarian or primary peritoneal carcinoma: a study of the Gynecologic Oncology Group. J Clin Oncol. 2009;27:2686–91.
Ramus SJ, Song H, Dicks E, Tyrer JP, Rosenthal AN, Intermaggio MP, Fraser L, Gentry-Maharaj A, Hayward J, Philpott S, et al. Germline Mutations in the BRIP1, BARD1, PALB2, and NBN Genes in Women With Ovarian Cancer. J Natl Cancer Inst. 2015;107:2.
Millstein J, Budden T, Goode EL, Anglesio MS, Talhouk A, Intermaggio MP, Leong HS, Chen S, Elatre W, Gilks B, et al. Prognostic gene expression signature for high-grade serous ovarian cancer. Ann Oncol. 2020;9:56.
Sonego M, Pellarin I, Costa A, Vinciguerra GLR, Coan M, Kraut A, D’Andrea S, Dall’Acqua A, Castillo-Tong DC, Califano D, et al. USP1 links platinum resistance to cancer cell dissemination by regulating Snail stability. Sci Adv. 2019;5:3235.
Gralewska P, Gajek A, Marczak A, Rogalska A. Participation of the ATR/CHK1 pathway in replicative stress targeted therapy of high-grade ovarian cancer. J Hematol Oncol. 2020;13:39.
Oza AM, Estevez-Diz MDP, Grischke EM, Hall M, Marme F, Provencher DM, Uyar DS, Weberpals JI, Wenham RM, Laing N, et al: A biomarker-enriched, randomized Phase II trial of adavosertib (AZD1775) plus paclitaxel and carboplatin for women with platinum-sensitive TP53-mutant ovarian cancer. Clin Cancer Res. 2020.
Ahmed AA, Etemadmoghadam D, Temple J, Lynch AG, Riad M, Sharma R, Stewart C, Fereday S, Caldas C, Defazio A, et al. Driver mutations in TP53 are ubiquitous in high grade serous carcinoma of the ovary. J Pathol. 2010;221:49–56.
Bast RC Jr, Hennessy B, Mills GB. The biology of ovarian cancer: new opportunities for translation. Nat Rev Cancer. 2009;9:415–28.
Karlberg S, Tiitinen A, Lipsanen-Nyman M. Failure of sexual maturation in Mulibrey nanism. N Engl J Med. 2004;351:2559–60.
Marusyk A, Almendro V, Polyak K. Intra-tumour heterogeneity: a looking glass for cancer? Nat Rev Cancer. 2012;12:323–34.
Brock A, Chang H, Huang S. Non-genetic heterogeneity–a mutation-independent driving force for the somatic evolution of tumours. Nat Rev Genet. 2009;10:336–42.
Puram SV, Tirosh I, Parikh AS, Patel AP, Yizhak K, Gillespie S, Rodman C, Luo CL, Mroz EA, Emerick KS, et al. Single-Cell Transcriptomic Analysis of Primary and Metastatic Tumor Ecosystems in Head and Neck Cancer. Cell. 2017;171:1611–24.
Lambrechts D, Wauters E, Boeckx B, Aibar S, Nittner D, Burton O, Bassez A, Decaluwe H, Pircher A, Van den Eynde K, et al. Phenotype molding of stromal cells in the lung tumor microenvironment. Nat Med. 2018;24:1277–89.
Zhang L, Li Z, Skrzypczynska KM, Fang Q, Zhang W, O’Brien SA, He Y, Wang L, Zhang Q, Kim A, et al. Single-Cell Analyses Inform Mechanisms of Myeloid-Targeted Therapies in Colon Cancer. Cell. 2020;181:442–59.
Yost KE, Satpathy AT, Wells DK, Qi Y, Wang C, Kageyama R, McNamara KL, Granja JM, Sarin KY, Brown RA, et al. Clonal replacement of tumor-specific T cells following PD-1 blockade. Nat Med. 2019;25:1251–9.
Shih AJ, Menzin A, Whyte J, Lovecchio J, Liew A, Khalili H, Bhuiya T, Gregersen PK, Lee AT. Identification of grade and origin specific cell populations in serous epithelial ovarian cancer by single cell RNA-seq. PLoS ONE. 2018;13:e0206785.
Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012;28:882–3.
Butler A, Hoffman P, Smibert P, Papalexi E, Satija R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nat Biotechnol. 2018;36:411–20.
Korsunsky I, Millard N, Fan J, Slowikowski K, Zhang F, Wei K, Baglaenko Y, Brenner M, Loh PR, Raychaudhuri S. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat Methods. 2019;16:1289–96.
Winterhoff BJ, Maile M, Mitra AK, Sebe A, Bazzaro M, Geller MA, Abrahante JE, Klein M, Hellweg R, Mullany SA, et al. Single cell sequencing reveals heterogeneity within ovarian cancer epithelium and cancer associated stromal cells. Gynecol Oncol. 2017;144:598–606.
Schelker M, Feau S, Du J, Ranu N, Klipp E, MacBeath G, Schoeberl B, Raue A. Estimation of immune cell content in tumour tissue using single-cell RNA-seq data. Nat Commun. 2017;8:2032.
Young MD, Mitchell TJ, Vieira Braga FA, Tran MGB, Stewart BJ, Ferdinand JR, Collord G, Botting RA, Popescu DM, Loudon KW, et al. Single-cell transcriptomes from human kidneys reveal the cellular identity of renal tumors. Science. 2018;361:594–9.
Chen YP, Yin JH, Li WF, Li HJ, Chen DP, Zhang CJ, Lv JW, Wang YQ, Li XM, Li JY, et al. Single-cell transcriptomics reveals regulators underlying immune cell diversity and immune subtypes associated with prognosis in nasopharyngeal carcinoma. Cell Res. 2020;30:1024–42.
Patel AP, Tirosh I, Trombetta JJ, Shalek AK, Gillespie SM, Wakimoto H, Cahill DP, Nahed BV, Curry WT, Martuza RL, et al. Single-cell RNA-seq highlights intratumoral heterogeneity in primary glioblastoma. Science. 2014;344:1396–401.
Trapnell C, Cacchiarelli D, Grimsby J, Pokharel P, Li S, Morse M, Lennon NJ, Livak KJ, Mikkelsen TS, Rinn JL. The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells. Nat Biotechnol. 2014;32:381–6.
Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16:284–7.
Hanzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics. 2013;14:7.
Vento-Tormo R, Efremova M, Botting RA, Turco MY, Vento-Tormo M, Meyer KB, Park JE, Stephenson E, Polanski K, Goncalves A, et al. Single-cell reconstruction of the early maternal-fetal interface in humans. Nature. 2018;563:347–53.
Sadik A, Somarribas Patterson LF, Ozturk S, Mohapatra SR, Panitz V, Secker PF, Pfander P, Loth S, Salem H, Prentzell MT, et al. IL4I1 Is a Metabolic Immune Checkpoint that Activates the AHR and Promotes Tumor Progression. Cell. 2020;182:1252–70.
Molinier-Frenkel V, Prevost-Blondel A, Castellano F. The IL4I1 Enzyme: A New Player in the Immunosuppressive Tumor Microenvironment. Cells. 2019;8:8.
Aubatin A, Sako N, Decrouy X, Donnadieu E, Molinier-Frenkel V, Castellano F. IL4-induced gene 1 is secreted at the immune synapse and modulates TCR activation independently of its enzymatic activity. Eur J Immunol. 2018;48:106–19.
Mazzoni A, Capone M, Ramazzotti M, Vanni A, Locatello LG, Gallo O, De Palma R, Cosmi L, Liotta F, Annunziato F, Maggi L. IL4I1 Is Expressed by Head-Neck Cancer-Derived Mesenchymal Stromal Cells and Contributes to Suppress T Cell Proliferation. J Clin Med. 2021;10:2356.
Zhang B, Miao T, Shen X, Bao L, Zhang C, Yan C, Wei W, Chen J, Xiao L, Sun C, et al. EB virus-induced ATR activation accelerates nasopharyngeal carcinoma growth via M2-type macrophages polarization. Cell Death Dis. 2020;11:742.
Gomez V, Eykyn TR, Mustapha R, Flores-Borja F, Male V, Barber PR, Patsialou A, Green R, Panagaki F, Li CW, et al. Breast cancer-associated macrophages promote tumorigenesis by suppressing succinate dehydrogenase in tumor cells. Sci Signal. 2020;13:9.
Liu M, Zhong YB, Shao J, Zhang C, Shi C. Tumor-associated macrophages promote human hepatoma Huh-7 cell migration and invasion through the Gli2/IGF-II/ERK1/2 axis by secreting TGF-beta1. Cancer Biol Ther. 2020;9:1–10.
Zheng P, Luo Q, Wang W, Li J, Wang T, Wang P, Chen L, Zhang P, Chen H, Liu Y, et al. Tumor-associated macrophages-derived exosomes promote the migration of gastric cancer cells by transfer of functional Apolipoprotein E. Cell Death Dis. 2018;9:434.
Feng Q, Chang W, Mao Y, He G, Zheng P, Tang W, Wei Y, Ren L, Zhu D, Ji M, et al. Tumor-associated Macrophages as Prognostic and Predictive Biomarkers for Postoperative Adjuvant Chemotherapy in Patients with Stage II Colon Cancer. Clin Cancer Res. 2019;25:3896–907.
Guan Y, Du Y, Wang G, Gou H, Xue Y, Xu J, Li E, Chan DW, Wu D, Xu P, et al. Overexpression of PLXDC2 in Stromal Cell-Associated M2 Macrophages Is Related to EMT and the Progression of Gastric Cancer. Front Cell Dev Biol. 2021;9:673295.
Han C, Yang Y, Sheng Y, Wang J, Li W, Zhou X, Guo L. The mechanism of lncRNA-CRNDE in regulating tumour-associated macrophage M2 polarization and promoting tumour angiogenesis. J Cell Mol Med. 2021;25:4235–47.
Baek SH, Lee HW, Gangadaran P, Oh JM, Zhu L, Rajendran RL, Lee J, Ahn BC. Role of M2-like macrophages in the progression of ovarian cancer. Exp Cell Res. 2020;395:112211.
Ono K, Sogawa C, Kawai H, Tran MT, Taha EA, Lu Y, Oo MW, Okusha Y, Okamura H, Ibaragi S, et al. Triple knockdown of CDC37, HSP90-alpha and HSP90-beta diminishes extracellular vesicles-driven malignancy events and macrophage M2 polarization in oral cancer. J Extracell Vesicles. 2020;9:1769373.
Wu A, Zhang S, Liu J, Huang Y, Deng W, Shu G, Yin G. Integrated Analysis of Prognostic and Immune Associated Integrin Family in Ovarian Cancer. Front Genet. 2020;11:705.
Fan L, Lin Y, Lei H, Shu G, He L, Yan Z, Rihan H, Yin G. A newly defined risk signature, consisting of three m(6)A RNA methylation regulators, predicts the prognosis of ovarian cancer. Aging (Albany NY). 2020;12:18453–75.
Presset M, Djordjevic D, Dupont A, Le Gall E, Molinier-Frenkel V, Castellano F. Identification of inhibitors of the immunosuppressive enzyme IL4I1. Bioorg Chem. 2020;94:103463.
Liu M, Pan Q, Xiao R, Yu Y, Lu W, Wang L. A cluster of metabolism-related genes predict prognosis and progression of clear cell renal cell carcinoma. Sci Rep. 2020;10:12949.
Ramspott JP, Bekkat F, Bod L, Favier M, Terris B, Salomon A, Djerroudi L, Zaenker KS, Richard Y, Molinier-Frenkel V, et al. Emerging Role of IL-4-Induced Gene 1 as a Prognostic Biomarker Affecting the Local T-Cell Response in Human Cutaneous Melanoma. J Invest Dermatol. 2018;138:2625–34.
Puiffe ML, Dupont A, Sako N, Gatineau J, Cohen JL, Mestivier D, Lebon A, Prevost-Blondel A, Castellano F, Molinier-Frenkel V. IL4I1 Accelerates the Expansion of Effector CD8(+) T Cells at the Expense of Memory Precursors by Increasing the Threshold of T-Cell Activation. Front Immunol. 2020;11:600012.
Carbonnelle-Puscian A, Copie-Bergman C, Baia M, Martin-Garcia N, Allory Y, Haioun C, Cremades A, Abd-Alsamad I, Farcet JP, Gaulard P, et al. The novel immunosuppressive enzyme IL4I1 is expressed by neoplastic cells of several B-cell lymphomas and by tumor-associated macrophages. Leukemia. 2009;23:952–60.
Copie-Bergman C, Boulland ML, Dehoulle C, Moller P, Farcet JP, Dyer MJ, Haioun C, Romeo PH, Gaulard P, Leroy K. Interleukin 4-induced gene 1 is activated in primary mediastinal large B-cell lymphoma. Blood. 2003;101:2756–61.
We thank Mengcheng Yao supported the bioinformatic analysis.
This work was supported by the National Natural Science Foundation of China (Grant/Award No. 81672838 China), the Beijing Municipal Administration of Hospitals Clinical Medicine Development of Special Funding (Grant No. XMLX201705) and the Beijing Municipal Science and Technology Commission (Grant No. Z181100001718193).
Ethics approval and consent to participate
Consent for publication
The authors declare that there are no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
About this article
Cite this article
Zhao, H., Teng, Y., Hao, W. et al. Single-cell analysis revealed that IL4I1 promoted ovarian cancer progression. J Transl Med 19, 454 (2021). https://doi.org/10.1186/s12967-021-03123-7