A nine-gene signature related to tumor microenvironment predicts overall survival with ovarian cancer

Qi Ding; Shanshan Dong; Ranran Wang; Keqiang Zhang; Hui Wang; Xiao Zhou; Jing Wang; Kee Wong; Ying Long; Shuai Zhu; Weigang Wang; Huayi Ren; Yong Zeng

doi:10.18632/aging.102914

Research Paper Volume 12, Issue 6 pp 4879—4895

A nine-gene signature related to tumor microenvironment predicts overall survival with ovarian cancer

Qi Ding^{1,2,
,} , Shanshan Dong^{1,2,
,} , Ranran Wang^1,2, , Keqiang Zhang^4, , Hui Wang^3, , Xiao Zhou^1,2, , Jing Wang^4, , Kee Wong^2, , Ying Long^1, , Shuai Zhu^1, , Weigang Wang^4, , Huayi Ren^1, , Yong Zeng^1,2, ,

¹ Translational Medicine Center, The Affiliated Cancer Hospital of Xiangya School of Medicine, Central South University/Hunan Cancer Hospital, Changsha, China
² Engineering Technology Research Center for Diagnosis-Treatment and Application of Tumor Liquid Biopsy, Changsha, China
³ Key Laboratory of Radiation Oncology, Department of Radiation Oncology, Hunan Cancer Hospital and The Affiliated Cancer Hospital of Xiangya School of Medicine, Central South University, Changsha, China
⁴ The Fifth Department of Gynecological Oncology, The Affiliated Cancer Hospital of Xiangya School of Medicine, Central South University/Hunan Cancer Hospital, Changsha, China

* Equal contribution

Received: October 1, 2019 Accepted: March 2, 2020 Published: March 24, 2020

https://doi.org/10.18632/aging.102914
How to Cite

Copyright © 2020 Ding et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY 3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Abstract

Mounting evidence suggests that immune cell infiltration within the tumor microenvironment (TME) is a crucial regulator of carcinogenesis and therapeutic efficacy in ovarian cancer (OC). In this study, 593 OC patients from TCGA were divided into high and low score groups based on their immune/stromal scores resulting from analysis utilizing the ESTIMATE algorithm. Differential expression analysis revealed 294 intersecting genes that influencing both the immune and stromal scores. Further Cox regression analysis identified 34 differentially expressed genes (DEGs) as prognostic-related genes. Finally, the nine-gene signature was derived from the prognostic-related genes using a Least Absolute Shrinkage and Selection Operator (LASSO) and Cox regression. This nine-gene signature could effectively distinguish the high-risk patients in the training (TCGA database) and validation (GSE17260) cohorts (all p < 0.01). A time-dependent receiver operating characteristic (ROC) analysis showed that the nine-gene signature had a reasonable predictive accuracy (AUC = 0.707, AUC =0.696) in both cohorts. In addition, this nine-gene signature is associated with immune infiltration in TME by Gene Set Variation Analysis (GSVA), and can be used to predict the survival of patients with OC.

Introduction

Ovarian cancer (OC) is one of the highest mortality rate malignant tumors of the female reproductive system [1]. There are more than 239,000 new cases, and about 152,000 deaths worldwide from OC every year [2]. The standard treatment plan for this disease is tumor cytoreductive surgery combined with platinum-based chemotherapy [3]. In this treatment mode, more than two-thirds of patients have a total survival of less than 10 years, and the survival rate with advanced (III-IV) stage patients is less than 20% [1].

The occurrence, development, and therapeutic efficacy of OC are closely related to many factors such as disease pathological type, TNM stage, treatment timing, and endocrine level [2, 4–6]. Most recent studies revealed genetic changes are notably linked with the occurrence and the treatment efficacy of OC. For example, from the perspective of the disease occurrence, the low-grade serous, endometrioid, clear cell, and mucinous subtypes are characterized as genetically stable, showing local invasive growth; therefore the patient has a better prognosis. About 75% of OC are high-grade serous type and genetically unstable. Majority of the patients carry p53 mutation and possible BRCA1 and BRCA2 mutations. Clinical observations of this type of patients are usually accompanied by metastatic lesions and poor prognosis [7]. From the perspective of treatment efficacy, the global loss of 5-Hydroxymethylcytosine is associated with platinum drug resistance, shortened progression-free survival (PFS), and shortened overall survival (OS) in patients with high grade serous OC [8]. OC patients with high expression of Cyclin-dependent kinase 9 (CDK9) in relapsed and metastatic lesions have a worse prognosis than patients with low expression of CDK9 [9]. A number of studies in transcription and epigenetics have confirmed that the occurrence, development, and therapeutic efficacy of OC are influenced by the dynamic changes of multiple oncogenes and tumor suppressor genes [10–15].

Existing research shows that tumor cell and host cell interaction is an important factor in promoting tumor growth and disease progression [16]. Immune cells (T lymphocytes and tumor-associated macrophages), stromal cells (fibroblasts, etc.), and extracellular matrix together form a tumor microenvironment (TME) in cancer patients [17, 18]. This TME plays a role in disease progression and formation of metastatic lesions. For example, cancer-associated fibroblasts (CAFs) facilitate OC metastasis by promoting angiogenesis, lymphangiogenesis, and tumor cell invasion [19]. CAFs induce the upregulation of Lipoma-preferred partner (LPP) in microvascular endothelial cells that can lead to chemoresistance in OC [20]. Matrix Metallopeptidase 1 (MMP1) mRNA in extracellular vesicles (EVs) secreted by OC cells can induce apoptosis of peritoneal mesothelial cells, thereby destroying the peritoneal mesothelial barrier and promoting the transfer of tumor cells to the peritoneum in OC patients [21].

Large and complicated biological data has been generated with the advent of high-throughput detection technology and bioinformatics development. The Cancer Genome Atlas (TCGA) database is one of the largest cancer genome program that provides researchers with multi-omics and standardized clinical data that can be used to design basic bioinformatics research [22, 23]. The ESTIMATE algorithm can predict tumor purity by calculating the immune and stromal scores based on specific molecular biomarker expression in both immune and stromal cells [24]. Subsequently, ESTIMATE has been applied to many neoplasms, such as prostate cancer [25], glioblastoma [26], and clear cell renal cell carcinoma [27]. However, the immune/stromal scores of OC have not been investigated in detail.

In the present study, 593 OC patients were obtained from TCGA, and their immune/stromal scores were derived from ESTIMATE algorithm. The patients were divided into high and low immune/stromal score groups with the immune/stromal score median value as the cut-off value. Differential expression analysis revealed 294 intersecting genes that influencing both the immune and stromal scores. Further univariate Cox analysis narrowed the list down to 34 genes. A final prognostic nine-gene signature was derived with a Lasso-Cox regression analysis. The prognostic nine-gene signature was trained and validated on the TCGA and Gene Expression Omnibus (GEO) datasets respectively. Time-dependent receiver operating characteristic (ROC) analysis was used to evaluate the performance of the nine-gene signature. Functional enrichment analysis and GSVA as well as Tumor Immune Estimation Resource (TIMER) were used to elucidate the valuable gene-related functions in the TME. The findings indicate that the prognostic nine-gene signature could be used as a predictive tool to assess the survival rate of patients with OC and provide novel strategies for future immunotherapy.

