Identification of MMP1 as a Potential Prognostic Biomarker and Correlating with Immune Infiltrates in Cervical Squamous Cell Carcinoma

Abstract

Cervical squamous cell carcinoma (CESC) is a human papillomavirus-driven tumor that the underlying molecular mechanisms are unclear. This study aimed to elucidate the key candidate genes and potential mechanism in CESC by bioinformatics analysis. A total of 132 common differentially expressed genes (DEGs) were identified based on three expression profile data sets. A multivariate Cox proportional regression model was used to develop a four-gene prognostic signature. Mechanistically, the correlationship between MMP1 and tumor-infiltrating lymphocytes was further analyzed. Furthermore, annotations were investigated by Gene Ontology (GO) and The Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses. Finally, potential drugs for CESC treatment were predicted by Connectivity Map. We profiled four genes (EFNA1, ANLN, MMP1, and ZWINT) with significant prognostic values for CESC. Multiple public available data sets were used for mRNA expression and prognostic characterization. Subsequently, GO and KEGG pathway analyses showed DEGs were mainly enriched in cell cycle, immunity, and metabolic-related pathway. We then conducted an integrated analysis of MMP1, and the expression of MMP1 showed significantly inverse association with the amount of CD8+ T cells, CD4+ T cells, and macrophages infiltration. Our findings suggest the four-gene signature may be associated with prognosis. We further revealed that MMP1 may be a novel biomarker for immunotherapy, and prognostic judgment of patients with cervical cancer.

Introduction

Cervical cancer remains one of the most common gynecologic malignancies, with >500,000 new cases and 300,000 deaths worldwide each year (Bray et al., 2018). The squamous cell carcinoma subtype is a main type, which accounts for 80–90% of cervical cancer cases (Printz, 2019). Ninety-five percent of the cases are caused by persistent infections with human papillomavirus (HPV), and HPV16 and HPV18 are the dominant subtypes (Schiffman et al., 2011). Although there were improvements in the diagnosis and therapy in the past decades, a significant number of patients were still diagnosed in late stage and could not get satisfactory clinical outcomes. Therefore, it is of great importance to further explore more effective prognostic biomarkers and therapeutic targets.

As HPV-driven cancer, cervical cancer appears at least partly mediated by the immune system (Komdeur et al., 2017). Checkpoint immunotherapy has shown significant effectiveness in lung, bladder, renal, and head and neck cancers (Otter et al., 2019). Cervical cancer will hopefully be within the expansion of the use of immunotherapeutic agents. However, the roles and mechanisms of how the interactions between tumor and the immune system impacting patients outcome remain unclear. Tumor-infiltrating lymphocytes (TILs) play an essential role in tumor-associated immune responses (Stanton et al., 2016). According to the specific antigens in the cell membrane, TILs could be divided into several subsets, and each of TIL subsets has different roles in the progression of the tumor. Many researchers have reported that TILs improved prognosis in various types of tumor, including colon, ovarian, lung, and breast cancer (Bremnes et al., 2016; Morris et al., 2018; Zhang et al., 2018; Truntzer et al., 2019). Of note, some studies have indicated a distinct correlation between different TIL intensities and different prognosis values (Wang et al., 2014; Vacchelli et al., 2015; Noble et al., 2016). The prognostic role of TILs is controversial currently; therefore, TILs have gained increasing attention in the treatment and prognosis prediction of tumors. Considering the considerable heterogeneity of TILs, the clinical impact of other immune cells in cervical squamous cell carcinoma (CESC) remains poorly understood.

This study aimed to identify useful signature genes to indicate poor clinical outcomes and investigate the correlations of tumor immunity with prognosis of CESC, which may assist us in evaluating therapeutic decisions and the prognosis for CESC patients. A signature of four genes was discovered and validated across the GEO cohorts, The Cancer Genome Atlas (TCGA) cohort, and public databases. The signature genes were not only predictive of clinical outcomes of CESC but also significantly associated with tumor immune. Furthermore, the most significant gene MMP1 was selected for further analysis. We investigated the correlation of MMP1 with tumor-infiltrating immune cells through Tumor Immune Estimation Resource (TIMER) and TISIDB. The results in this report revealed the critical role of MMP1 in CESC as well as provide an underlying mechanism and the potential relationship between MMP1 and tumor-immune interactions.

Materials and Methods

Gene Expression Omnibus data set and data processing

Gene expression raw microarray cell intensity (CEL) profiles of cervical cancer were evaluated in three independent data sets from the Gene Expression Omnibus database (accession nos. GSE63514, GSE6791, and GSE7410) (Pyeon et al., 2007; Biewenga et al., 2008; den Boon et al., 2015), which included 28 tumor tissue samples and 24 cancerous cervix samples; 20 tumor tissue samples and 8 cancerous cervix samples; 35 tumor tissue samples and 5 cancerous cervix samples, respectively. When more than one probe matched the same gene ID, the mean expression value of the gene was used for our study.

TCGA data set

The CESC data set containing 255 tumor samples and two normal samples, which included raw counts of RNAseq expression data and corresponding clinical information were obtained from TCGA data set. The TCGA data set was divided into three parts, two parts were selected randomly as a training cohort, and the complete data set was chosen as a validation cohort.

Identification of common differential expression genes

The GSE63514, GSE7410, and GSE6791 expression profiles were normalized and analyzed using the limma package (Smyth et al., 2005). The differential expression genes (DEGs) overlapping in the three cohorts were named as common DEGs. Furthermore, the featured genes were chosen based on the following criteria: (1) paired t-test analysis with p-value <0.05; (2) gene median expression in cervical cancer >0, and in adjacent normal tissue >0; (3) a two-tailed t-test was used to identify DEGs with the limma R package, and the |log2foldchange (log2FC)| >1 and adjusted p < 0.05 was set as the threshold to define the differences genes. All the data processing and normalization were performed using the R scripts.

Functional analysis of common DEGs

To identify their associated pathways and function annotations, Gene Ontology (GO) and The Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were conducted using DAVID Bioinformatics Resources (Huang da et al., 2009). GO analyses classified the common DEGs into three categories, such as biological process (BP), cellular component (CC), and molecular function (MF). KEGG analyses were conducted to determine the pathways in which the DEGs were significantly enriched in p < 0.05 and enrichment score (ES) >1.0 were used as the cutoff significant criteria.

Identification and selection of prognosis-related genes

Cox's proportional hazards model is commonly employed in survival analysis. To assess the association between the expression of common DEGs and patients' overall survival, univariable Cox regression analysis was applied by the survival package and 10 significant related genes were selected out (p < 0.05). Therefore, as candidate variables, the 10 prognosis-related genes were further analyzed through stepwise multivariate Cox proportional hazards regression analysis. By fitting all candidate predictive models into Akaike information criterions (AICs), an optimal predictive model was selected as with the lowest AIC value, and the best performance efficacy predictive model was acquired. Subsequently, four genes (EFNA1, ANLN, MMP1, and ZWINT) were selected, and a four-gene-based prognostic model was established to evaluate the survival risk of each patient as follows:

