Molecular classification and immunologic characteristics of immunoreactive high‐grade serous ovarian cancer

Abstract High‐grade serous ovarian cancer (HGS‐OvCa) is one of the most lethal gynaecological malignancies. Molecular classification identified an immunoreactive subtype of HGS‐OvCa; however, the immunologic characteristics of immunoreactive HGS‐OvcA remain unclear. In this study, 121 immunoreactive HGS‐OvCa samples were identified from a meta‐analysis of 5 large transcriptome profiling data sets using a cross‐platform immunoreactive HGS‐OvCa subgroup‐specific classifier. By comparing the gene expression profiles of immunoreactive HGS‐OvCa samples and normal tissues, 653 differentially expressed genes (DEGs) were identified. KEGG pathway analysis revealed that the leukocyte transendothelial migration pathways were significantly enriched in the immunoreactive HGS‐OvCa. Protein‐protein interaction analysis identified a module that showed strong involvement of the immune‐related chemokine signalling pathway. Moreover, the GSEA enrichment analysis showed a T‐cell subgroup and M1 macrophages were significantly enriched in immunoreactive OvCa compared with normal samples. Macrophage infiltration levels were significantly elevated in immunoreactive HGS‐OvCa compared with other OvCa subtypes. In addition, expression of immune checkpoint molecules VTCN1 and IDO1 was significantly increased in immunoreactive HGS‐OvCa. In summary, our results suggest that the immunoreactive HGS‐OvCa has unique molecular characteristics and a tumour‐associated immune microenvironment featured by increased infiltration of macrophages, rather than lymphocytes. VTCN1 could be potential targets for the treatment of immunoreactive HGS‐OvCa.


| BACKG ROU N D
High-grade serous ovarian cancer (HGS-OvCa) is one of the most common malignancies of the female reproductive system, ranking third in morbidity among all gynaecological cancers. 1 Although surgery combined with chemotherapy is used as a standard of care, 75% of treated patients may experience drug resistance, relapse and short survival times. Only a fraction of the treatment-sensitive patients have a long disease-free survival. 2 The clinical application of immune checkpoint inhibitors is one of the great successes of anticancer treatment in recent years. 3 The immune checkpoints mediate the balance between immune surveillance and immune escape. 4 Immune checkpoint inhibitors such as PD-1/PD-L1 inhibitors have shown promising antitumour activity and limited adverse effects in treating several types of cancers. 5 Although HGS-OvCa had a strong immune recognition susceptibility, 6 the response of HGS-OvCa to anti-PD-1/PD-L1 monotherapy was minimal in some patients. The lack of therapeutic efficacy may be because of the insufficient and heterogeneous expression of PD-1 in the tumour-associated microenvironment of HGS-OvCa. 7 However, there still lacks an effective way to distinguish the immunotherapy-sensitive cohort in HGS-OvCa patients.
Cancer heterogeneity is directly related to disease progression and patient prognosis. Take breast cancer as an example-it was classified into five main subtypes: luminal A, luminal B, HER2overexpression, basal-like and normal-like. The prognosis of HER2positive breast cancer patients who use trastuzumab as a targeted treatment is significantly better than that of patients with basal-like breast cancers. 8 HGS-OvCa also exhibits high heterogeneity among its molecular characteristics. Similar to breast cancer, the classifications of different subtypes of HGS-OvCa might be useful to predict the drug sensitivity and the treatment outcome of patients. 9 The application of high-throughput gene expression profiling methods enables accurate identification of the molecular subtypes of HGS-OvCa. Tothill et al 10 identified four distinct molecular subtypes of HGS-OvCa. Among them, C1 (mesenchymal) subtype correlates with the high stromal response, the C2 (immunoreactive) subgroup exhibits high expression of immune cell-related genes, the C4 (differentiated) subtype shares some common features with serous borderline tumours, and the C5 (proliferative) subtype demonstrates low expression levels of differentiation markers. The Cancer Genome Atlas (TCGA) Research Network showed that immunoreactive HGS-OvCa exhibits high expression levels of the T-cell chemokine ligands CXCL11 and CXCL10, and its receptor CXCR3. They suggested that the C2 subset has a unique immune microenvironment. Thus, this subset may benefit from immune checkpoint-targeted treatments.
However, controversy exists on whether the prognosis of patients could be affected by the molecular subtypes of HGS-OvCa.
Zhang et al 11 classified the survival pattern of four subtypes of HGS-OvCa using the data from TCGA. They concluded that it was the tumour-associated stroma, not HGS-OvCa subtypes, that was associated with the patient's prognosis. 11 On the other hand, Shilpi et al 12 developed a novel classification system using the exon array and RNA sequencing data of HGS-OvCa from TCGA.
They showed that the molecular subtypes of HGS-OvCa could stratify patients into different survival patterns. 12 More importantly, they revealed that the prognosis of the immunoreactive subtype was not the best, which was unexpected. 12  The aim of the present study was to elucidate the unique molecular and immune characteristics of the immunoreactive HGS-OvCa and to identify the potential immune checkpoint inhibitors for its treatment. Here, we established a cross-platform classifier to distinguish the immunoreactive subtypes based on gene expression profiles. We identified the differentially expressed genes (DEGs) between immunoreactive HGS-OvCa tissues and normal tissues according to the classifier. A series of bioinformatic analyses were conducted to investigate the distinct molecular characteristics of immunoreactive HGS-OvCa. After that, we compared the enrichment status of immune cells between immunoreactive OvCa tissue and the normal tissue. Then, the comparison of the immune cell abundance and fractions between subtypes was also performed. The immune checkpoints expression patterns in subgroups were evaluated and validated by multiple methods using varies datasets. 14 Finally, we further explored the possible mechanism that maintains the immune-balanced microenvironment in immunoreactive HGS-OvCa.