Results

Clinical characteristics of the study patients

Figure 1 shows the workflow for the identification, validation, and functional analysis of the prognostic nine-gene signature. Four hundred sixty-five OC patients from the TCGA database were included as training cohort. One hundred and nine OC patients from the GEO dataset GSE17260 were used as a validation cohort. The detailed clinical characteristics of the training and validation cohort were summarized in Table 1.

The overall design of the study. TCGA-OC: TCGA-ovarian serous adenocarcinoma; ESTIMATE: Estimation of STromal and Immune cells in Malignant Tumor tissues using Expression data; LASSO: least absolute shrinkage and selection operator; GSVA: Gene Set Variation Analysis.

Figure 1. The overall design of the study. TCGA-OC: TCGA-ovarian serous adenocarcinoma; ESTIMATE: Estimation of STromal and Immune cells in Malignant Tumor tissues using Expression data; LASSO: least absolute shrinkage and selection operator; GSVA: Gene Set Variation Analysis.

Table 1. Baseline characteristics of study patients.

Variables	Training cohort No. (%)	Validation cohort No. (%)
No. of patients	465	109
Age (years)	59.68±11.49 (mean ± SD)
Vital status
Alive	207(44.5%)	63 (58.2%)
Dead	258(55.5%)	46 (41.8%)
FIGO stage
Stage II	24(5.2%)
Stage III	362(77.8%)	92(84.5%)
Stage IV	17(16.1%)	17(15.5%)
Unknown	4(0.9%)
Grade
GB	1(0.2%)
G1		26(23.6%)
G2	56(12.0%)	40(37.3%)
G3	397(85.4%)	43(39.1%)
G4	1(0.2%)
Unknown	10(2.2%)
Venous invasion
NO	52(11.2%)
YES	70(15.1%)
Unknown	343(73.8%)
Lymphatic invasion
NO	60(12.9%)
YES	112(24.1%)
Unknown	293(63.0%)
Tumor residual disease
No Macroscopic disease	88(18.9%)
1-10 mm	214(46.0%)
11-20 mm	28(6.0%)
>20 mm	85(18.3%)
Unknown	50(10.8%)

Analysis of differential gene expression profile with immune and stromal scores in OC

By comparing the gene expression profiles of patients with high immune scores against those with low immune scores, a total of 480 (438 upregulated, 42 downregulated) DEGs were identified (Figure 2A). Four hundred thirty-two (414 upregulated, 18 downregulated) DEGs were identified by comparing the high and low stromal score groups (Figure 2B). A fold-change > 1.5 and normalized p values < 0.05 were used as criterions for screening DEGs. A total of 281 DEGs were in common among the high immune/stromal score groups. A total of 13 DEGs were in common among the low immune/stromal score groups. (Figure 2E and 2F).

Differentially expressed genes based on immune scores and stromal scores. (A) The volcano plot showed that 438 genes were up-regulated and 42 genes down-regulated in the high immune scores group compared with the low scores group. (B) In a similar way, 414 upregulated genes and 18 downregulated genes were identified by comparing stromal scores. (C, D) Significantly enriched gene sets of the immune or stromal score group. (E, F) A total of 281 DEGs were in common among the high immune/stromal score groups and 13 DEGs in low immune/stromal score groups.

Figure 2. Differentially expressed genes based on immune scores and stromal scores. (A) The volcano plot showed that 438 genes were up-regulated and 42 genes down-regulated in the high immune scores group compared with the low scores group. (B) In a similar way, 414 upregulated genes and 18 downregulated genes were identified by comparing stromal scores. (C, D) Significantly enriched gene sets of the immune or stromal score group. (E, F) A total of 281 DEGs were in common among the high immune/stromal score groups and 13 DEGs in low immune/stromal score groups.

DEGs functional enrichment analysis

To dissect the underlying biological function of DEGs, we performed a functional enrichment analysis utilizing the R package clusterProfiler. Immune- related 480 DEGs were significantly enriched in Human T-cell leukemiavirus 1 infection, Human immunodeficiency virus 1, Human cytomegalovirus, Herpes simplex virus 1, and cytokine-cytokine receptor interaction (Figure 2C). The statistically significant pathways of 432 stromal-related DEGs are as follows: focal adhesion, human papillomavirus infection, PI3K-Akt signaling pathway, proteoglycans in cancer, and the cytokine-cytokine receptor interaction (Figure 2D).

Derivation of prognostic DEGs and construction of a gene risk score model

In the process of screening for prognostic-related biomarkers, the 294 DEGs in common among the high immune/stromal score and low immune/stromal score groups were subjected to a univariate Cox proportional hazard regression analysis. Out of which, 34 DEGs were found to be significantly (p < 0.05) correlated with the OS of the 465 OC patients from the TCGA database. Subsequently, these 34 candidate markers were used to construct a prognostic model with a Lasso-Cox proportional hazards regression. The resulting optimal prognostic signature for predicting the OS consists of nine genes: Ubiquitin D (UBD), V-Set And Immunoglobulin Domain Containing 4 (VSIG4), C-X-C Motif Chemokine Ligand 11 (CXCL11), Guanylate Binding Protein 2 (GBP2), C-X-C Motif Chemokine Ligand 13 (CXCL13), C-X3-C Motif Chemokine Receptor 1 (CX3CR1), Complement C5a Receptor 1 (C5AR1), Tissue Factor Pathway Inhibitor 2 (TFPI2), and DNA segment on chromosome 4 (unique) 234 expressed sequence (D4S234E). The Cox proportional hazard assumption was examined and validated through a Schoenfeld residuals test (P = 0.2259). The detailed information regarding the nine genes is provided in Table 2. The following is the formula for calculating the prognosis risk score:

Table 2. Nine prognostic genes significantly associated with OS in the training cohort.

Name	Coefficient	Type	Down/up-regulated	HR	95%CI	P value
UBD	-0.033	Protective	Up	0.90	0.85 - 0.96	<0.001
VSIG4	0.066	Risky	Up	1.14	1.05 - 1.25	0.002
CXCL11	-0.049	Protective	Up	0.86	0.80 - 0.93	<0.001
GBP2	-0.035	Protective	Up	0.90	0.82 – 1.00	0.048
CXCL13	-0.003	Protective	Up	0.81	0.73 - 0.91	<0.001
CX3CR1	0.009	Risky	Up	1.15	1.06 - 1.25	<0.001
C5AR1	0.009	Risky	Up	1.21	1.04 - 1.40	0.012
TFPI2	0.032	Risky	Up	1.08	1.01 - 1.15	0.028
D4S234E	-0.050	Protective	Down	0.90	0.83 - 0.98	0.011
Abbreviations: OS, overall survival; CI, confidence interval.

Risk score = (-0.033 × expressionUBD) + (0.066 × expressionVSIG4) + (-0.049 × expressionCXCL11) + (-0.035 × expressionGBP2) + (-0.003 × expressionCXCL13) + (0.009 × expressionCX3CR1) + (0.009 × expressionC5AR1) + (0.032 × expressionTFPI2) + (-0.050×expressionD4S234E). Each patient was assigned a risk score based on the formula and divided into either high-risk group or low-risk group according to the best cut-off in two cohorts. The distribution of the gene-based risk scores, OS, OS status, and the nine-gene expression profile of the patients in the training and validation cohorts are presented in Figure 3. The heat map showed that the five protective genes (UBD, CXCL11, GBP2, CXCL13, and D4S234E) exhibit low expression in the high-risk group. In contrast, the four risk genes (VSIG4, CX3CR1, CA5R1, and TFP12) have high expression in the high-risk group. Moreover, Kaplan-Meier curves were used to compare the OS of the two groups and the analysis showed that the OS of the high-risk group was substantially shorter than the low-risk group (p < 0.001; Figure 4).