Consequently, each of the four prognosis-related genes weighted by their estimated regression coefficients in the multivariable Cox regression analysis developed the risk score formula in both the training set and the validation set.

Risk score formula establishment and validation

To validate the four-gene risk signature in the internal validation data sets, we calculated the risk score for each patient in the complete TCGA cohort. The patients were then divided into a high-risk and a low-risk group based on the corresponding median risk score. The prediction accuracy of this risk model was determined by a time-dependent receiver operating characteristic (ROC) analysis. To compare the survival time difference between the low-risk and high-risk group, Kaplan–Meier curve was plotted by the survival package using the log-rank test.

Signature genes selected and validation

Furthermore, prognosis-related genes in CESC was validated using public database PROGgeneV2 (Goswami and Nakshatri, 2014), Gene Expression Profiling Interactive Analysis (GEPIA) database (Tang et al., 2017), and Oncomine database (Rhodes et al., 2004). PROGgeneV2 contains data from 134 cohorts from 21 cancer types, was used to test the survival differences in combined signature genes expression. The data deposited in GEPIA are retrieved from the UCSC Xena server, including 9736 tumor samples and 8587 normal samples, whereas ONCOMINE is a cancer microarray database that contains 65 gene expression data sets comprising nearly 48 million gene expression measurements from >4700 microarray experiments. These two databases were used to verify the prognostic signature of the four genes.

Gene set enrichment analysis identifies a high-risk group-related signaling pathway

Samples of the complete cohort were divided into a high-risk group and low-risk group depending on the expression of the four-signature genes. To identify the potential mechanism of four-gene prognosis signature, gene set enrichment analysis (GSEA) was performed by GSEA software version 3.0 (Subramanian et al., 2005). Annotated gene sets c2.cp.kegg.v5.2.symbols.gmt were chosen as the reference gene sets. false discovery rate (FDR) <0.25, and |ES| >0.8 were used to determine statistical significance, and set as the cutoff criteria.

Mining the immune-related mechanism of MMP1

The TISIDB database integrates 988 reported immune-related anti-tumor genes, high-throughput screening techniques, molecular profiling, and paracancerous multi-omics data, as well as various resources for immunological data retrieved from seven public databases (Ru et al., 2019). In this study, the TISIDB was employed to analyze the relationship between the clinical outcome and the expression levels of MMP1. Then, the correlations between MMP1 expression and TILs were investigated further. TIMER is a web resource for systematical evaluations of the clinical impact of different immune cells in diverse cancer types. The TIMER was used to explore the correlations between expression of MMP1 and each TILs (B cells, CD4+ T cells, CD8+ T cells, neutrophils, macrophages, and dendritic cell types) (Li et al., 2017). Purity-corrected partial Spearman method was selected. REACTOME is function analysis module in TISIDB, resource for understanding functions of genes from the molecular level, was used to exhibit the enrichment analysis pathways that involved of MMP1. p-Value <0.05 was set as the threshold value.

Related small molecule drugs screening

The common DEGs between CESC samples and adjacent normal samples were used to query the CMap database, which was used to predict potential small molecule drugs (Lamb et al., 2006). A positive connectivity score suggests that a drug can induce the input signature, whereas a negative connectivity score suggests that a drug can reverse the input signature in human cell lines, and the drugs with negative connectivity scores were identified as candidate molecules with potential therapeutic value. After sorting all instances, the connectivity score of various instances was filtered by the number of instances (N > 10) and p-value < 0.05, and seven small molecules were selected. The ES representing similarity were calculated (ranging from −1 to 1). As a result, a negative connectivity score drug was selected for it means that the drug can reverse the input signature in human cell lines and identified as candidate molecules with potential therapeutic value. Besides, tomographs of these candidate small molecular medicine were observed in the PubChem database.

Results

Identification of common DEGs and functional annotation

The overall flowchart of this study is summarized in Figure 1. After the analyses of GSE63514, GSE6791, and GSE7410 data sets, a total of 120 CESC samples (83 tumor tissues and 37 normal tissues) were explored. There were 3201 DEGs screened from the GSE63514 data set, including 2340 upregulated and 861 downregulated genes; 2811 DEGs screened from the GSE6791 data set, including 312 upregulated and 2499 downregulated genes; 843 DEGs selected from the GSE7410 data set, including 491 upregulated and 352 downregulated genes. The volcano plots of DEGs among each data set are shown in Figure 2A–C. Cluster analyses of DEGs showed two significantly different distribution patterns between the tumor tissue and normal samples, suggesting crucial roles of DEGs in tumorigenesis and progression of CESC (Fig. 2D–F). By Venn diagram analysis, 132 common DEGs in the intersection of the three data sets were identified and selected for further analysis, and the top 16 genes are shown in Figure 3B.

FIG. 1.

Flowchart for this study. AIC, Akaike information criterion; DEGs, differential expression genes; GEPIA, Gene Expression Profiling Interactive Analysis; GO, Gene Ontology; GSEA, gene set enrichment analysis; KEGG, Kyoto Encyclopedia of Genes and Genomes; KM, Kaplan–Meier; TCGA, The Cancer Genome Atlas; TIMER, Tumor Immune Estimation Resource. Color images are available online.

FIG. 2.

DEGs in three data sets. (A–C) The volcano plots visualize the DEGs in GSE63514, GSE7410, and GSE6791, respectively. The red nodes represent upregulated genes. The blue nodes represent downregulated genes. (D–F) Heatmap of the top 100 DEGs according to the value of |logFC| > 1 and p < 0.05. The color in heat maps from green to red shows the progression from low expression to high expression. logFC, log fold change. Color images are available online.

FIG. 3.

(A) GO analysis of common DEGs. (B) Genes in three data sets and the top 16 genes are shown. BP, biological process; CC, cellular component; MF, molecular function. Color images are available online.

Functional annotation and KEGG pathway enrichment of common DEGs

GO analysis and KEGG pathway enrichment analysis were performed for the common DEGs to explore the potential BPs and function pathways associated with CESC. A total of 95 GO terms were significantly enriched (FDR <25%) (Supplementary Table S1). The top 10 GO terms in each functional category, including MF, BPs, and CC, are signified in Figure 3A. The GO terms majorly clustered in cell division, cell proliferation, cell cycle, immunity, and metabolic-related pathway, which are essential for rapid tumor growth (Supplementary Table S2). These results suggest that changes in immune and metabolic are necessary for tumorigenesis and development.

Moreover, the results of KEGG pathway analyses demonstrated that most of these target genes were enriched in cancer-related pathways, including p53 signaling pathway, transcriptional misregulation in cancer, and chemokine signaling pathway, indicating a role of the common DEGs in affecting tumorigenesis and progression of CESC (Table 1). We found that one gene may participate in multiple pathway, therefore, an intersection visualization plot was generated through UpSet plot (Fig. 4).