| Establishment of a cross-platform classifier of immunoreactive HGS-OvCa
The mRNA profile of ovarian cancer and related data (ovarian cancer, RNA Seq V2, Illumina GA-DNASeq) in the TCGA database was acquired from the UCSC cancer genome database (http://xena.ucsc. edu/). The expression data of genes shared across multiple data sets including the TCGA dataset and 5 Gene Expression Omnibus (GEO) datasets (GSE06008, GSE18520, GSE26712, GSE27651, GSE9891) were extracted (Table 1). These gene expression values were normalized and scaled using the scikit-learn library in R. The 299 HGS-OvCa samples from TCGA were sorted into four subtypes, including differential, immunoreactive, proliferative and mesenchymal types, according to TCGA classification (Additional file S1). 15 The top 50 expressed genes were selected as feature genes based on the filter methods for feature selection 16 and an expression matrix was formed. Then, the 299 well-characterized samples from the TCGA mRNA expression data set were randomly divided into a 250-case training cohort and a 49-case validation cohort, respectively. Samples in the training cohort were classified as either the immunoreactive type or non-immunoreactive type by applying cluster analysis using the BP neural network model in Python. The same method was applied to the validation cohort to verify the accuracy of the model.

| A meta-analysis based on the crossplatform model
GEO2R (http://www.ncbi.nlm.nih.gov/geo/geo2r) were used to computed the P-value, adjust. P value, t value, LogFC value and B value of each sequencing site of the corresponding probe. For the same gene, the most significant sequencing site remains and others were excluded. Then, the commonly expressed genes among the four GEO data sets were obtained using Perl and the Merge package in R software. Two meta-analyses were performed on four GSE data sets, including cancer tissues, and the control group, using the MAMA and RankProd package. Z scores (which had |7| as a cut-off value) and the pval-test (which had |5| as a cut-off value) were used to filter the DEGs. Genes that met the above criteria were regarded as the final selected DEGs.

| GO annotations and KEGG pathway enrichment analysis
The Gene Ontology (GO) 17

| Identification of hub genes and significant modules
Hub genes were identified based on eigenvector centrality using CentiScaPe 2.1. 22 Eigenvector centrality is a measurement to evaluate the influence of a node in a certain network. The most significant cluster gene module was chosen with the condition that the degree cut-off = 2, node score cut-off = 0.2, k-core = 2, and max. depth = 100 using Molecular Complex Detection (MCODE) software. Moreover, the DEGs in each module with an FDR less than 0.1 were subjected to the WEB-based GEne SeT AnaLysis Toolkit to perform GO and KEGG analyses.

| Kaplan-Meier (KM) survival analysis
KM plotter (http://kmplot.com/analy sis/) is an effective online program that contains 1816 HGS-OvCa patients. 23 The gene expression file of the hub genes and relative survival information was extracted from the GEO, European Genome-phenome Archive (EGA) and TCGA databases. Survival analysis was performed to explore the relationships between the expression levels of selected genes and the prognosis of patients with HGS-OvCa.

| Analysis of enrichment and abundance of tumour-infiltrating immune cells in immunoreactive HGS-OvCa
Gene Set Enrichment Analysis (GSEA) was performed to evaluate the differential expression between immunoreactive HGS-OvCa and normal tissue based on a pre-defined leukocyte gene signature matrix using the WEB-based GEne SeT AnaLysis Toolkit. 24 To investigate the immune cell abundance in immunoreactive HGS-OvCa, we used the Tumor Immune Estimation Resource (TIMER; cistrome.shinyapps.io/timer) to estimate the abundance of immune cells (B cells, CD4 T cells, CD8 T cells, macrophages and dendritic cells) in TCGA samples. The abundance data were analysed and verified using pathological estimations. 25

| Characterization of the expression patterns of the selected genes in bulk expression and single-cell RNA-sequencing (scRNAseq) data sets
The expression levels of specific genes in four subtypes of HGS-OvCa were acquired from 5 GSE data sets. The comparison between subtypes was applied by using the beanplot package in R 3.5.1.
The raw data of scRNA-seq of HGS-OvCa were obtained from Shih et al 26 (GSE118828). Seurat, another R package, was used to analyse the scRNA-seq data. The cell population that expressed the selected genes were identified after data normalization, scaling, linear dimensional reduction and visualization using UMAP. 27

| Establishment of a cross-platform classifier of HGS-OvCa subtypes
To build a comprehensive cross-platform classifier, we first ex- The ROC curve suggested that the accuracy of the model applied to GSE98981 was 80.2%, with a 95% CI of 0.718-0.870 ( Figure 1B).
The expression profiles of marker genes selected by TCGA 15 further confirmed the accuracy of our classifier. As shown in Additional file S2, the expression levels of these genes among the different HGS-OvCa subtypes were consistent with the original conclusion.
For example, CXCL11 and CXCR3 were highly expressed in immunoreactive HGS-OvCa, while the expression of CXCL12 was low.

| Exploration of the molecular characteristics of immunoreactive HGS-OvCa
We compared the expression levels of 12,818 genes in the 5 GEO and CX3CL1 were associated with prognosis (Additional file S10).

| Investigation of the immune characteristics in immunoreactive HGS-OvCa
To understand the involvement of immune cell subsets in immunoreactive HGS-OvCa, GSEA was performed to investigate the

| The exploration of the relationships of macrophages and VTCN1 expression based on singlecell analysis
A previous study revealed that VTCN1 was mainly expressed in tumour-associated macrophages, which had negative effects on the T-cell response. 29 Our analysis above found that immunoreactive HGS-OvCa presents high abundant macrophage infiltration and up-regulated expression of VTCN1. We analysed a published single-cell RNA-seq data set to explore the relationship between VTCN1 and tumour-infiltrating macrophages. Macrophage signature CD68 was used to locate the macrophage cluster in the GESE118828 single-cell data set ( Figure 6A,B). 26 The results confirmed that the VTCN1 is indeed mainly expressed in the tumour-infiltrating macrophages ( Figure 6C), but not in tumour cells. Finally, we reported the overall survival based on VTCN1 expression in each subgroup of HGS-OvCa based on the TCGA cohort. It should be mentioned that although it did not reach the statistical significance, the high expression level of VTCN1 has the potential to be the risk factor ( Figure 6D). This result suggests

HGS-OvCa exhibits high levels of heterogeneity in its molecular
and histological characteristics. The immunoreactive subgroup has a distinct immune cell-infiltrated microenvironment and a generally favourable outcome. 30 The aim of the present study was to gain a deeper understanding of the underlying mechanism and to iden- HGS-OvCa patients. 41 Our results confirmed that VTCN1 was significantly up-regulated in the macrophages in immunoreactive HGS-OvCa ( Figure 6). found that although it did not reach the statistical significance, the high expression level of VTCN1 has the potential to be the risk factor.
In the present study, we built a cross-platform HGS-OvCa subgroup classifier and applied a set of integrated bioinformatic tools to systematically analyse the characteristics of immunoreactive HGS-OvCa. We clarified that TEM and the chemokine signalling pathway were specifically involved in the formation of the tumour-associated microenvironment. Moreover, macrophages showed significantly high expression levels in immunoreactive HGS-OvCa either compared with normal tissue or other subtypes. Two immune checkpoints were up-regulated in immunoreactive HGS-OvCa. Their differential expression contributes to the distinct immunoreactive HGS-OvCa tumour-associated microenvironment.

| CON CLUS IONS
In conclusion, the present study built an effective HGS-OvCa subtype classifier to sort immunoreactive HGS-OvCa samples in TCGA and GEO data sets. The unique molecular and immune characteristics of immunoreactive HGS-OvCa are revealed. Tang, #20161ACB20022 to to F. Fu).

CO N FLI C T O F I NTE R E S T
The authors declare that they have no competing interests.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are available from the corresponding author upon reasonable request.