GRAP2 : A Novel Immune-Related Prognosis Biomarker in Cervical Cancer

Background : Immune infiltration of the tumor microenvironment offers unlimited possibilities for therapeutic strategies in cervical cancer, where GRAP2 is an adaptor protein engaged in diverse signal activations. However, uncertainty exists regarding GRAP2 ’s prognostic significance and its relationship to immune infiltration. Methods : The data on cervical cancer cases were downloaded from The Cancer Genome Atlas (TCGA) database. The ESTIMATE computational technique was utilized to calculate the amount of immunological and stromal components, which helped us to identify the differential expression genes (DEGs). Among them, GRAP2 was considered to be related to overall survival based on a protein-protein interaction network and a univariate Cox regression analysis. Thus, based on the Gene Expression Omnibus (GEO) and TCGA databases, we evaluated GRAP2 ’s influence on clinical prognosis. Furthermore, GRAP2 expression was analyzed by Gene Set Enrichment Analysis (GSEA). Finally, we used CIBERSORTx analysis to assess the proportion of tumor-infiltrating immune cells (TICs) and the connection between GRAP2 and the tumor immune microenvironment. Results : ES-TIMATEScore was associated with cervical cancer patient’s prognosis. There are 791 DEGs and 11 potential key genes were identified including GRAP2 . In survival analyses with clinical information, We found that the GRAP2 high expression group exhibited a significantly longer overall survival (OS) than the low expression group and that the gene expression gradually declined as the Federation of International of Gynecologists and Obstetricians (FIGO) stage and M classification increased. GRAP2 was strongly linked with immunity and metabolism, according to GSEA. Finally, we discovered that 11 different TIC types and GRAP2 expressions were linked. Conclusions : GRAP2 may be a novel immune-related prognosis biomarker in cervical cancer.


Introduction
The term cervical cancer is used to define the malignant neoplasm that arises from the uterine cervix, whereby it poses a threat to female health.Persistent high-risk human papillomavirus (HPV) infection is a primary etiological factor in cervical cancer.Although such an important mechanism has been identified (which makes the issue of cervical cancer one of the most preventable types of cancer), cervical cancer is persistently the fourth most common female cancer worldwide [1].In the United States, which is the highest quality medical resource owner, cervical cancer took the lives of 4152 women in 2019, half of whom were under the age of 50 [2].Furthermore, due to low prophylactic HPV vaccination rates, the number is larger in most low-to middle-income countries, wherein cervical cancer has a higher incidence and mortality [3,4].
Although cervical cancer patients can be categorized according to the International Federation of Gynecology and Obstetrics (FIGO) to estimate their tumor burden and prognosis, the therapeutic strategy varies according to the tumor stage and the status of the patient.Chemotherapy and radiation therapy are the standard therapeutic strategies.For women with low-risk, early-stage disease, conservative, fertility-preserving surgical techniques are now considered the standard of care [5].However, for advancedstage patients with metastatic or recurrent disease, the overall prognosis remains poor [6].Revolutionary immunotherapy has been able to provide hope for patients in this situation, and this therapy aims to recruit and "educate" immune cells to target cancer cells more effectively and specifically.The major strategy includes adoptive cell therapy (ACT) and reversing immunosuppression or effector T-cell suppression.Nevertheless, solid tumors have generally been less responsive to ACT [7], whereas only pembrolizumab (which is a PD1-blocking antibody) has been approved by the Food and Drug Administration (FDA) for PD-L1positive metastatic or recurrent cervical cancer; moreover, even the response rate for these patients is only ~15% [8].
A sophisticated tumor microenvironment (TME) offers the foundation for immunotherapy, which consists of extracellular matrix, immune cells, inflammatory cells, fibroblasts, and surrounding blood vessels.Several special biomarkers in the TME have also been indicated to improve the prediction of several cancer prognoses.Approximately 1 decade ago, immune infiltration in the TME was observed as a prognostic factor that should not be ignored [9].According to recent research, tumor-infiltrating immune cells (TICs) and stromal components play an important role in the development of breast cancer [10,11], colorectal cancer [12], and cutaneous melanoma [13].
In this study, we screened differentially expressed genes (DEGs) related to the TME and clinical characteristics, through The Cancer Genome Atlas (TCGA) database and Gene Expression Omnibus (GEO), as well as several bioinformatic algorithms and tools, such as ESTIMATE algorithms and CIBERSORTx.ESTIMATE is a tool to predict tumor purity by using expression data to estimate stromal and immune cells in tumor; and it generates three scores: StromalScore (records the presence of stroma in the tumor tissue) ImmuneScore (represents the infiltration of immune cells in the tumor tissue) ESTIMATEScore (infers tumor purity) [14].CIBERSORTx is a tool for the deconvolution of the expression matrix of human immune cell subtypes based on linear support vector regression [15].Further analysis confirmed that GRAP2 may be a novel and potential immune-related prognosis biomarker for cervical cancer.