FIG. 4.

(A) UpSet plots of KEGG enrichment pathways in CESC. The horizontal axis represents the gene number of pathway types. The vertical axis represents the number of genes for one or several pathway types (B) Venn diagram of KEGG enrichment pathways in CESC. CESC, cervical squamous cell carcinoma. Color images are available online.

Table 1.

The Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analysis of Common Differential Expression Genes

Category	Term	Count	%	p	Genes	List total	Pop hits	Pop total	Fold enrichment	Bonferroni	Benjamini	FDR
KEGG_PATHWAY	hsa04110:cell cycle	8	6.060606061	1.13E-04	CCNB1, CDC7, CCNE1, CDKN2A, DBF4, PCNA, CHEK1, MCM6	63	124	6879	7.044546851	0.011549328	0.011549328	0.125797836
KEGG_PATHWAY	hsa03030:DNA replication	5	3.787878788	2.84E-04	RFC5, PCNA, RNASEH2A, FEN1, MCM6	63	36	6879	15.16534392	0.028830298	0.014520572	0.316495542
KEGG_PATHWAY	hsa04115:p53 signaling pathway	5	3.787878788	0.00300068	CCNB1, CCNE1, CDKN2A, RRM2, CHEK1	63	67	6879	8.148542999	0.266211642	0.098033781	3.298495596
KEGG_PATHWAY	hsa05161:hepatitis B	6	4.545454545	0.009586817	DDX58, CCNE1, MMP9, PCNA, TLR2, STAT1	63	145	6879	4.518226601	0.629242055	0.219680219	10.19374869
KEGG_PATHWAY	hsa04620:Toll-like receptor signaling pathway	5	3.787878788	0.014994424	TLR2, CXCL9, CXCL11, STAT1, CXCL10	63	106	6879	5.15049416	0.789047363	0.267450474	15.51714822
KEGG_PATHWAY	hsa04623:cytosolic DNA-sensing pathway	4	3.03030303	0.019675268	DDX58, POLR2H, AIM2, CXCL10	63	64	6879	6.824404762	0.870846605	0.289030499	19.89137572
KEGG_PATHWAY	hsa05200:pathways in cancer	9	6.818181818	0.023840409	CCNE1, CKS1B, CDKN2A, MMP9, LAMC2, BIRC3, STAT1, MMP1, FZD6	63	393	6879	2.500545256	0.916699989	0.298857267	23.6092257
KEGG_PATHWAY	hsa05162:measles	5	3.787878788	0.031447925	DDX58, CCNE1, TLR2, MX1, STAT1	63	133	6879	4.10490512	0.962789414	0.337275061	29.99691553
KEGG_PATHWAY	hsa03410:base excision repair	3	2.272727273	0.035284211	PCNA, PARP1, FEN1	63	33	6879	9.926406926	0.97527497	0.337082251	33.03004168
KEGG_PATHWAY	hsa05222:small cell lung cancer	4	3.03030303	0.040865219	CCNE1, CKS1B, LAMC2, BIRC3	63	85	6879	5.13837535	0.986397995	0.349330733	37.22931152
KEGG_PATHWAY	hsa05206:microRNAs in cancer	7	5.303030303	0.042965198	KIF23, CCNE1, CDKN2A, MMP9, CDCA5, BRCA1, PLAU	63	286	6879	2.672494172	0.989146838	0.337151608	38.74627012
KEGG_PATHWAY	hsa05219:bladder cancer	3	2.272727273	0.052353222	CDKN2A, MMP9, MMP1	63	41	6879	7.989547038	0.996068227	0.369697484	45.1282018
KEGG_PATHWAY	hsa05202:transcriptional misregulation in cancer	5	3.787878788	0.063111845	MMP9, GZMB, MMP3, PLAU, MLF1	63	167	6879	3.269175934	0.998787063	0.403404667	51.6936538
KEGG_PATHWAY	hsa04668:TNF signaling pathway	4	3.03030303	0.071582149	MMP9, BIRC3, MMP3, CXCL10	63	107	6879	4.08188696	0.999524031	0.420992497	56.35024083
KEGG_PATHWAY	hsa05145:toxoplasmosis	4	3.03030303	0.076397049	TLR2, LAMC2, BIRC3, STAT1	63	110	6879	3.970562771	0.999721396	0.42057245	58.81129009
KEGG_PATHWAY	hsa04062:chemokine signaling pathway	5	3.787878788	0.086118463	CXCL13, CXCL9, CXCL11, STAT1, CXCL10	63	186	6879	2.935227855	0.999906317	0.439947484	63.39938473

FDR, false discovery rate; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Signature genes showed the greatest prognostic value in patients with cervical tumors

Given the significant association of common DEGs with tumorigenesis and development, we further investigate the potential roles in clinical outcomes. Univariate Cox proportional hazards regression revealed that EFNA1, ANLN, MMP1, TFRC, ZWINT, APOBEC3B, TRIM59, IFT80, FZD6, and PLOD2 were significant prognostic factors (Table 2). Multivariate Cox proportional hazards regression analysis of the 10 genes revealed a significant prognostic model comprising four different genes (EFNA1, ZWINT, ANLN, and MMP1). Among the four genes EFNA1, ANLN, and MMP1 were risk factors (hazard ratio (HR)>1), whereas ZWINT was identified as protective factor (HR <1). More information about these four genes is given in Table 3. The two ROCs exhibited high predictive abilities in both training cohort (n = 160) and validation cohort (n = 255) (Fig. 5A). The Kaplan–Meier analysis showed that patients in the high-risk group had a shorter overall survival time than patients in the low-risk group (p < 0.01; Fig. 5B, C).

FIG. 5.

Validation the prognostic value of the four-gene signature. (A) Five-year survival time-dependent ROC curve in training cohort and validation cohort, separately. The area under the ROC curve (AUC) was 0.878 and 0.782, which indicated modest diagnostic value. (B, C) Kaplan–Meier survival curve of overall survival between the high-risk group and low-risk group. (D) Validation of the survival difference between the high-risk group and low-risk group in PROGgeneV2 database. AUC, area under the curve; ROC, receiver operating characteristic. Color images are available online.

Table 2.

The Detailed Information of the 10 Genes Significantly Associated with Overall Survival in Univariable Analyses