Figure 3. The nine-gene signature predicts overall survival with ovarian cancer. (A, B) The distribution of risk score, overall survival, vital status, and the heat map of the nine gene expression profile in the training cohort and validation cohort.

Kaplan-Meier curves to compare overall survival of high-risk and low-risk groups based on the nine-gene signature in the training cohort (A) and validation cohort (B).

Figure 4. Kaplan-Meier curves to compare overall survival of high-risk and low-risk groups based on the nine-gene signature in the training cohort (A) and validation cohort (B).

The results of the univariate and multivariate Cox proportional hazard regression analyses identified that the nine-gene signature, age, and tumor residual disease as independent prognostic variables for the OS (Table 3).

Table 3. Univariate and multivariate Cox proportional hazards regression analyses in the training cohort.

Variables	Univariate analysis		Multivariate analysis
Variables	Hazard ratios (95%CI)	P-value	Hazard ratios (95%CI)	P-value
Age	1.022(1.012-1.033)	<0.001	1.017(1.006-1.028)	<0.001
Grade
G2	Referent
G3	1.165(0.820-1.654)	0.394
Unknown	1.194(0.531-2.684)	0.669
FIGO stage
II	Referent
III	2.355(1.109-5.001)	0.026
IV	2.961(1.350-6.495)	0.007
Unknown	3.923(0.814-18.91)	0.089
Venous invasion
No	Referent
Yes	0.967(0.560-1.671)	0.905
Unknown	1.249 (0.867-1.934)	0.318
Lymphatic invasion
No	Referent
Yes	1.264(0.798-2.001)	0.127
Unknown	1.094 (0.732-1.636)	0.374
Tumor residual disease
No Macroscopic disease	Referent		Referent
1-10 mm	1.899(1.324-2.722)	<0.001	1.469(1.006-2.144)	0.046
11-20 mm	2.191(1.259-3.814)	<0.001	2.034(1.140-3.629)	0.016
>20 mm	2.313(1.536-3.483)	<0.001	1.803(1.177-2.762)	0.007
Unknown	0.975(0.595-1.597)	0.919	1.034(0.626-1.709)	0.896
Nine-mRNA signature	21.48(10.15-45.42)	<0.001	15.60(6.963-34.96)	<0.001

Risk score model accuracy assessment

The time-dependent ROC curve analysis was conducted and the area under the curve (AUC) value was used to evaluate the predictive effect of the nine-gene signature. In the training cohort, the three-year AUC was 0.684 and the five-year AUC was 0.707. In the validation cohort, the three-year AUC was 0.606 and the five-year AUC was 0.696 (Figure 5A and 5B). By comparing the nine-gene signature against other prognostic factors and single gene individually, the nine-gene signature demonstrated a higher prognostic accuracy (Figure 5C and 5D).

Time-dependent ROC curves were generated to evaluate the nine-gene signature performance. (A, B) Three-years or five-years ROC curves of the nine-gene signature in the training cohort and validation cohort. © Five-years ROC curves for nine-gene signature and single gene. (D) Five-years ROC curves for nine-gene signature and clinical risk factor.

Figure 5. Time-dependent ROC curves were generated to evaluate the nine-gene signature performance. (A, B) Three-years or five-years ROC curves of the nine-gene signature in the training cohort and validation cohort. © Five-years ROC curves for nine-gene signature and single gene. (D) Five-years ROC curves for nine-gene signature and clinical risk factor.

Comparing the immune infiltration between the high- and low-risk groups

To provide novel insight into the biological role of each of the risk groups, we performed immune infiltration analysis using the GSVA method. Among the low-risk group, the level of immune infiltration (e.g., “Activated B cell”, “Activated CD4 cell”, “Activated CD8 cell”, “Effector memory CD8 T cell”, and “Immature B cell”) were found to be significantly higher than that of the high-risk group (Figure 6A). In contrast, “Central memory CD8 T cells”, “Immature dendritic cells”, and “Plasmacytoid dendritic cells” were significantly enriched in the high-risk group (Figure 6A).

The relationship between the nine-gene signature and immune infiltration. (A) Comparison of relative immune cell abundance based on GSVA score in high-risk and low-risk groups (B) Partial Spearman's correlation of nine genes expression and immune infiltrates. *: Statistically significant p

Figure 6. The relationship between the nine-gene signature and immune infiltration. (A) Comparison of relative immune cell abundance based on GSVA score in high-risk and low-risk groups (B) Partial Spearman's correlation of nine genes expression and immune infiltrates. *: Statistically significant p < 0.05, **: Statistically significant p < 0.01.

Since tumor-infiltrating lymphocytes are an independent prognostic predictor of survival in various tumors [28–30], we performed correlation analysis between the nine genes and immune infiltration level for OC. The results showed that the association between the nine genes and the immune microenvironment is significant, as each of the nine genes had a significant correlation with tumor purity (Figure 6B). Among these genes, D4S234E was positively correlated with tumor purity, whereas the other eight genes were negatively correlated with tumor purity. The most relevant genes among the nine gene signature associated with immune infiltration included: CX3CR1 (related B cell, cor = 0.311), GBP2 (related CD8⁺ T cell, dendritic cell, cor = 0.403, 0.495), CXCL13 (related CD4+ T cell, cor = 0.308), and VSIG4 (related macrophage and neutrophil cor = 0.51, 0.605).

Stratification analysis based on clinical information

Risk stratification analysis was performed to test whether the nine-gene signature could predict OS regardless of tumor residual disease. The results shown that patients with low-risk scores had significantly longer OS than patients with high-risk scores in no macroscopic disease (p=0.0014, 1-10 mm (p =0.00026), >20 mm (p=0.0001)) (Supplementary Figure 1).

Discussion

OC represents one of the diseases with the highest mortality rate of the female reproductive system [1]. Due to the lack of early and effective detection methods, most OC patients were diagnosed at an advanced stage who subsequently missed the optimal treatment period and resulting in a poor clinical outcome. Recent studies have shown that TME played a vital role during OC progression [31–33]. Moreover, TME-related molecular markers can be used as predictors to precisely assess patients’ immunotherapy response, thereby enhancing their clinical outcome [34–36]. However, immune infiltration and its molecular mechanisms have not been thoroughly explained in OC.