Data Preparation
The transcriptome profiling data of patients with cervical cancer that we analyzed were retrieved from the TCGA database (https://portal.gdc.cancer.gov/) on May 30, 2022.The data included 309 RNA-seq files.The tumorcorresponding clinical characteristics of cervical cancer were also downloaded from TCGA, and the data included 304 cases of tumors.The GSE52903 cohort was downloaded from the GEO database (https://www.ncbi.nlm.nih.g ov/geo/), which includes 55 cervical cancer RNA-seq files and corresponding clinical characteristics.

Estimation of ImmuneScores, StromalScores, and ESTIMATEScores
All of the transcriptome data were analyzed by using the "ESTIMATE" R package (https://www.r-project.org,v1.0.13) via Rstudio (https://www.rstudio.com,v2022.02.2+485).The calculated immune, stromal, and ESTIMATE scores were observed to display a compo-nent ratio in TME immune cell infiltration.ESTIMATES scores consisted of both ImmuneScores and StromalScores.Higher scores indicated a larger number of components in the TME.

Survival Analysis and Clinical Relevance
The samples were divided into the "high scores" group and "low scores" group according to the medians of the ImmuneScores, StromalScores, and ESTIMATEScores.A Kaplan-Meier survival analysis was performed in two groups via the R package "survival" (https://www.bioconductor.org/,v3.3-1) and "survminer" (https://www.bioconductor.org/,v0.4.9).The primary prognostic endpoint was overall survival (OS).In addition, some clinical characteristic data were permitted to assess the associations with the scores, which included age, stage (defined by FIGO), histopathological grades, and Tumour, Nodes, and Metastasis (TNM) stage.These comparisons were carried out using the Kruskal-Wallis rank-sum test, with p < 0.05 considered statistically significant.

Protein and Protein Interaction Network and COX Regression Analysis
The protein and protein >0.4 and homo protein interaction (PPI) networks were constructed on the website STRING (https://cn.string-db.org/).All of the DEGs were screened with confidence scores, and interactions based only on experimental evidence were selected as the organism to build the network.The network was visualized using Cytoscape (https://cytoscape.org/,v3.8.0).Moreover, the modular analysis was performed using the Cytoscape tool cytoHubba (https://www.bioconductor.org/,v2.5.3), and the top 10 genes with at least 15 adjacent nodes were selected according to the multinetwork clustering method (PPI key genes).We also performed a COX regression analysis with the R package "survival" to explore the re-lationship between each DEG and OS.Similarly, p < 0.05 was considered the cutoff criteria, wherein we obtained 147 genes related to patient survival status (the COX key genes).Subsequently, the R package "VennDiagram" was utilized to portray the intersection of key COX and PPI genes.

Survival Analyses of GRAP2
The concerned gene GRAP2 was selected.Afterward, we summarized the correlation between clinical information and GRAP2 expression of patients in the TCGA database, Pearson's Chi-squared test, Fisher's exact test, and Wilcoxon rank sum test were used to test the pvalue.The expression of each gene was matched to clinical features using "ggpubr" (https://www.bioconductor.org/,v0.4.0) and "survival".For these comparisons, a Wilcoxon rank-sum was employed, and p < 0.05 was considered statistically significant.To validate the reliability of the results, we extracted the GSE52903 data of cervical cancer patients from the GEO database.Due to the fewer number of samples (only 55 samples), we used the R package "maxstat" (https://www.bioconductor.org/,v0.7-25) to calculate the best cutoff value of |log2(fold change)|, wherein the minimum grouping sample size was set greater than 25%, and the maximum sample size grouping was less than 75%.Based on the best cutoff value of 5.81137, the patients were divided into two high (>5.81137)and low (≤5.81137)groups.By using the R package "survival" we performed a Kaplan-Meier survival analysis to estimate survival differences between the two groups.

Gene Set Enrichment Analysis for GRAP2
To further explore the role of GRAP2 in cervical cancer, GSEA, a computational method that determines whether an a priori defined set of genes shows statistically, was performed by using GSEA ( http://www.gsea-msigd b.org/gsea/index.jsp,v4.2.3) software to analyze the gene matrix of GO and KEGG associated with the gene GRAP2 expression level.Moreover, p < 0.05 was considered to be significantly enriched.

Flowchart of the Study
The current study was performed by using the following analysis process (Fig. 1).First, we downloaded the RNA-seq profiles of 304 cervical cancer cases and the corresponding clinical information from the TCGA and GEO databases.Afterward, calculations of TIC ratios and immune and stromal proportions were carried out separately using the CIBERSORT and ESTIMATE algorithms.We identified DEGs between groups with low and high Im-muneScore and StromalScore.Simultaneously, clinical information and TME scores were integrated for univariate Cox regression analyses, and we found that 147 genes showed a significant association with prognosis (the Cox key genes).Subsequently, the DEGs were subjected to PPI network analysis and univariate Cox regression analyses, and 30 genes were selected by the cytoHubba plugin (the PPI key genes).Finally, an intersection analysis of Cox key genes and PPI key genes was used to obtain 12 prognosis key genes.The gene GRAP2 attracted our attention; therefore, GRAP2 was further analyzed, including expression, survival, and clinicopathological characteristics, as well as GSEA and correlations with TICs proportions.

ESTIMATE Score is Associated with Cervical Cancer Patient Prognosis
Cervical cancer patients' immune infiltration scores were evaluated for prognoses (Fig. 2).For each of the TME scores (ImmuneScore, StromalScore, and ESTI-MATEScore), half of the 304 tumor samples (with complete clinical data) were classified into the high score group (152 cases), and the other half of the samples were classified into the low score group (152 cases) based on each median score.Thereafter, Kaplan-Meier survival curves were constructed.
OS was better for patients with a high ESTI-MATEScore (log-rank test p = 0.024, hazard ratio (HR) = 0.58) than the low scores (Fig. 2C).However, the Stroma- HR, hazard ratio; CI, confidence interval.lScore and ImmuneScore alone were not significantly associated with OS (Fig. 2A,B).According to these findings, the entire cervical TME might play an important role in cervical cancer prognosis.

Associations of TME-Scores with Clinical Characteristics among Cervical Cancer Patients
We also analyzed the relationship between TME scores and clinical characteristics (Fig. 3).The latter parameter included age, FIGO stage, pathologic grade, and the TNM stage.Along with the progression of the M classification, the StromalScore and ESTIMATEScore notably declined (except for the ImmuneScore) (Fig. 3N,O).This result is different from that of a previous report [16] due to more samples being included and more cases being assessed by our team.These findings clarified that the TME of cervical cancer plays a crucial role in the invasion and metastasis of cancer.

Identification and Analyses of DEGs
The median TME score was used to divide patients into low and high TME groups, which can help us to identify a correlation between the TME and gene expression (Fig. 4).1131 upregulated DEGs and 632 downregulated DEGs were acquired from ImmuneScore.1299 upregulated DEGs and 51 downregulated DEGs were obtained from the StromalScore group.The intersection of the two groups, including upregulation and downregulation, was considered the final representation of the DEGs, which contained a total of 791 genes (749 upregulated genes and 42 downregulated genes).Moreover, these genes were used for the GO and KEGG enrichment analyses (Fig. 5).Based on the GO analysis, these DEGs are associated with immune functions, such as plasma membrane structure and leukocyte mediation (Fig. 5A-C).According to KEGG, DEGs are mainly involved in chemokine signaling pathways and cytokinecytokine receptor interactions (Fig. 5D).It appears that immune factors play a pivotal role in cervical cancer TME.

Identification of Prognostic Key Genes
PPI networks, which were based on DEGs established in the STRING database and shown by Cytoscape software, were used to explore potential associations between these factors (Fig. 6).Subsequently, the top 25 key PPI genes (namely, PPI-key genes) were screened via the cytoHubba plugin in Cytoscape (Fig. 6A,B).147 genes were significantly associated with the 791 DEGs using univariate Cox proportional hazard regression analysis (namely, COX-key genes).Thereafter, 11 genes existing in both groups of PPIand COX-key genes were defined as being prognosis key genes and were visualized in the Venn diagram (Fig. 6C).Finally, we selected GRAP2, which was the best hazard ratio gene, for further analysis (Fig. 6D).

GRAP2 Expression was Related to Survival and Clinicopathological Characteristics in Cervical Cancer
The GRAP2 gene encodes GRB2-related adapter protein 2, which is a 37 kDa protein in humans and is also known as a GRB2-related adaptor downstream of Shc (GADS).GADS is involved in Ras signaling [17].The latest research suggests that GADS is recruited to phosphorylated CD28 to activate T cells [18].Increasing evidence has demonstrated that GRAP2 seems to serve as a protective factor in lung cancer [19].Furthermore, in this study, we outlined the relationship between GRAP2 expression and clinical information in the TCGA-CESC database (Table 1) and observed that GRAP2 expression is closely associated with survival time and M stage in cervical cancer patients.Besides, the GRAP2 high-expression group possessed a relatively longer OS than the low-expression group (Fig. 7A).To avoid biased results due to a single database source, we selected 55 samples from GEO (GSE:52903) to verify this result.Similarly, we found that high GRAP2 expression predicted better OS (Fig. 7B).After conducting an analysis of GRAP2 combined with clinicopathological characteristics, we found that GRAP2 expression gradually decreased following the advanced FIGO stage and M classification    (Fig. 7C-H).It is widely recognized that GRAP2 expression is associated with improved prognosis, including survival rates and metastasis rates, in cervical cancer.

GRAP2 is likely Involved in the Immune Activation and Metabolism Regulation of Cervical Cancer
According to the abovementioned results, we believe that GRAP2 expression may play a role in the malignant progression of cervical cancer.To further verify our theory, GSEA was used to explore the enrichment of GRAP2 expression in function and pathway.Genes within the GRAP2 high-expression tended to be associated with immune processes, such as cell activation involved in immune responses, cytokine-mediated signaling pathways, and mononuclear migration (Fig. 8A).Moreover, in the KEGG pathway analysis, the group mainly included immune-related pathways, including cytokine receptor interaction, chemokine signaling pathway, antigen processing and presentation, and T-cell receptor signaling pathway (Fig. 8C).Synchronously, the GRAP2 lowexpression cohort was enriched in several different biological processes, including binding of sperm to the zona pellucida, establishment of mitochondrion localization, isoprenoid biosynthetic process, and negative regulation of nuclear division (Fig. 8B).Similarly, a wide variety of metabolic function pathways were enriched in low expression, which included the TCA cycle, glycerolipid metabolism, GPI anchor biosynthesis, and glyoxylate and dicarboxylate metabolism (Fig. 8D).In this figure, only the top ten items of pathways are shown in each section.All of the results illustrated that GRAP2 may serve as a promising indicator of immune and metabolic statuses in cervical cancer.

Relationship between GRAP2 and the Proportion of TICs
The CIBERSORT algorithm was applied to further investigate the interplay between GRAP2 expression and TICs.First, Among the 22 immune cell types, 11 TICs had significant differences between the GRAP2 expression groups (Fig. 9A).Coincidentally, 11 TICs were strongly associated with GRAP2 expression from the correlation and difference analyses (Fig. 9B).Surprisingly, the intersection of the analyses showed a perfect coincidence of the two analyses (Fig. 9C).These TICs contain CD8 T cells, resting memory CD4 T cells, activated memory CD4 T cells, follicular helper T cells, M0 macrophages, M1 macrophages, activated dendritic cells, resting dendritic cells, resting mast cells, activated mast cells, eosinophils.According to these results, GRAP2 is essential to the immune TME of cervical cancer.

Discussion
In this study, we found there is a close correlation between the microenvironment and the prognosis of cervical cancer patients.We also identified several immune-related prognostic genes from the TCGA database.Furthermore, an important role for GRAP2 was found in immune-related biological processes.
An essential element of tumor biology is the immune resistance of the tumor microenvironment, which refers to the mutual resistance of malignant tumors and tumor cells.[20].Herein, we used the ESTIMATE algorithm to estimate the purity of the tumor microenvironment's components based on gene expression to calculate stromal scores, immune scores, and estimate scores, thereby aiding to identify candidate immune TME-related biomarkers.The results suggested that a high EstimateScore led to better clin-ical outcomes, including survival rate and M classification.Specifically, some TME-related genes were potentially able to predict prognosis.The same conclusions could be drawn from other studies [10,21].Immediately thereafter, a total of 749 upregulated and 42 downregulated genes were detected between groups with high and low ImmuneScores and StromalScores, respectively.GO and KEGG analyses suggested that the vast majority of these DEGs were involved in immune-related processes, such as leukocytemediated immunity and the chemokine signaling pathway.To further investigate TME-related genes associated with prognosis, we constructed a PPI network by using Cytoscape based on the DEGs, and we selected the top 25 key PPI genes via cytoHubba plugin scores.Thereafter, 11 genes were defined as key prognostic genes in both groups of PPI-and COX-key genes.In the survival analyses, all 11 genes indicated a favorable prognosis.Interestingly, GRAP2 (log-rank p < 0.001, HR = 0.499, 95% confidence interval (CI) = 0.289-0.862)was identified as the best prognostic gene among the 11 genes.For cervical cancer, this seems to be the first report describing the idea that GRAP2 may be a prognostic biomarker.Therefore, further investigation was framed around GRAP2.
The GRAP2 gene encodes the GRB2-related adaptor protein 2, a protein involved in leukocyte-specific signaling by protein-tyrosine kinases [22].This protein can interact with other proteins, such as GRB2-associated binding pro-tein 1 (GAB1) and the SLP-76 leukocyte protein (LCP2), to activate T cells via multiple signaling pathways [23].However, uncertainty surrounds the process of GRAP2 in cervical cancer.In this study, we found that patients with cervical cancer who expressed a high level of GRAP2 had an earlier FIGO and M stage, and a better survival ending.This suggests that GRAP2 may be used as a new prognostic biomarker to guide clinical treatment strategies.Similar findings were confirmed in another study [19], wherein GRAP2 was shown to predict survival status in lung adenocarcinoma patients.Furthermore, it has been claimed that GRAP2 is linked to metastasis in cervical cancer.Unfortunately, there is no research in this field to support this idea.
In the current study, we analyzed the components associated with GRAP2 in the cervical cancer microenvironment.Both GO and KEGG analyses, as well as GSEA, showed that GRAP2 was associated with immune function.This is consistent with the known function of GRAP2.For instance, GRAP2's carboxy-terminal SH3 structural domain constitutively binds SLP-76, whereas the SH2 structural domain of GRAP2 can bind to the tyrosine-phosphorylated linker for activation of T lymphocytes (LAT).A heterotrimeric complex made up of LAT, GRAP2, and SLP-76 mediates TCR signaling, thus linking proximal membrane signals to downstream signaling pathways and activating T cells as a result [23,24].Interestingly, we also found that the low expression group of GRAP2 sig- nificantly enriched several metabolism-related pathways in GSEA, such as glycerolipid metabolism.It is not difficult to understand this phenomenon, as immune cells are also key to the microenvironment and tumor metabolic reprogramming [25,26].Although there is little proof that GRAP2 is linked to tumor metabolic reprogramming, this could be a fruitful area of future investigation.This characteristic of GRAP2 provides new insights into immunotherapy for cervical cancer.There is a critical clinical need for new medicines or therapeutic approaches for advanced cervical cancer, notwithstanding the significant contribution of HPV vaccinations to cervical cancer prevention.Regardless of the second-line drug that is employed, the outcomes are quite dismal for patients who relapse after receiving first-line platinum-based chemotherapy [27].Immunotherapy is an attractive treatment strategy for cervical cancer, which is closely associated with human papillomavirus infection and is often accompanied by high PD-L1 expression and T-cell infiltration.Nevertheless, response rates with immune monotherapy remain under 15%, despite mounting evidence in favor of immunotherapy for cervical cancer [28].The TME is a critical factor limiting the efficacy of immunotherapy in cervical cancer [20].Another reason for the suboptimal therapeutic efficacy of adoptive T-cell therapies (ACT) in cervical cancer is the influence of "cold tumors", which are cancers that have poor immune activation and lymphocyte infiltration, thus making it challenging to kill tumor cells by activating immune effector cells to accomplish this feat [29].Herein, we found that GRAP2 plays an essential role in cervical cancer, particularly in the immunological microenvironment.We have reason to believe that GRAP2 has a promising future as an immunotherapy target and prognostic biomarker in cervical cancer.

Conclusions
GRAP2 demonstrates immunity-related activity in the TME of cervical cancer, thus suggesting that it may influence the biological behavior of cervical cancer, particularly in regard to metastasis.As a result, GRAP2 may be a unique molecular type of indicator, and it can be combined with the FIGO stage to enhance the prediction of clinical outcomes.An improvement in the identification of immunotherapysensitive cervical cancer patients will aid in the advancement of precision medicine.

Fig. 3 .
Fig. 3. Associations of ImmuneScore, StromalScore, and ESTIMATEScore with age, FIGO stage, grade, and TNM staging.(A-C) Kruskal-Wallis rank sum test of distribution of ImmuneScore, StromalScore, and ESTIMATEScore in the age classification, (D-F) in the FIGO stage, (G-I) in grades, (J-L) in T classifications, (M-O) in N classifications, (P-R) in M classifications.

Fig. 4 .
Fig. 4. Distribution of DEGs: Volcano plots, Venn diagrams, and Heatmap.Volcano plots showing the significance of DEGs between groups with high and low (A) ImmuneScore and (D) StromalScore.Venn diagrams showing the DEGs that were matched by the ImmuneScore and StromalScore assessments and were (B) upregulated and (E) downregulated.Heatmap for DEGs produced by contrasting the (C) ImmuneScore and (F) StromalScore groups with high and low scores.The gene name appears in the heatmap's row, while the sample IDs that aren't displayed in the plot appear in its column.

Fig. 5 .
Fig. 5. GO and KEGG enrichment analyses for 791 DEGs.(A-C) GO and (D) KEGG enrichment analyses, terms with p and q < 0.05 were thought to be considerably enriched.

Fig. 6 .
Fig. 6.PPI network visualization with univariate Cox regression analyses.(A) An experimental evidence type was selected to build the PPI network in the STRING database.The size and color of the node represent the degree value, the color of the edge represents the total score value, and the edge represents the interaction between the target gene.The degree value that corresponds to the color from yellow to blue increases with the size of the node.(B) Bar-plot diagrams of the top 25 genes identified by Cytoscape's cytoHubba plugin.(C) A Venn diagram showing where the important COX and PPI genes overlap.(D) The GRAP2 gene is highlighted in red on the Forest plot, which shows the substantial deferentially altered genes between the two cohorts.

Fig. 8 .
Fig. 8. GSEA for GRAP2 expression.Up-regulated genes are situated on the left, which approaches the origin of the coordinates; in contrast, down-regulated genes are located on the right of the x-axis.Each line represents one specific gene set with a distinct color.The only gene sets that were deemed significant have NOM p < 0.05 and FDR q < 0.25.(A) Top 10 significant GOBP associated with high GRAP2 expression.(B) Top 10 significant pathways associated with low GRAP2 expression.(C) The enriched gene sets in KEGG with high GRAP2 expression.(D) The enriched gene sets in KEGG with the low GRAP2 expression.

Fig. 9 .
Fig. 9.The proportions of TICs in groups with high and low GRAP2 expression and their correlation with GRAP2 expression.(A) By using the Wilcoxon rank sum test, 22 immune cell types in tumor tissues with low (green) or high (red) GRAP2 expression were compared.(B) An illustration of Pearson's correlation between the number of TICs and GRAP2 expression, with blue lines denoting the best fit linear models.(C) An intersection between analyses in (A,B) shows 11 TICs shared by both.