Ensemble ID	Gene symbol	HR	p	Genomic locations (GRCh38/hg38)
ENSG00000169242	EFNA1	1.965778251	0.002510695	chr1:155,127,873-155,134,910
ENSG00000011426	ANLN	1.89300833	0.007164518	chr7:36,389,803-36,453,791
ENSG00000196611	MMP1	1.198677	0.011758633	chr11:102,789,910-102,798,235
ENSG00000072274	TFRC	1.415421416	0.025752334	chr3:196,027,183-196,082,161
ENSG00000122952	ZWINT	0.687348177	0.026615596	chr10:56,357,227-56,361,302
ENSG00000179750	APOBEC3B	0.801839856	0.03394592	chr22:38,982,347-38,992,804
ENSG00000213186	TRIM59	0.62506258	0.035688667	chr3:160,432,445-160,485,773
ENSG00000068885	IFT80	0.504393715	0.041878078	chr3:160,256,986-160,399,880
ENSG00000164930	FZD6	1.627645866	0.049183351	chr8:103,298,433-103,332,866
ENSG00000152952	PLOD2	1.429278727	0.049192392	chr3:146,069,436-146,163,653

HR, hazard ratio.

Table 3.

The Four Prognosis-Associated Genes to Establish the Risk Score System

Gene symbol	Description	Coefficient	HR	95% CI	p
EFNA1	Ephrin A1	0.4988	1.6467	0.8789–3.0853	0.11949
ZWINT	ZW10 Interacting Kinetochore Protein	−0.4720	0.6238	0.4088–0.9518	0.02860
ANLN	Anillin Actin Binding Protein	0.7490	2.1149	1.3060–3.4245	0.00232
MMP1	Matrix Metallopeptidase 1	0.2603	1.2974	1.1069–1.5206	0.00131

Bold font indicates p-value is less than 0.05.

CI, confidence interval.

Validation and analysis of signature genes

We further compared the survival difference between the high-risk group and the low-risk group in external database PROGgeneV2. The consistent results was obtained using from external verification that the high-risk group had had a shorter overall survival time (Fig. 5D, HR: 1.68 [1.27–2.22], p = 0.0002777). In addition, we analyzed the differential mRNA expression levels of EFNA1, ZWINT, ANLN, and MMP1 between cervical cancer and normal cervical tissues using the GEPIA database, which are generally consistent with those obtained from the Oncomine database (Fig. 6). Both databases showed that EFNA1, ZWINT, ANLN, and MMP1 mRNA expression levels were upregulated in cervical cancer compared with normal tissues, and indicated that EFNA1, ZWINT, ANLN, and MMP1 play important roles during the tumorigenesis and development of CESC.

FIG. 6.

Validation of the gene expression levels of EFNA1, ZWINT, ANLN, and MMP1 between CESC tissues and normal cervical in GEPIA database. (A–D) EFNA1, ZWINT, ANLN, and MMP1 are significantly upregulated in CESC tissues compared with normal tissues (p < 0.01). The red asterisk represents p < 0.01. (E) Gene expression level of the four genes (based on Oncomine database). Color images are available online.

Gene set enrichment analysis

Owing to such different prognostic outcomes between the high-risk and low-risk groups, we sought to investigate possible differences and molecular mechanisms by GSEA, which verified an enhanced antiapoptosis and epithelial-mesenchymal transition (EMT) signature in the high-risk group (Fig. 7). Therefore, we propose the existence of an intense regulatory role for the development and progression in high-risk cervical cancer patients.

FIG. 7.

GSEA results for high- vs. low-risk group (A) apoptosis, (B) epithelial-mesenchymal-transition. Color images are available online.

Regulation of immune molecules by MMP1

The KEGG pathway enrichment analysis showed that MMP1 is a hub gene of two cancer-related pathways. Besides, TISIDB-Reactome pathway analysis support suggested that MMP1 was involved in the cytokine signaling in the immune system (Table 4). Hence, we conducted an integrated analysis to predict the potential biological roles of MMP1 in tumor-immune of CESC. The results show that the expression of MMP1 was correlated with the abundance of TILs, whereas copy number does not appear to correlate with PD-L1 expression in cervical cancer (Fig. 8). More important, we found that patients with low expression of MMP1 exhibited significant treatment response of checkpoint immunotherapy (Table 5). To identify signatures and mechanisms of the relationship between MMP1 and response to checkpoint blockade (e.g., anti-PDL1 and anti-PD1), correlation analysis was performed, indicating that the expression of MMP1 is positively correlated with infiltration of macrophages, while negatively correlated with infiltration of CD8+ memory T and CD4+ memory T, which may be the reason why high MMP1 expression leads to better immune checkpoint treatment response (Figs. 8 and 9).

FIG. 8.

(A) Prognostic values of MMP1 in CESC. (B) Associations between MMP1 expression and overall survival across human cancers. (C) Correlation between MMP1 expression and immunity in CESC. (D) The “SCNA” module provides the comparison of the abundance of TIICs among tumors with different somatic copy number aberrations for a given gene. SCNAs are defined by GISTIC 2.0 (14), including deep deletion, shallow deletion, diploid/normal, low-level gain, and high amplification. Box plots are presented to show the distributions of each TIIC subset for each copy number status in selected cancer types using the same statistical tests as in the “mutation” module. Color images are available online.

FIG. 9.

(A) Pearman's correlation of ANXA8 with lymphocytes and immunomodulators (TISIDB). (B) Relations between the abundance of macrophages and MMP1 expression. Color images are available online.

Table 4.

Reactome-MMP1

Stable identifier	Pathway
R-HSA-1592389	Activation of Matrix Metalloproteinases
R-HSA-210991	Basigin interactions
R-HSA-202733	Cell surface interactions at the vascular wall
R-HSA-1442490	Collagen degradation
R-HSA-1280215	Cytokine signaling in immune system
R-HSA-1474228	Degradation of the extracellular matrix
R-HSA-1474244	Extracellular matrix organization
R-HSA-109582	Hemostasis
R-HSA-168256	Immune system
R-HSA-6785807	Interleukin-4 and -13 signaling
R-HSA-392499	Metabolism of proteins
R-HSA-381426	Regulation of insulin-like growth factor transport and uptake by insulin-like growth factor binding proteins

Table 5.

Points in the Scatter Plot Represent the Expression Difference of MMP1 in Various Data Sets (TISIDB)

No.	PMID	Cancer type	Group	Drug	Log2 (fold change)	p
1	26997480	Melanoma	All	Anti-PD-1	−2.485	0.0458
2	26997480	Melanoma	MAPKi	Anti-PD-1	−4.534	0.0355
3	28552987	Urothelial cancer	All	Anti-PD-L1	−2.367	0.0168
4	28552987	Urothelial cancer	Smoking	Anti-PD-L1	−3.955	0.0132

See scatter plot, Fig. 8C.

CMap analysis

By CMap analyzed, we identified seven drugs with positive CMap ES of >0 and PCMap-value of <0.001 the most related small molecule drug trichostatin A with highly significant correlations and the maximum number of instances were identified. Trichostatin A, a histone deacetylase (HDAC) inhibitor, showed a significantly negative correlation and the tomography of trichostatin A was found in the PubChem database and shown in Figure 10. The screened small molecules are listed in Table 6.

FIG. 10.