To our knowledge, our work is the first to use the ESTIMATE algorithm combined with LASSO-Cox to explore molecular markers associated with OC prognosis. Firstly, we derived a series of TME-associated DEGs by comparing the transcriptional expression profiles in 593 OC patients with high versus low stromal/immune scores based on TCGA data. The DEGs functional enrichment analysis indicated that the main pathways were associated with immune response and cancer (e.g., cytokine-cytokine receptor interaction, Human immunodeficiency virus 1 infection, focal adhesion, proteoglycans in cancer, and PI3K-Akt signaling pathway), which are in agreement with findings that immune response and cancer progression exhibit crosstalk and interact with each other [37]. Based on the obtained DEGs, we built a nine-gene signature that was notably related to the OS in OC patients in both the training and validation cohorts. The patients could be divided into high-risk and low-risk groups with distinct differences in the five-year OS with this nine-gene signature. The results of the GSVA analysis showed that “Activated CD8 T cells”, “Effector memory CD8 T cell”, “Activated B cells”, and “Activated CD4 T cells” were associated with significantly lower infiltration in the high-risk group. In contrast, “Central memory CD8 T cells”, “Immature dendritic cells”, and “Plasmacytoid dendritic cells” were associated with significantly higher infiltration in the high-risk group.

Numerous studies have documented that CD8⁺ T cell high infiltration in the TME is associated with positive anti-tumor effects in various cancer [38–40]. Natural killer (NK)-dendritic cell (DC) cross talk results in upregulation of Chemokine (C-X-C motif) ligand 9 (CXCL9), Chemokine (C-X-C motif) ligand 10 (CXCL10), and Chemokine (C-C motif) ligand 5 (CCL5) on DCs leading to CD8⁺ effector T cells recruitment into TME, thereby promote antitumor immune response in OC [41]. On the other hand, the expression of inhibitory molecules such as CTLA4, PD-1, and LAG3 on CD8⁺ T cells are promoted by IL-6 and IL-10, that produced by tumor cells and tumor-associated macrophages, in turn inhibit CD8⁺ T cells infiltration [42–44]. Another example of CD8⁺ T cells inhibition shown C-C motif chemokine 22 (CCL22) can promote CTLA4⁺ FOXP3⁺ GITR⁺ Tregs and CCR4+ Tregs infiltration in TME, thereby inhibiting CD8⁺ T cells activation [45, 46].

Curdin et al. reported that plasmacytoid DCs can induce immunosuppression in OC by providing ICOS⁺ Treg cells with Inducible T Cell Costimulator Ligand (ICOS-L) stimulation, thereby enhancing the capability of Treg cell impairments of T cell proliferation [47]. Immature DCs express low levels of Major histocompatibility complex (MHC) and co-stimulatory molecules, therefore T cells activation by immature DCs is inefficient [48, 49]. These results suggest that a low density of activated T cells or high infiltration of immature/plasmacytoid DCs may be the cause of poor clinical outcome of cancer patients.

Tumor progression and metastasis typically occur in adipose tissue-rich areas such as omentum, one of the main metastatic sites in OC. Adipose tissue is composed of a variety of cells including adipocytes, adipose stem cells, endothelial cells, and infiltrating immune cells that secrete diverse soluble tumor-promoting factors such as hormones, cytokines, reactive oxygen species, extracellular matrix, and lipid metabolites. These secreted factors not only directly promote tumor progression but can also reduce the anti-tumor immune response by altering the TME. For example, chemokines (e.g., TNF-α, IL-6, and IL-1b) recruit immunosuppressive neutrophils and M2 macrophages to the TME, thereby inhibiting anti-tumor cell activity (e.g., TCD 8+ lymphocytes and NKTs) [50]. Besides, adipose tissue associated PD-L1 is found to attenuate T cell activation which also contributes to an immune suppressive microenvironment [51]. A large number of adipose tissue infiltrating M1 macrophages can lead to adipocyte death. Moreover, the release of intracellular substances of dead adipocyte not only aggravates inflammation but also provides energy required for tumor cell growth. These factors all provide a favorable microenvironment for tumor growth [52].

Cancer is a heterogeneous disease for which the identification of dysregulated genes involved in tumorigenesis and progression might aid in improving prognostic and treatment strategies. In this study, we identified a group of nine genes (CX3CR1, UBD, GBP2, D4S234E, CXCL11, CXCL13, VSIG4, TFPI2 and C5AR1) that can effectively predict the OS in OC. Among these genes, CX3CR1 and UBD can promote tumor metastasis and the epithelial to mesenchymal transition [53–55]. As p53-related genes, GBP2 and D4S234E have been previously shown to regulate mitochondrial fission and apoptosis of cancer cells [56–58]. CXCL11 and CXCL13 are associated with CD8⁺ T cell and B cell infiltration, which act as a protective factor inhibiting tumor [59, 60]. VSIG4 inhibits T cell proliferation and IL-2 production as well as regulates Treg differentiation and stability leading to immune tolerance [61, 62]. TFPI2 is a serum biomarker for the detection of ovarian clear cell adenocarcinoma, and its predictive values are superior to that of CA125 [63]. Great promise for immunotherapies has been achieved by the Conduct Phase II and III clinical trials for the discovery of a drug that targets the C5a-C5aR1 pathway [64]. Therefore, our nine-gene signature could potentially be used as a predictive tool for risk assessment and might offer potential targets for immunotherapy in the clinical management of OC.

There are also some limitations associated with this research that should be addressed. First, the biological function of the nine identified genes should be validated in wet lab experiments, particularly regarding the association with immune infiltration. Second, missing information in OC clinical characteristics (contains many patients with “unknown” information in venous invasion and lymphatic invasion) in TCGA limited us in building a nomogram for incorporating clinical characteristics to improve the predicted accuracy of the model. Third, the risk score model requires further validation in multiple cohorts to evaluate the model generalization ability.

In conclusion, the gene expression profile and clinical characteristics of the TCGA database were analyzed by ESTIMATE and a Lasso-Cox algorithm to obtain the nine gene prognostic signature related to TME in OC. This molecular signature can effectively distinguish high-risk populations from OC patients in TCGA and GSE17260 datasets. In addition, the expression of each gene in the model is significantly correlated with the TME components, which further supports the important role of the TME in the occurrence and development of OC.

Materials and Methods

Data source