The tomographs of the candidate small molecule drug trichostatin A for CESC. Color images are available online.

Table 6.

The Top 10 CESC-Related Small Molecules with Highly Significant Correlations in Results of CMap Analysis

Rank	CMap name	Mean	n	Enrichment	p	Specificity	% Non-null
1	Phenoxybenzamine	−0.745	4	−0.975	0	0	100
2	Resveratrol	−0.755	9	−0.9	0	0	100
3	Vorinostat	−0.648	12	−0.82	0	0.0442	100
4	Trichostatin A	−0.388	182	−0.451	0	0.2673	64
5	LY-294002	−0.253	61	−0.333	0	0.2577	54
6	Nortriptyline	−0.664	4	−0.932	0.00002	0	100
7	Clotrimazole	−0.611	5	−0.87	0.00008	0	100
8	15-Delta prostaglandin J2	−0.417	15	−0.559	0.00008	0.0602	60
9	Trifluoperazine	−0.434	16	−0.502	0.00022	0.1154	62
10	Puromycin	−0.68	4	−0.892	0.00026	0.0373	100
11	Monobenzone	−0.692	4	−0.892	0.00026	0	100
12	Chlorcyclizine	−0.398	6	−0.724	0.00085	0	66
13	Fluphenazine	−0.368	18	−0.449	0.00099	0.0588	61

Drug in bold type is the most promising treatment small molecule drug identified in this study.

Discussion

In this study, we used comprehensive bioinformatics analyses to identify four prognosis-related signature genes (EFNA1, ZWINT, ANLN, and MMP1), the expression of which was significantly related to OS in CESC patients. In addition, the risk score of the four-gene signature demonstrated superior ability to divide CESC patients into high-risk and low-risk groups with significantly different OS in training data set and validation data set, as well as in Pan Cancer Prognostics Database (PROGgeneV2).

Furthermore, ROC curve analysis also suggested that these genes had an excellent value for predicting the OS for CESC patients. In the further GSEAs, we found that high-risk group showed an enhanced antiapoptosis and epithelial-mesenchymal pathway, implying that the four-gene signature may be involved in regulation of the prognosis of cervical cancer through antiapoptosis and the EMT pathways. Activation of EMT pathway promotes local tumor invasion and tumor cell metastasis in various cancers (De Craene and Berx, 2013; Jung et al., 2015). Qureshi et al. (2015) have reported that the EMT plays a vital role in cervical cancer progression and metastasis. In addition, previous studies have demonstrated that activation of EMT pathway promotes lymph node (LN) metastasis in cervical cancer (Noordhuis et al., 2011). Taken together, these results suggest that the high-risk populations may have a worse prognosis because EMT pathway activation promotes cervical cancer progression, which contributed to the poor prognosis of CESC.

The signature gene MMP1 was identified as biomarkers in lung adenocarcinoma, esophageal cancer, bladder cancer (Wieczorek et al., 2015; Liu et al., 2016; Ni and Sun 2019). ANLN was reported as a prognostic marker for head and neck squamous cell carcinoma (Shimizu et al., 2007). However, no previous study has considered the possible role of EFNA1, ZWINT in cancer. We first proposed two new prognostic markers for cervical cancer. Therefore, the signature may be a potential good predictor of prognosis for cervical cancer patients. This research provides more deep insights into possible molecular mechanisms of prognosis for CESC. However, all of these four genes need further experimental verification to demonstrate their roles in CESC.

Significantly, MMP1 was chose for subsequent in-depth investigations for its correlations with tumor immunity and possible mechanism. Finally, we selected a small molecule drug that has the potential values to treat CESC by CMap database.

The mRNA expression level of the four genes was significantly upregulated in cervical cancer compared with normal, both revealed in silico and Oncomine database analyses. Patients with high EFNA1, ZWINT, ANLN, and MMP1 expression had significantly shorter survival than patients with lower expression of these genes. Moreover, the signature was confirmed to harbor high predictive value of CESC both in the validation cohort and the online databases PROGgeneV2. EFNA1 has been found as key serum markers for the detection of hepatocellular carcinoma development and progression, but no reporter in literature in CESC to date (Tian et al., 2010). Previously, Suman and Mishra (2018) reported that ZWINT were independent prognostic factors. Mechanistically, ANLN may be responsible for progressing cervical neoplasia to cervical cancer at an early stage (Xia et al., 2018). MMP1 was considered having a pivotal role in the regulation of cervical tumor growth and LN metastasis through EMT (Tian et al., 2019). These results are consistent with ours, and suggest that the four signature genes had an excellent value for predicting the 5-year survival prognosis for CESC.

To further explore the potential mechanism of the four genes, we observe that high-risk phenotype was associated with apoptosis pathway and EMT by GSEA. Several researchers had reported that EMT is the most important mechanism that plays a crucial role in the metastasis and invasion of cervical cancer (Feng et al., 2019; Li et al., 2019; Zeng et al., 2019) These results represent candidate biomarkers for predicting prognosis of CESC patients.

Given the prognostic role of MMP1 expression in CESC has been reported in several studies, its role in cervical cancer the immune modulation remains unstudied (Zhang and Zhao, 2016; Roszik et al., 2018; Tian et al., 2018). The KEGG pathway enrichment analysis showed that MMP1 is a hub gene of two cancer-related pathways. Besides, TISIDB-Reactome pathway analysis support suggested that MMP1 was involved in the cytokine signaling in immune system. Hence, we conducted an integrated analysis to predict the potential biological roles of MMP1 in tumor-immune of CESC.

In further study, we found out that there was a significant difference in the expression level of MMP1 between responders and nonresponders of PD-1 in four clinical immunotherapy researches. Indeed, low expression of MMP1 with better immunotherapy response rate. In addition, correlation analysis reported that the expression of MMP1 was correlated with the abundance of TILs, whereas copy number does not appear to correlate with PD-L1 expression in cervical cancer. PD-L1 immunohistochemistry alone is not yet an adequate biomarker for routine clinical use in deciding which patients to offer anti-PD-1 or anti-PD-L1 therapy. It has recently confirmed that baseline tumor-infiltrating lymphocyte status could also serve as a predictive biomarker for checkpoint inhibitor immunotherapy (Otter et al., 2016). Otter et al. (2019) suggested that the best predictive marker of response to PD-1 blocking agents in pretreatment biopsies was CD8+ memory T cell density at the invasive margin. In fact, the CD8+ memory T cell subset, displayed an activated and proliferative state and shown to be crucial for effective in vivo PD-1 blockade (Heeren et al., 2019). Therefore, CD8+ T cell density was more predictive of the response to treatment of PD-1 and PD-L1 blocking agents. However, no report has demonstrated a correlation between MMP1 expression and CD8+ memory in CESC. In this study, it is suggested that low expressed MMP1 may be crucial for evoking the enrichment of the effective (active) CD8+ memory T. In this study, we have for the first time revealed this inverse correlation between MMP1 and CD8+ memory T cells.

In contrast, CD4+ memory T cells have been shown to predict prognosis better than TNM (tumor-node-metastasis) classification system in colorectal cancer (Mlecnik et al., 2011). Moreover, recent reports have suggested that also showed a significant correlation between circulating highly differentiated memory CD4 T cells and objective clinical responses. Patients with high percentages of these memory subsets within the CD4 cell population also had about 50% response rates (Zuazo et al., 2018, 2019). Interestingly, patients with high memory CD4 T cells that also had lower percentages of PD-L1+ CD11b+ cells did not show objective responses. These results suggest that combining key systemic cellular parameters could be used to stratify patients with a high probability of responding to PD-L1/PD-1 blockade. It is likely that PD-L1 blockade in these cells may enhance their antigen presentation capacities, and together with memory CD4 T cells trigger strong antitumor activities. Consistently, PD-L1 expression in systemic immune cell populations was identified as a potential predictive biomarker of responses to PD-L1/PD-1 blockade therapy in lung cancer (Bocanegra et al., 2019). These results suggest that a combined memory CD4 T cells and expression of PD-L1 may be beneficial biomarkers. In fact, our results complemented the immune environment regulation values of MMP1 in CESC with multidimensional analysis.

Macrophages are one of the most abundant immune cells in the tumor microenvironment of solid tumors, and their presence correlates with reduced survival in most cancers (Nielsen and Schmid, 2017). Tumor-associated macrophages (TAMs) display several protumorigenic functions that have essential roles in cancer development and progression such as the ability to provide cytokines and induce tumor angiogenesis (Grivennikov et al., 2010). TAMs are a source of tumor-promoting interleukin (IL)-6 in several murine tumor models. Promotes colon tumor cell proliferation and protection from apoptosis through activation of signal transducers and activator of 66 transcription-3 (STAT3) (Bollrath et al., 2009; Grivennikov et al., 2009; Shrivastava et al., 2013). The poor prognosis of CESC caused by high expression of MMP1 may be related to the massive infiltration of TAM. Furthermore, through such a STAT3 activation process, CESC patients with high MMP1 expression may lead to advanced disease and present a higher risk for poor outcome, which becomes the basis for developing prognosis prediction biomarkers.

Furthermore, GSVA (TISIDB-Reactome) analysis demonstrated that the BPs of the genes were mainly associated with immune suppression, such as infection by HPV, cytokine signaling in immune system, IL-4 and IL-13 signaling. Previous studies indicated that STAT3, one of STAT family members, has been reported to be activated by inflammatory cytokines such as IL-4 and IL-13 in 68 integral signaling pathways (Lesina et al., 2011). IL-4 and IL-13 were involved in the STAT3 pathway, indicating that MMP1 may control macrophages to activate the STAT3 pathway to promote tumor progression. Collectively, these results demonstrated that MMP1 is a high candidate as a biomarker for immunotherapy and prognostic judgments of cervical cancer.

Revealing the relationship between MMP1 and immune infiltrating cells and identifying MMP1 as an essential prognostic factor are the highlights of this study. We are currently expanding these analyses during immunotherapy treatments to evaluate the dynamic changes of these cervical cancer populations in responders and progressors. Considering both the immune and gene expression profiling data provided by our study, we feel MMP1 is not only a potential prognostic biomarker but also a high priority consideration for predicting the efficacy of PD-1/PD-L1 targeted therapy. Furthermore, to identify more potential drugs for CESC treatment based on the common DEG profiles, CMap analysis was conducted and 18 candidate drugs (p < 0.001) were obtained from the CMap data set.

Moreover, we identified several small molecules with potential therapeutic efficacy against CESC. Trichostatin A (connectivity score of −0.788) was particularly interesting and it is considered to be the most important small molecule in our screening for the most research data supporting. Trichostatin A is a HDAC inhibitor (HDACi). HDACi, a novel class of anticancer drug, has been recently reported to suppress host immunity. Many studies demonstrated that trichostatin A may be a useful therapeutic tool for the reversion of EMT in various types of cancers, including colon cancer (Huang et al., 2019) and breast cancer (Vigushin et al., 2001). Reasonably it is presumed that the potential anticancer effect of trichostatin A could be further characterized and validated in further study.

Limitation

We would like to acknowledge some of the limitations of this study. First, this is a retrospective study; all the data of this study were obtained from publicly available databases. Second, further studies, including in vivo and in vitro experiments, are required to confirm the results. Third, the exact mechanisms of action of the four genes in CESC remain to be fully elucidated, and the molecular mechanism of MMP1 regulating the TILs infiltration will be another research point for us. Despite this limitation, these findings present an opportunity to rationally approach future combination immunotherapy trials in the treatment of recurrent cervical cancer. These discoveries suggest that a comprehensive investigation of tumor and immune cell interplay would assist both the elucidation of cervical cancer pathogenesis and the development of novel immunotherapy strategies.

Conclusion

In this study we identified a four-gene signature that may prove to be associated with progression and prognosis, which can provide the favorable prognostic value based on these four gene expression levels, and a candidate small molecule drug was identified, which provides direction for CESC targeted therapy. With further exploration, our data revealed that MMP1 was significantly highly expressed in CESC tissues and that high MMP1 expression was significantly correlated with poor prognosis TILs. Therefore, MMP1 may prove to be a novel biomarker for immunotherapy and prognostic judgment of patients with cervical cancer.

Footnotes

Acknowledgments

We would like to acknowledge the KEGG database developed by Kanehisa Laboratories. We also would like to acknowledge the Oncomine, GEPIA, TISIDB, TIMER, and TCGA databases for free use.

Authors' Contributions

M.Y. designed the study, as well as designed the figures and tables. S.Z. contributed to the statistical analysis, as well as wrote and corrected the article. All authors read and approved the final article.

Ethics Approval and Consent to Participate

Not applicable.

Patient Consent for Publication

Not applicable.

Disclosure Statement

The authors declare that they have no competing interests.

Funding Information

This study was supported by Natural Science Foundation of Liaoning Province (Grant No. 20180550781).

Supplementary Material

Supplementary Table S1

Supplementary Table S2

References

Biewenga

, Buist

M.R.

, Moerland

P.D.

, Ver Loren van Themaat

, van Kampen

A.H.

, ten Kate

F.J.

, et al. (2008). Gene expression in early stage cervical cancer. Gynecol Oncol, 108, 520–526.

Bocanegra

, Fernandez-Hinojal

, Zuazo-Ibarra

, Arasanz

, Garcia-Granda

M.J.

, Hernandez

, et al. (2019). PD-L1 expression in systemic immune cell populations as a potential predictive biomarker of responses to PD-L1/PD-1 blockade therapy in lung cancer. Int J Mol Sci, 2, 20.

Bollrath

, Phesse