The available Level 3 gene expression profiles of the OC patients were downloaded from the TCGA database (https://tcga-data.nci.nih.gov/tcga/). RNA expression detection of ovarian serous adenocarcinoma was performed on Affymetrix HT Human Genome U133a microarray. Clinicopathological characteristics including age, histological type, FIGO stage, venous invasion, lymphatic invasion, tumor residual disease, survival, and outcome were also downloaded from the TCGA data coordination center. ESTIMATE algorithm was used to calculate immune and stromal scores using "estimate" package (https://r-forge.r-project.org; repos= rforge, dependencies=TRUE) [24].

There were 593 ovarian cancer samples with gene expression data in the TCGA. The data from 576 primary solid tumors patients were retained after removing 17 recurrent solid tumors samples. The samples were further trimmed to 465 to include only patients with corresponding clinical information. The distribution of the values for the samples we have selected was viewed graphically as a box plot and no outliers were identified. All 465 serous adenocarcinoma samples contained complete survival information (survival time and outcome). All 12,042 coding genes of each patient did not contain missing values. Gene expression values were normalized by log2 (Affy RMA). The 465 OC patients were grouped as a training cohort. To confirm whether our nine-gene signature could predict prognosis of OC patients in an independent dataset, we selected to use publicly available data (GSE17260) with clinical informations such as OS time, samples were diagnosed as serous adenocarcinoma, obtained from primary lesion, and gene expression profiling was measured by microarray. The raw data for the validation set GSE17260 was downloaded from the GEO database, which include 109 OC samples (a low quality sample was removed) and its corresponding clinical information (Grade, FIGO stage). R package limma was used for quality control and normalization. Gene expression value for genes with multiple probes was calculated as the average of the probes.

Identification of prognosis-related genes

We divided the patients into two groups based on the median value of the immune scores and stromal scores. The differential expression analysis was performed with R package limma [65], and the fold-change (> 1.5) and an eBayes test p-value (< 0.05) were used as criteria to screen for DEGs between the high and low groups. We subsequently performed a functional enrichment analysis with the R package “clusterprofiler” [66] to identify the potential biological function of the DEGs. The intersection of immune and stromal DEGs was subjected to a univariate Cox regression analysis to identify the OS-related signature. A threshold of p < 0.05 was deemed significant.

Construction of risk score system of OC

We performed the least absolute shrinkage and selection operator (LASSO) on the Cox regression model using R package glmnet [67]. The 10-fold cross-validation approach and “one-standard error (1se)” was used to identify the optimum parameter λ. After the Lasso-Cox analysis, we obtained the corresponding regression coefficients and signatures of the nine genes and constructed the following formula:

Risk score = sum of each gene^’s (regression coefficients × level of gene expression)

Each patient was assigned a risk score according to this formula. The optimal cut-off of the risk score was determined by the “surv_cutpoint” function of the “survminer” R package (https://www.r-project.org/) and used to stratify the OC patients into high- and low-risk groups. A comparison of the survival between the two groups was analyzed with a Kaplan-Meier estimator and log-rank test.

Time-dependent ROC curve analysis

A time-dependent ROC Curve method can be implemented to estimate the three- and five-years prognostic model prediction performance in a training and validation cohort [68]. A stratified analysis was conducted to investigate whether the prognostic model was widely applicable to clinical characteristics. Furthermore, the AUC was used to determine whether the prognostic model was superior to that of other risk factors.

The correlation between gene expression and immune infiltration

Gene sets for 28 subpopulations of tumor-infiltrating lymphocytes resulting from the study by Charoentong et al. [69], containing cell types related to adaptive immunity (activated, central memory, effector memory CD4⁺ and CD8⁺ T cells, γδ T cells, type 1 helper T (TH1) cells, TH2 cells, TH17 cells, regulatory T cells, follicular helper T cells, as well as activated immature and memory B cells) and innate immunity (macrophages; monocytes; mast cells; eosinophils; neutrophils; activated, plasmacytoid and immature dendritic cells; natural killer cells; natural killer T cells, and myeloid-derived suppressor cells). The gene set parameter was subjected to GSVA analysis. GSVA transformed gene expression into an absolute enrichment score, which was represented as relative immune cell abundance in each sample [70]. A t-test was performed to compare the GSVA score between the high- and low-risk groups. The Tumor Immune Estimation Resource (TIMER, https://cistrome.shinyapps.io/timer/) was used to investigate the correlation between gene expression and tumor-infiltrating immune cells, including B cells, CD4⁺ T cells, CD8⁺ T cells, neutrophils, macrophages, and dendritic cells in OC [71]. A heat map was generated with partial Spearman's correlation and p < 0.05 was regarded as statistically significant.

Supplementary Materials

Supplementary Figure 1

Author Contributions

Q.D., S.D. and R.W. designed the study and performed data analysis; Y.Z., Y.L., X.Z, and H.R. revised the paper; S.Z. and W.W. performed literature search and data collection; K.W. constructed figures and improved the language; H.W., J.W., and K.Z. directed the overall project. All authors reviewed the manuscript.

Acknowledgments

We are very grateful to TCGA and GEO databases for providing valuable data resources to enable us to conduct this research.

Conflicts of Interest

The authors declare that there is no conflicts of interest to disclose.

Funding

National Natural Science Foundation (No. 81201730 and No. 81703005), Key Research and Development Project of Hunan Science and Technology Department (No. 2018SK2126), Key Project of Changsha Science and Technology Bureau (No. kq1706046), Research Project of National Cancer Center Cancer (No.NCC2017A21), Research Project of Health and Family Planning Commission of Hunan Province (No. B2017098), The Natural Science Foundation of Hunan Province (No. 12JJ5073, No. 2017JJ3195, and No. 2018JJ3311), The Provincial Key Research and Development Program of Hunan Province(2018SK2123), The Provincial Key Clinical Specialist Construction Projects of Hunan Province.

References

1. Menon U, Karpinskyj C, Gentry-Maharaj A. Ovarian Cancer Prevention and Screening. Obstet Gynecol. 2018; 131:909–27. https://doi.org/10.1097/AOG.0000000000002580 [PubMed]
2. Reid BM, Permuth JB, Sellers TA. Epidemiology of ovarian cancer: a review. Cancer Biol Med. 2017; 14:9–32. https://doi.org/10.20892/j.issn.2095-3941.2016.0084 [PubMed]
3. Jayson GC, Kohn EC, Kitchener HC, Ledermann JA. Ovarian cancer. Lancet. 2014; 384:1376–88. https://doi.org/10.1016/S0140-6736(13)62146-7 [PubMed]
4. Rosen DG, Yang G, Liu G, Mercado-Uribe I, Chang B, Xiao XS, Zheng J, Xue FX, Liu J. Ovarian cancer: pathology, biology, and disease models. Front Biosci. 2009; 14:2089–102. https://doi.org/10.2741/3364 [PubMed]
5. Prat J, and FIGO Committee on Gynecologic Oncology. Staging classification for cancer of the ovary, fallopian tube, and peritoneum. Int J Gynaecol Obstet. 2014; 124:1–5. https://doi.org/10.1016/j.ijgo.2013.10.001 [PubMed]
6. McGee J, Bookman M, Harter P, Marth C, McNeish I, Moore KN, Poveda A, Hilpert F, Hasegawa K, Bacon M, Gatsonis C, Brand A, Kridelka F, et al, and participants of the 5th Ovarian Cancer Consensus Conference. Fifth Ovarian Cancer Consensus Conference: individualized therapy and patient factors. Ann Oncol. 2017; 28:702–10. https://doi.org/10.1093/annonc/mdx010 [PubMed]
7. El Bairi K, Kandhro AH, Gouri A, Mahfoud W, Louanjli N, Saadani B, Afqir S, Amrani M. Emerging diagnostic, prognostic and therapeutic biomarkers for ovarian cancer. Cell Oncol (Dordr). 2017; 40:105–18. https://doi.org/10.1007/s13402-016-0309-1 [PubMed]
8. Tucker DW, Getchell CR, McCarthy ET, Ohman AW, Sasamoto N, Xu S, Ko JY, Gupta M, Shafrir A, Medina JE, Lee JJ, MacDonald LA, Malik A, et al. Epigenetic Reprogramming Strategies to Reverse Global Loss of 5-Hydroxymethylcytosine, a Prognostic Factor for Poor Survival in High-grade Serous Ovarian Cancer. Clin Cancer Res. 2018; 24:1389–401. https://doi.org/10.1158/1078-0432.CCR-17-1958 [PubMed]
9. Wang J, Dean DC, Hornicek FJ, Shi H, Duan Z. Cyclin-dependent kinase 9 (CDK9) is a novel prognostic marker and therapeutic target in ovarian cancer. FASEB J. 2019; 33:5990–6000. https://doi.org/10.1096/fj.201801789RR [PubMed]
10. Watson ZL, Yamamoto TM, McMellen A, Kim H, Hughes CJ, Wheeler LJ, Post MD, Behbakht K, Bitler BG. Histone methyltransferases EHMT1 and EHMT2 (GLP/G9A) maintain PARP inhibitor resistance in high-grade serous ovarian carcinoma. Clin Epigenetics. 2019; 11:165. https://doi.org/10.1186/s13148-019-0758-2 [PubMed]
11. Chen LY, Huang RL, Chan MW, Yan PS, Huang TS, Wu RC, Suryo Rahmanto Y, Su PH, Weng YC, Chou JL, Chao TK, Wang YC, Shih IM, Lai HC. TET1 reprograms the epithelial ovarian cancer epigenome and reveals casein kinase 2α as a therapeutic target. J Pathol. 2019; 248:363–76. https://doi.org/10.1002/path.5266 [PubMed]
12. Tan M, Asad M, Heong V, Wong MK, Tan TZ, Ye J, Kuay KT, Thiery JP, Scott C, Huang RY. The FZD7-TWIST1 axis is responsible for anoikis resistance and tumorigenesis in ovarian carcinoma. Mol Oncol. 2019; 13:757–780. https://doi.org/10.1002/1878-0261.12425 [PubMed]
13. Zeng M, Kwiatkowski NP, Zhang T, Nabet B, Xu M, Liang Y, Quan C, Wang J, Hao M, Palakurthi S, Zhou S, Zeng Q, Kirschmeier PT, et al. Targeting MYC dependency in ovarian cancer through inhibition of CDK7 and CDK12/13. eLife. 2018; 7:7. https://doi.org/10.7554/eLife.39030 [PubMed]
14. Park SM, Park SH, Ryu KJ, Kim IK, Han H, Kim HJ, Kim SH, Hong KS, Kim H, Kim M, Cho BI, Heo JD, Kim NH, et al. Downregulation of CHIP promotes ovarian cancer metastasis by inducing Snail-mediated epithelial-mesenchymal transition. Mol Oncol. 2019; 13:1280–95. https://doi.org/10.1002/1878-0261.12485 [PubMed]
15. Ma Y, Zhang H, Xiong C, Liu Z, Xu Q, Feng J, Zhang J, Wang Z, Yan X. CD146 mediates an E-cadherin-to-N-cadherin switch during TGF-β signaling-induced epithelial-mesenchymal transition. Cancer Lett. 2018; 430:201–14. https://doi.org/10.1016/j.canlet.2018.05.016 [PubMed]
16. Grivennikov SI, Greten FR, Karin M. Immunity, inflammation, and cancer. Cell. 2010; 140:883–99. https://doi.org/10.1016/j.cell.2010.01.025 [PubMed]
17. Klemm F, Joyce JA. Microenvironmental regulation of therapeutic response in cancer. Trends Cell Biol. 2015; 25:198–213. https://doi.org/10.1016/j.tcb.2014.11.006 [PubMed]
18. Swartz MA, Iida N, Roberts EW, Sangaletti S, Wong MH, Yull FE, Coussens LM, DeClerck YA. Tumor microenvironment complexity: emerging roles in cancer therapy. Cancer Res. 2012; 72:2473–80. https://doi.org/10.1158/0008-5472.CAN-12-0122 [PubMed]
19. Zhang Y, Tang H, Cai J, Zhang T, Guo J, Feng D, Wang Z. Ovarian cancer-associated fibroblasts contribute to epithelial ovarian carcinoma metastasis by promoting angiogenesis, lymphangiogenesis and tumor cell invasion. Cancer Lett. 2011; 303:47–55. https://doi.org/10.1016/j.canlet.2011.01.011 [PubMed]
20. Leung CS, Yeung TL, Yip KP, Wong KK, Ho SY, Mangala LS, Sood AK, Lopez-Berestein G, Sheng J, Wong ST, Birrer MJ, Mok SC. Cancer-associated fibroblasts regulate endothelial adhesion protein LPP to promote ovarian cancer chemoresistance. J Clin Invest. 2018; 128:589–606. https://doi.org/10.1172/JCI95200 [PubMed]
21. Yokoi A, Yoshioka Y, Yamamoto Y, Ishikawa M, Ikeda SI, Kato T, Kiyono T, Takeshita F, Kajiyama H, Kikkawa F, Ochiya T. Malignant extracellular vesicles carrying MMP1 mRNA facilitate peritoneal dissemination in ovarian cancer. Nat Commun. 2017; 8:14470. https://doi.org/10.1038/ncomms14470 [PubMed]
22. Zheng M, Hu Y, Gou R, Wang J, Nie X, Li X, Liu Q, Liu J, Lin B. Integrated multi-omics analysis of genomics, epigenomics, and transcriptomics in ovarian carcinoma. Aging (Albany NY). 2019; 11:4198–215. https://doi.org/10.18632/aging.102047 [PubMed]
23. Wu J, Xu WH, Wei Y, Qu YY, Zhang HL, Ye DW. An Integrated Score and Nomogram Combining Clinical and Immunohistochemistry Factors to Predict High ISUP Grade Clear Cell Renal Cell Carcinoma. Front Oncol. 2018; 8:634. https://doi.org/10.3389/fonc.2018.00634 [PubMed]
24. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, Treviño V, Shen H, Laird PW, Levine DA, Carter SL, Getz G, Stemke-Hale K, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013; 4:2612. https://doi.org/10.1038/ncomms3612 [PubMed]
25. Shah N, Wang P, Wongvipat J, Karthaus WR, Abida W, Armenia J, Rockowitz S, Drier Y, Bernstein BE, Long HW, Freedman ML, Arora VK, Zheng D, Sawyers CL. Regulation of the glucocorticoid receptor via a BET-dependent enhancer drives antiandrogen resistance in prostate cancer. eLife. 2017; 6:6. https://doi.org/10.7554/eLife.27861 [PubMed]
26. Jia D, Li S, Li D, Xue H, Yang D, Liu Y. Mining TCGA database for genes of prognostic value in glioblastoma microenvironment. Aging (Albany NY). 2018; 10:592–605. https://doi.org/10.18632/aging.101415 [PubMed]
27. Xu WH, Xu Y, Wang J, Wan FN, Wang HK, Cao DL, Shi GH, Qu YY, Zhang HL, Ye DW. Prognostic value and immune infiltration of novel signatures in clear cell renal cell carcinoma microenvironment. Aging (Albany NY). 2019; 11:6999–7020. https://doi.org/10.18632/aging.102233 [PubMed]
28. Fu Q, Chen N, Ge C, Li R, Li Z, Zeng B, Li C, Wang Y, Xue Y, Song X, Li H, Li G. Prognostic value of tumor-infiltrating lymphocytes in melanoma: a systematic review and meta-analysis. Oncoimmunology. 2019; 8:1593806. https://doi.org/10.1080/2162402X.2019.1593806 [PubMed]
29. Ohtani H. Focus on TILs: prognostic significance of tumor infiltrating lymphocytes in human colorectal cancer. Cancer Immun. 2007; 7:4. [PubMed]
30. Azimi F, Scolyer RA, Rumcheva P, Moncrieff M, Murali R, McCarthy SW, Saw RP, Thompson JF. Tumor-infiltrating lymphocyte grade is an independent predictor of sentinel lymph node status and survival in patients with cutaneous melanoma. J Clin Oncol. 2012; 30:2678–83. https://doi.org/10.1200/JCO.2011.37.8539 [PubMed]
31. Singel KL, Emmons TR, Khan AN, Mayor PC, Shen S, Wong JT, Morrell K, Eng KH, Mark J, Bankert RB, Matsuzaki J, Koya RC, Blom AM, et al. Mature neutrophils suppress T cell immunity in ovarian cancer microenvironment. JCI Insight. 2019; 4:4. https://doi.org/10.1172/jci.insight.122311 [PubMed]
32. Curtis M, Mukherjee A, Lengyel E. The Tumor Microenvironment Takes Center Stage in Ovarian Cancer Metastasis. Trends Cancer. 2018; 4:517–19. https://doi.org/10.1016/j.trecan.2018.06.002 [PubMed]
33. Au K, Peterson N, Truesdell P, Reid-Schachter G, Khalaj K, Ren R, Francis JA, Graham CH, Craig AW, Koti M. CXCL10 alters the tumour immune microenvironment and disease progression in a syngeneic murine model of high-grade serous ovarian cancer. Gynecol Oncol. 2017; 145:436–45. https://doi.org/10.1016/j.ygyno.2017.03.007 [PubMed]
34. Isomoto K, Haratani K, Hayashi H, Shimizu S, Tomida S, Niwa T, Yokoyama T, Fukuda Y, Chiba Y, Kato R, Tanizaki J, Tanaka K, Takeda M, et al. Impact of EGFR-TKI treatment on the tumor immune microenvironment in EGFR mutation-positive non-small cell lung cancer. Clin Cancer Res. 2020. [Epub ahead of print]. https://doi.org/10.1158/1078-0432.CCR-19-2027 [PubMed]
35. Ma Z, Shuai Y, Gao X, Wen X, Ji J. Circular RNAs in the tumour microenvironment. Mol Cancer. 2020; 19:8. https://doi.org/10.1186/s12943-019-1113-0 [PubMed]
36. Cai D, Li J, Liu D, Hong S, Qiao Q, Sun Q, Li P, Lyu N, Sun T, Xie S, Guo L, Ni L, Jin L, Dong C. Tumor-expressed B7-H3 mediates the inhibition of antitumor T-cell functions in ovarian cancer insensitive to PD-1 blockade therapy. Cell Mol Immunol. 2020; 17:227–36. https://doi.org/10.1038/s41423-019-0305-2 [PubMed]
37. Yu YR, Ho PC. Sculpting tumor microenvironment with immune system: from immunometabolism to immunoediting. Clin Exp Immunol. 2019; 197:153–60. https://doi.org/10.1111/cei.13293 [PubMed]
38. Piersma SJ, Jordanova ES, van Poelgeest MI, Kwappenberg KM, van der Hulst JM, Drijfhout JW, Melief CJ, Kenter GG, Fleuren GJ, Offringa R, van der Burg SH. High number of intraepithelial CD8+ tumor-infiltrating lymphocytes is associated with the absence of lymph node metastases in patients with large early-stage cervical cancer. Cancer Res. 2007; 67:354–61. https://doi.org/10.1158/0008-5472.CAN-06-3388 [PubMed]
39. Kim PS, Ahmed R. Features of responding T cells in cancer and chronic infection. Curr Opin Immunol. 2010; 22:223–30. https://doi.org/10.1016/j.coi.2010.02.005 [PubMed]
40. Kmiecik J, Poli A, Brons NH, Waha A, Eide GE, Enger PO, Zimmer J, Chekenya M. Elevated CD3+ and CD8+ tumor-infiltrating immune cells correlate with prolonged survival in glioblastoma patients despite integrated immunosuppressive mechanisms in the tumor microenvironment and at the systemic level. J Neuroimmunol. 2013; 264:71–83. https://doi.org/10.1016/j.jneuroim.2013.08.013 [PubMed]
41. Wong JL, Berk E, Edwards RP, Kalinski P. IL-18-primed helper NK cells collaborate with dendritic cells to promote recruitment of effector CD8+ T cells to the tumor microenvironment. Cancer Res. 2013; 73:4653–62. https://doi.org/10.1158/0008-5472.CAN-12-4366 [PubMed]
42. Hodi FS, Butler M, Oble DA, Seiden MV, Haluska FG, Kruse A, Macrae S, Nelson M, Canning C, Lowy I, Korman A, Lautz D, Russell S, et al. Immunologic and clinical effects of antibody blockade of cytotoxic T lymphocyte-associated antigen 4 in previously vaccinated cancer patients. Proc Natl Acad Sci USA. 2008; 105:3005–10. https://doi.org/10.1073/pnas.0712237105 [PubMed]
43. Matsuzaki J, Gnjatic S, Mhawech-Fauceglia P, Beck A, Miller A, Tsuji T, Eppolito C, Qian F, Lele S, Shrikant P, Old LJ, Odunsi K. Tumor-infiltrating NY-ESO-1-specific CD8+ T cells are negatively regulated by LAG-3 and PD-1 in human ovarian cancer. Proc Natl Acad Sci USA. 2010; 107:7875–80. https://doi.org/10.1073/pnas.1003345107 [PubMed]
44. Huang RY, Eppolito C, Lele S, Shrikant P, Matsuzaki J, Odunsi K. LAG3 and PD1 co-inhibitory molecules collaborate to limit CD8+ T cell signaling and dampen antitumor immunity in a murine ovarian cancer model. Oncotarget. 2015; 6:27359–77. https://doi.org/10.18632/oncotarget.4751 [PubMed]
45. Landskron J, Helland Ø, Torgersen KM, Aandahl EM, Gjertsen BT, Bjørge L, Taskén K. Activated regulatory and memory T-cells accumulate in malignant ascites from ovarian carcinoma patients. Cancer Immunol Immunother. 2015; 64:337–47. https://doi.org/10.1007/s00262-014-1636-6 [PubMed]
46. Curiel TJ, Coukos G, Zou L, Alvarez X, Cheng P, Mottram P, Evdemon-Hogan M, Conejo-Garcia JR, Zhang L, Burow M, Zhu Y, Wei S, Kryczek I, et al. Specific recruitment of regulatory T cells in ovarian carcinoma fosters immune privilege and predicts reduced survival. Nat Med. 2004; 10:942–49. https://doi.org/10.1038/nm1093 [PubMed]
47. Conrad C, Gregorio J, Wang YH, Ito T, Meller S, Hanabuchi S, Anderson S, Atkinson N, Ramirez PT, Liu YJ, Freedman R, Gilliet M. Plasmacytoid dendritic cells promote immunosuppression in ovarian cancer via ICOS costimulation of Foxp3(+) T-regulatory cells. Cancer Res. 2012; 72:5240–49. https://doi.org/10.1158/0008-5472.CAN-12-2271 [PubMed]
48. Manicassamy S, Pulendran B. Dendritic cell control of tolerogenic responses. Immunol Rev. 2011; 241:206–27. https://doi.org/10.1111/j.1600-065X.2011.01015.x [PubMed]
49. Lutz MB, Schuler G. Immature, semi-mature and fully mature dendritic cells: which signals induce tolerance or immunity? Trends Immunol. 2002; 23:445–49. https://doi.org/10.1016/S1471-4906(02)02281-0 [PubMed]
50. Xu YX, Mishra S. Obesity-Linked Cancers: Current Knowledge, Challenges and Limitations in Mechanistic Studies and Rodent Models. Cancers (Basel). 2018; 10:10. https://doi.org/10.3390/cancers10120523 [PubMed]
51. Wu B, Sun X, Gupta HB, Yuan B, Li J, Ge F, Chiang HC, Zhang X, Zhang C, Zhang D, Yang J, Hu Y, Curiel TJ, Li R. Adipose PD-L1 Modulates PD-1/PD-L1 Checkpoint Blockade Immunotherapy Efficacy in Breast Cancer. Oncoimmunology. 2018; 7:e1500107. https://doi.org/10.1080/2162402X.2018.1500107 [PubMed]
52. Corrêa LH, Heyn GS, Magalhaes KG. The Impact of the Adipose Organ Plasticity on Inflammation and Cancer Progression. Cells. 2019; 8:8. https://doi.org/10.3390/cells8070662 [PubMed]
53. Singh SK, Mishra MK, Singh R. Hypoxia-inducible factor-1α induces CX3CR1 expression and promotes the epithelial to mesenchymal transition (EMT) in ovarian cancer cells. J Ovarian Res. 2019; 12:42. https://doi.org/10.1186/s13048-019-0517-1 [PubMed]
54. Shen F, Zhang Y, Jernigan DL, Feng X, Yan J, Garcia FU, Meucci O, Salvino JM, Fatatis A. Novel Small-Molecule CX3CR1 Antagonist Impairs Metastatic Seeding and Colonization of Breast Cancer Cells. Mol Cancer Res. 2016; 14:518–27. https://doi.org/10.1158/1541-7786.MCR-16-0013 [PubMed]
55. Yuan R, Wang K, Hu J, Yan C, Li M, Yu X, Liu X, Lei J, Guo W, Wu L, Hong K, Shao J. Ubiquitin-like protein FAT10 promotes the invasion and metastasis of hepatocellular carcinoma by modifying β-catenin degradation. Cancer Res. 2014; 74:5287–300. https://doi.org/10.1158/0008-5472.CAN-14-0284 [PubMed]
56. Zhang J, Zhang Y, Wu W, Wang F, Liu X, Shui G, Nie C. Guanylate-binding protein 2 regulates Drp1-mediated mitochondrial fission to suppress breast cancer cell invasion. Cell Death Dis. 2017; 8:e3151. https://doi.org/10.1038/cddis.2017.559 [PubMed]
57. Guimarães DP, Oliveira IM, de Moraes E, Paiva GR, Souza DM, Barnas C, Olmedo DB, Pinto CE, Faria PA, De Moura Gallo CV, Small IA, Ferreira CG, Hainaut P. Interferon-inducible guanylate binding protein (GBP)-2: a novel p53-regulated tumor marker in esophageal squamous cell carcinomas. Int J Cancer. 2009; 124:272–79. https://doi.org/10.1002/ijc.23944 [PubMed]
58. Kudoh T, Kimura J, Lu ZG, Miki Y, Yoshida K. D4S234E, a novel p53-responsive gene, induces apoptosis in response to DNA damage. Exp Cell Res. 2010; 316:2849–58. https://doi.org/10.1016/j.yexcr.2010.06.025 [PubMed]
59. Gao Q, Wang S, Chen X, Cheng S, Zhang Z, Li F, Huang L, Yang Y, Zhou B, Yue D, Wang D, Cao L, Maimela NR, et al. Cancer-cell-secreted CXCL11 promoted CD8⁺ T cells infiltration through docetaxel-induced-release of HMGB1 in NSCLC. J Immunother Cancer. 2019; 7:42. https://doi.org/10.1186/s40425-019-0511-6 [PubMed]
60. Bindea G, Mlecnik B, Tosolini M, Kirilovsky A, Waldner M, Obenauf AC, Angell H, Fredriksen T, Lafontaine L, Berger A, Bruneval P, Fridman WH, Becker C, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity. 2013; 39:782–95. https://doi.org/10.1016/j.immuni.2013.10.003 [PubMed]
61. Vogt L, Schmitz N, Kurrer MO, Bauer M, Hinton HI, Behnke S, Gatto D, Sebbel P, Beerli RR, Sonderegger I, Kopf M, Saudan P, Bachmann MF. VSIG4, a B7 family-related protein, is a negative regulator of T cell activation. J Clin Invest. 2006; 116:2817–26. https://doi.org/10.1172/JCI25673 [PubMed]
62. Yuan X, Yang BH, Dong Y, Yamamura A, Fu W. CRIg, a tissue-resident macrophage specific immune checkpoint molecule, promotes immunological tolerance in NOD mice, via a dual role in effector and regulatory T cells. eLife. 2017; 6:6. https://doi.org/10.7554/eLife.29540 [PubMed]
63. Arakawa N, Miyagi E, Nomura A, Morita E, Ino Y, Ohtake N, Miyagi Y, Hirahara F, Hirano H. Secretome-based identification of TFPI2, a novel serum biomarker for detection of ovarian clear cell adenocarcinoma. J Proteome Res. 2013; 12:4340–50. https://doi.org/10.1021/pr400282j [PubMed]
64. Hawksworth OA, Li XX, Coulthard LG, Wolvetang EJ, Woodruff TM. New concepts on the therapeutic control of complement anaphylatoxin receptors. Mol Immunol. 2017; 89:36–43. https://doi.org/10.1016/j.molimm.2017.05.015 [PubMed]
65. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015; 43:e47. https://doi.org/10.1093/nar/gkv007 [PubMed]
66. Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012; 16:284–87. https://doi.org/10.1089/omi.2011.0118 [PubMed]
67. Friedman J, Hastie T, Tibshirani R. Regularization Paths for Generalized Linear Models via Coordinate Descent. J Stat Softw. 2010; 33:1–22. https://doi.org/10.18637/jss.v033.i01 [PubMed]
68. Kamarudin AN, Cox T, Kolamunnage-Dona R. Time-dependent ROC curve analysis in medical research: current methods and applications. BMC Med Res Methodol. 2017; 17:53. https://doi.org/10.1186/s12874-017-0332-6 [PubMed]
69. Charoentong P, Finotello F, Angelova M, Mayer C, Efremova M, Rieder D, Hackl H, Trajanoski Z. Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep. 2017; 18:248–62. https://doi.org/10.1016/j.celrep.2016.12.019 [PubMed]
70. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics. 2013; 14:7. https://doi.org/10.1186/1471-2105-14-7 [PubMed]
71. Li B, Severson E, Pignon JC, Zhao H, Li T, Novak J, Jiang P, Shen H, Aster JC, Rodig S, Signoretti S, Liu JS, Liu XS. Comprehensive analyses of tumor immunity: implications for cancer immunotherapy. Genome Biol. 2016; 17:174. https://doi.org/10.1186/s13059-016-1028-7 [PubMed]

Research Paper Volume 12, Issue 6 pp 4879—4895