T.J.

, von Burstin

V.A.

, Putoczki

, Bennecke

, Bateman

, et al. (2009). gp130-mediated Stat3 activation in enterocytes regulates cell survival and cell-cycle progression during colitis-associated tumorigenesis. Cancer Cell, 15, 91–102.

Bray

, Ferlay

, Soerjomataram

, Siegel

R.L.

, Torre

L.A.

, and Jemal

(2018). Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin, 68, 394–424.

Bremnes

R.M.

, Busund

L.T.

, Kilvær

T.L.

, Andersen

, Richardsen

, and Paulsen

E.E.

(2016). The role of tumor-infiltrating lymphocytes in development, progression, and prognosis of non-small cell lung cancer. J Thorac Oncol, 11, 789–800.

Cui

X.D.

, Lee

M.J.

, Yu

G.R.

, Kim

I.H.

, Yu

H.C.

, Song

E.Y.

, and Kim

D.G.

(2010). EFNA1 ligand and its receptor EphA2: potential biomarkers for hepatocellular carcinoma. Int J Cancer, 126, 940–949.

De Craene

, and Berx

(2013). Regulatory networks defining EMT during cancer initiation and progression. Nat Rev Cancer, 13, 97–110.

Den Boon

J.A.

, Pyeon

, Wang

S.S.

, Horswill

, Schiffman

, Sherman

, et al. (2015). Molecular transitions from papillomavirus infection to cervical precancer and cancer: role of stromal estrogen receptor signaling. Proc Natl Acad Sci U S A, 112, 3255–3264.

Feng

, Liu

, Bai

, Pan

, Jia

, Zhang

, et al. (2019). LncRNA-CTS promotes metastasis and epithelial-to-mesenchymal transition through regulating miR-505/ZEB2 axis in cervical cancer. Cancer Lett, 465, 105–117.

10.

Gibney

G.T.

, Weiner

L.M.

, and Atkins

M.B.

(2016). Predictive biomarkers for checkpoint inhibitor-based immunotherapy. Lancet Oncol, 17, e542–e551.

11.

Goswami

C.P.

, and Nakshatri

(2014). PROGgeneV2: enhancements on the existing database. BMC Cancer, 17, 14–970.

12.

Grivennikov

, Karin

, Terzic

, Mucida

, Yu

G.Y.

, Vallabhapurapu

, et al. (2009). IL-6 and Stat3 are required for survival of intestinal epithelial cells and development of colitis-associated cancer. Cancer Cell, 15, 103–113.

13.

Grivennikov

S.I.

, Greten

F.R.

, and Karin

(2010). Immunity, inflammation, and cancer. Cell, 140, 883–899.

14.

Heeren

A.M.

, Rotman

, Stam

A.G.M.

, Pocorni

, Gassama

A.A.

, Samuels

, et al. (2019). Efficacy of PD-1 blockade in cervical cancer is related to a CD8+FoxP3+CD25+ T-cell subset with operational effector functions despite high immune checkpoint levels. J Immunother Cancer, 7, 43.

15.

Huang

, Wu

X.F.

, and Wang

X.L.

(2019). Trichostatin A inhibits phenotypic transition and induces apoptosis of the TAF-treated normal colonic epithelial cells through regulation of TGF-β pathway. Int J Biochem Cell Biol, 114, 105565.

16.

Huang da

, Sherman

B.T.

, and Lempicki

R.A.

(2009). Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc, 4, 44–57.

17.

Jung

H.Y.

, Fattet

, and Yang

(2015). Molecular pathways: linking tumor microenvironment to epithelial-mesenchymal transition in metastasis. Clin Cancer Res, 21, 962–968.

18.

Komdeur

F.L.

, Prins

T.M.

, van de Wall

, Plat

, Wisman

G.B.A.

, Hollema

, et al. (2017). CD103+ tumor-infiltrating lymphocytes are tumor-reactive intraepithelial CD8+ T cells associated with prognostic benefit and therapy response in cervical cancer. Oncoimmunology, 24, 6.

19.

Lamb

, Crawford

E.D.

, Peck

, Modell

J.W.

, Blat

I.C.

, Wrobel

M.J.

, et al. (2006). The Connectivity Map: using gene-expression signatures to connect small molecules, genes, and disease. Science, 313, 1929–1935.

20.

Lesina

, Kurkowski

M.U.

, Ludes

, Rose-John

, Treiber

, Klöppel

, et al. (2011). Stat3/Socs3 activation by IL-6 transsignaling promotes progression of pancreatic intraepithelial neoplasia and development of pancreatic cancer. Cancer Cell, 19, 456–469.

21.

, Zhang

, Yan

, Qiu

, Chen

, Sun

, et al. (2019). Nucleolar and spindle associated protein 1 promotes metastasis of cervical carcinoma cells by activating Wnt/β-catenin signaling. J Exp Clin Cancer Res, 38, 33.

22.

, Fan

, Wang

, Traugh

, Chen

, Liu

J.S.

, et al. (2017). TIMER: a web server for comprehensive analysis of tumor-infiltrating immune cells. Cancer Res, 77, e108–e110.

23.

Liu

, Hu

, Zhang

M.F.

, Luo

K.J.

, Xie

X.Y.

, Wen

, et al. (2016). MMP1 promotes tumor growth and metastasis in esophageal squamous cell carcinoma. Cancer Lett, 377, 97–104.

24.

Mlecnik

, Tosolini

, Kirilovsky

, Berger

, Bindea

, Meatchi

, et al. (2011). Histopathologic-based prognostic factors of colorectal cancers are associated with the state of the local immune reaction. J Clin Oncol, 29, 610–618.

25.

Morris

, Platell

, and Iacopetta

(2008). Tumor-infiltrating lymphocytes and perforation in colon cancer predict positive response to 5-fluorouracil chemotherapy. Clin Cancer Res, 14, 1413–1417.

26.

K.W.

, and Sun

G.Z.

(2019). The identification of key biomarkers in patients with lung adenocarcinoma based on bioinformatics. Math Biosci Eng, 16, 7671–7687.

27.

Nielsen

S.R.

, and Schmid

M.C.

(2017). Macrophages as key drivers of cancer progression and metastasis. Mediators Inflamm, 2017, 9624760.

28.

Noble

, Mellows

, McCormick Matthews

L.H.

, Bateman

A.C.

, Harris

, Underwood

T.J.

, et al. (2016). Tumour infiltrating lymphocytes correlate with improved survival in patients with oesophageal adenocarcinoma. Cancer Immunol Immunother, 65, 651–662.

29.

Noordhuis

M.G.

, Fehrmann

R.S.

, Wisman

G.B.

, Nijhuis

E.R.

, van Zanden

J.J.

, Moerland

P.D.

, et al. (2011). Involvement of the TGF-beta and beta-catenin pathways in pelvic lymph node metastasis in early-stage cervical cancer. Clin Cancer Res, 17, 1317–1330.

30.

Otter

S.J.

, Chatterjee

, Stewart

A.J.

, and Michael

(2019). The role of biomarkers for the prediction of response to checkpoint immunotherapy and the rationale for the use of checkpoint immunotherapy in cervical cancer. Clin Oncol (R Coll Radiol), 31, 834–843.

31.

Printz

(2019). Cervical cancer subtype less detectable by Papanicolaou testing. Cancer, 125, 3287.

32.

Pyeon

, Newton

M.A.

, Lambert

P.F.

, den Boon

J.A.

, Sengupta

, Marsit

C.J.

, et al. (2007). Fundamental differences in cell cycle deregulation in human papillomavirus-positive and human papillomavirus-negative head/neck and cervical cancers. Cancer Res, 15, 4605–4619.

33.

Qureshi

, Arora

, and Rizvi

M.A.

(2015). EMT in cervical cancer: its role in tumour progression and response to therapy. Cancer Lett, 356, 321–331

34.

Rhodes

D.R.

, Yu

, Shanker

, Deshpande

, Varambally

, Ghosh

, et al. (2004). ONCOMINE: a cancer microarray database and integrated data-mining platform. Neoplasia, 6, 1–6.

35.

Roszik

, Ring

K.L.

, Wani

K.M.

, Lazar

A.J.

, Yemelyanova

A.V.

, Soliman

P.T.

, et al. (2018). Gene expression analysis identifies novel targets for cervical cancer therapy. Front Immunol, 9, 2102.

36.

, Wong

C.N.

, Tong

, Yi Zhong

, Wa Zhong

S.S.

, Wu

W.C.

, et al. (2019). TISIDB: an integrated repository portal for tumor-immune system interactions. Bioinformatics, 35, 4200–4202.

37.

Schiffman

, Wentzensen

, Wacholder

, Kinney

, Gage

J.C.

, and Castle

P.E.

(2011). Human papillomavirus testing in the prevention of cervical cancer. J Natl Cancer Inst, 103, 368–383.

38.

Shimizu

, Seki

, Sugimoto

, Horiguchi

, Tanzawa

, Hanazawa

, et al. (2007). Identification of molecular targets in head and neck squamous cell carcinomas based on genome-wide gene expression profiling. Oncol Rep, 18, 1489–1497.

39.

Shrivastava

, Llovera

, Recasens

, Chertoff

, Giménez-Llort

, Gonzalez

, et al. (2013). Temporal expression of cytokines and signal transducer and activator of transcription factor 3 activation after neonatal hypoxia/ischemia in mice. Dev Neurosci, 35, 212–225.

40.

Smyth

G.K.

, Michaud

, and Scott

H.S.

(2005). Use of within-array replicate spots for assessing differential expression in microarray experiments. Bioinformatics, 21, 2067–2075.

41.

Stanton

S.E.

, and Disis

M.L.

(2016). Clinical significance of tumor-infiltrating lymphocytes in breast cancer. J Immunother Cancer, 18, 4–59.

42.

Subramanian

, Tamayo

, Mootha

V.K.

, Mukherjee

, Ebert

B.L.

, Gillette

M.A.

, et al. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A, 102, 15545–15550.

43.

Suman

, and Mishra

(2018). An interaction network driven approach for identifying biomarkers for progressing cervical intraepithelial neoplasia. Sci Rep, 8, 12927.

44.

Tang

, Li

, Kang

, Gao

, Li

, and Zhang

(2017). GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res, 45, 98–102.

45.

Tian

, Li

, Gao

, Li

, Yang

, and Wang

(2018). Identification and validation of the role of matrix metalloproteinase-1 in cervical cancer. Int J Oncol, 52, 1198–1208.

46.

Truntzer

, Isambert

, Arnould

, Ladoire

, and Ghiringhelli

(2019). Prognostic value of transcriptomic determination of tumour-infiltrating lymphocytes in localised breast cancer. Eur J Cancer, 120, 97–106.

47.

Vacchelli

, Semeraro

, Enot

D.P.

, Chaba

, Poirier Colame

, Dartigues

, et al. (2015). Negative prognostic impact of regulatory T cell infiltration in surgically resected esophageal cancer post-radiochemotherapy. Oncotarget, 28, 20840–20850.

48.

Vigushin

D.M.

, Ali

, Pace

P.E.

, Mirsaidi

, Ito

, Adcock

, and Coombes

R.C.

(2001). Trichostatin A is a histone deacetylase inhibitor with potent antitumor activity against breast cancer in vivo. Clin Cancer Res, 7, 971–976.

49.

Wang

, Jia

, Wang

, Zhang

, Tan

, Zhang

, and Cheng

(2014). The clinical significance of tumor-infiltrating neutrophils and neutrophil-to-CD8+ lymphocyte ratio in patients with resectable esophageal squamous cell carcinoma. J Transl Med, 12, 7.

50.

Wieczorek

, Jablonowski

, Tomasik

, Konecki

, Jablonska

, Gromadzinska

, et al. (2015). MMP, VEGF and TIMP as prognostic factors in recurring bladder cancer. Clin Biochem, 48, 1235–1240.

51.

Xia

, Su

, Shen

, Meng

, Yan

, Zhang

, et al. (2018). ANLN functions as a key candidate gene in cervical cancer as determined by integrated bioinformatic analysis. Cancer Manag Res, 10, 663–670.

52.

Zeng

Y.T.

, Liu

X.F.

, Yang

W.T.

, and Zheng

P.S.

(2019). REX1 promotes EMT-induced cell metastasis by activating the JAK2/STAT3-signaling pathway by targeting SOCS1 in cervical cancer. Oncogene, 38, 6940–6957.

53.

Zhang

A.W.

, McPherson

, Milne

, Kroeger

D.R.

, Hamilton

P.T.

, and Miranda

(2018). Interfaces of malignant and immunologic clonal dynamics in ovarian cancer. Cell, 173, 1755–1769.

54.

Zhang

Y.X.

, and Zhao

Y.L.

(2016). Pathogenic network analysis predicts candidate genes for cervical cancer. Comput Math Methods Med, 2016, 3186051.

55.

Zuazo

, Arasanz

, Fernández-Hinojal

, García-Granda

M.J.

, Gato

, Bocanegra

, et al. (2019). Functional systemic CD4 immunity is required for clinical responses to PD-L1/PD-1 blockade therapy. EMBO Mol Med, 11, e10293.

56.

Zuazo

, Arasanz

, Fernandez-Hinojal

, Gato-Canas

, Hernandez-Marin

, Martinez-Aguillo, et al. (2018). Highly differentiated CD4 T cells unequivocally identify primary resistance and risk of hyperprogression to PD-L1/PD-1 immune checkpoint blockade in lung cancer. bioRxiv. https://doi.org/10.1101/320176

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.01 MB

0.02 MB