MAPK Activated Protein Kinase 3 Is a Prognostic-Related Biomarker and Associated With Immune Infiltrates in Glioma

Glioma is the most common primary brain tumor that causes significant morbidity and mortality. MAPK activated protein kinase 3 (MAPKAPK3/MK3) is a serine/threonine protein kinase regulating various cellular responses and gene expression. However, the role of MK3 in tumor progress, prognosis, and immunity for glioma remains unclear. Here, we determined the expression and prognostic values of MK3. We further analyzed the correlation of MK3 expression with immune infiltrations by using the biochemical methods and bioinformatic approaches with available databases. We find that MK3 is aberrantly upregulated in glioma. In addition, the higher MK3 expression is closely linked to the poor clinicopathologic features of glioma patients. Importantly, MK3 expression is negatively correlated with the prognosis of patients with glioma. Mechanistically, we demonstrated that the correlated genes of MK3 were mainly enriched in pathways that regulate tumor immune responses. The MK3 level was significantly associated with tumor-infiltrating immune cells and positively correlated with the majority of tumor immunoinhibitors, chemokines, and chemokine receptors in glioma. Thus, these findings suggest the novel prognostic roles of MK3 and define MK3 as a promising target for glioma immunotherapy.


INTRODUCTION
Glioma is known as the most common primary brain and spinal cord tumor (1), with an estimated annual incidence rate of 6.57 per 100,000 persons in the United States (2). Despite their relatively rare rate, glioma usually causes significant morbidity and mortality because of the low treatment success rate and poor overall survival (OS) rate. Among the newly diagnosed glioma, glioblastoma (GBM) is the most malignant form of brain cancer accounting for about 50% (3), with a 5-year survival rate of approximately 5% (4). So far, the main therapeutic options for glioma patients include surgery resection (5), chemotherapy (6), radiotherapy (7), and immunotherapy (8). However, the glioma, especially GBM, still remains incurable with poor prognosis, which imposes enormous pressure on society, although some clinically relevant epigenetic biomarkers such as IDH status, O-6-methylguanine-DNA methyltransferase (MGMT) promoter methylation status, histone code, and chromatin organization have been used for classification of glioma and treatment decisions (9). Thus, searching for the novel specific biomarker for advancing the prognosis of glioma remains an imperative challenge.
A previous study has reported that MK3 could act as a reliable prognostic indicator for colorectal cancer patients (25) as well as regulate NK cell cytotoxicity and CD4 T-cell development (23). However, the expression, clinical significance, biological roles, and potential molecular mechanisms of MK3 in glioma have yet not been investigated. Here, we revealed that MK3 was aberrantly overexpressed in glioma tissues and cell lines. We reported that the MK3 expression was closely associated with the poor clinicopathologic features and prognosis of glioma. Through the function and pathway enrichment analyses, we demonstrated that the correlated genes of MK3 were mainly enriched in immune regulatory pathways. Finally, we recovered the close correlation of MK3 expression with immune infiltration, immune-related genes, and immune checkpoints in glioma. Together, our study deciphered the essential role of MK3 in glioma prognosis and tumor immunoregulation.

Human Tissues
We collected 92 glioma tissues that were identified by the pathologists according to the 2016 WHO classification criteria from the Department of Pathology of the Affiliated Hospital of Xuzhou Medical University between 2016 and 2017. We obtained the ethical review and approval from the institutional ethics committee of Affiliated Hospital of Xuzhou Medical University (ethical review no. XYFY2018-KL056-01).

Gene Expression Analysis
We analyzed the MK3 expression with the "Single Gene Analysis" module of Gene Expression Profiling Interactive Analysis (GEPIA) web (http://gepia.cancer-pku.cn/) and observed the differential expression of MK3 between tumor and adjacent normal tissues for the different tumors or specific tumor subtypes of The Cancer Genome Atlas (TCGA). RNAsequencing (RNA-seq) data of 689 GBM and low-grade glioma (LGG) tissues from TCGA dataset (https://www.cancer.gov/ about-nci/organization/ccg/research/structural-genomics/tcga) and 1,157 normal tissues from Genotype-Tissue Expression (GTEx) (https://www.genome.gov/Funded-Programs-Projects/ Genotype-Tissue-Expression-Project) were used for MK3 expression analyses. RNA-seq data of 413 glioma tissues (batch I) and 273 glioma tissues (batch II) datasets after deletion of incomplete data from the Chinese Glioma Genome Atlas (CGGA) dataset were also utilized for the expression analyses. In addition, we analyzed the expression of MK3 in glioma by using the Gene Expression Omnibus (GEO) database, GSE4290, and GSE7696 (https://www.ncbi.nlm.nih.gov/geo/). t-Test was used for the expression analyses by using R package "ggplot2."

Western Blotting
Cells were lysed with radioimmunoprecipitation assay (RIPA) extraction reagent [50 mM of Tris-HCl, pH 8.0, 150 mM of NaCl, 0.5% sodium deoxycholate, 1% NP-40, and 0.1% sodium dodecyl sulfate (SDS)] containing the protease inhibitors (Roche, Germany). Total protein was separated using 10% polyacrylamide gel electrophoresis and transferred to NC membrane (Millipore, USA). The membrane was blocked with 5% bovine serum albumin (BSA) in TBST for 1 h at room temperature and then incubated with the primary antibodies against MK3 (1:200, Santa Cruz, USA) and GAPDH (1:10,000, ProteinTech, USA) at 4°C overnight. Then the membrane was washed three times using 0.1% TBST buffer and incubated with the horseradish peroxidase (HRP)-conjugated secondary antibodies for 1 h at room temperature, followed by three times TBST buffer washing. Finally, the signal was detected with enhanced chemiluminescence reagent.

Immunohistochemistry
The 4-mm-thick glioma sections were stained with a specific primary antibody against MK3 (1:50, Santa Cruz, USA) at 4°C overnight. Then the sections were incubated with secondary antibody (ZSGB-BIO, China) for 1 h and DAB for 2-5 min at room temperature, followed by hematoxylin staining. The images were acquired by Olympus microscopy and scored by two experienced pathologists without knowing the patients' characteristics. Scores were calculated based on the intensity and percentage of positive tumor cells within the whole tissue, which were evaluated by using the German semiquantitative scoring method. The intensity score of cytoplasmic staining of each specimen was defined as follows: 0, negative; 1, weak; 2, moderate; and 3, strong. The percentage of positive cells was evaluated with 0.0% staining; 1, 1%-24 % staining; 2, 25%-49% staining; 3, 50%-74% staining; and 4, 75%-100% staining. The final immunoreactive score was calculated by the multiplication of the intensity scores and proportion scores.

Survival Analysis
We obtained the survival and clinical phenotype data of each sample from TCGA and CGGA datasets. The Kaplan-Meier (KM) curve analyses were conducted by the R package "survival" and "survminer." Cox analysis was used for survival analyses by using the R package "survival." The receiver operating characteristic (ROC) analysis was performed with the R package "pROC."

Univariate and Multivariate Cox Analyses
Univariate and multivariate Cox analyses were conducted by using the R package "survival." The WHO grade, IDH status, 1p/ 19q codeletion, primary therapy outcome, gender, age, and MK3 expression were included in these analyses.

Protein-Protein Interaction Network of MK3 Analysis
We searched the MK3 with the query of protein name ("MK3") and organism ("Homo sapiens") on the STRING website (https://string-db.org/). The basic setting parameters were the network type ("full STRING network"), meaning of network edges ("evidence"), active interaction sources ("experiments"), minimum required interaction score ("low confidence (0.150)"), and max number of interactors to show ("no more than 50 interactors in 1 st shell").

Gene Ontology and Kyoto Encyclopedia of Genes and Genome Analyses
Gene expression data of GBM and LGG in HTSeq-Counts were downloaded from TCGA website for further analysis. The correlated genes of MK3 were screened with Pearson's correlation coefficients (|r| >0.4 and p < 0.001) with R package "deseq2." Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genome (KEGG) analysis were conducted on MK3 correlated genes with R package "clusterProfiler" to identify the possible biological functions and signaling pathways affected by MK3. The biological process (BP), cellular component (CC), and molecular function (MF) were applied in GO analysis.

Gene Set Enrichment Analysis
The Gene Set Enrichment Analysis (GSEA) was performed to dissect the cancer-related pathways by using the MSigDB Collection (c2.cp.v7.2.symbols.gmt) of the clusterProfiler R package. The gene sets with normal p. adjust value and false discovery rate (FDR) q value both less than 0.05 were considered as significantly enriched.

Immune Cell Infiltration Analysis
The immune cell infiltration levels from tumors were quantified by the single-sample GSEA (ssGSEA) with the R package "GSVA" (27). The correlation between the MK3 expression and immunoinhibitors, chemokines, and chemokine receptors across human cancers was determined on the Tumor-Immune System Interaction Database (TISIDB) website (http://cis.hku. hk/TISIDB/index.php). The correlation analysis between the expression of MK3 and immunoinhibitors, chemokines, and chemokine receptor genes of the glioma was evaluated by using the R-package "ggplot2."

Tumor Immune Estimation Resource Database Analysis
The correlation of MK3 expression with immune checkpoint expression from glioma was analyzed on the Tumor IMmune Estimation Resource (TIMER) website (https://cistrome. shinyapps.io/timer/) and adjusted by the tumor purity.

Statistical Analysis
All gene expression data were normalized with log2 transformation. Statistical analyses were conducted by R (3.6.3). The survival data from the CGGA database were acquired by KM. The correlation analysis was evaluated by chi-square (c 2 ) test, Pearson's correlation, or Spearman's correlation analysis. Student's t-test was used to determine the statistical significance of differences between groups, and the differences between more than two groups were analyzed by one-way ANOVA or the Kruskal-Wallis test. A two-sided Mann-Whitney test was used for non-parametric data. p-Value < 0.05 was considered statistically significant. *p < 0.05, **p < 0.01, and ***p < 0.001.

The MK3 Expression Analysis in Glioma
To determine the MK3 expression in glioma, we analyzed the GEPIA database that contains the gene expression profile information across all tumor samples and paired normal tissues. The MK3 was aberrantly overexpressed in GBM (n = 163) and LGG (n = 518) tissues compared with the normal brain tissues (n = 207) that had a low level of MK3 mRNA expression among different normal human organs ( Figures 1A, B). We further validated the data by using the GEO database; a similar overexpression pattern was observed in the glioma tissues compared with the non-tumor (GSE4290) or normal (GSE7696) tissues ( Figures 1C, D). We also assessed the MK3 expression levels in glioma tissues by immunohistochemistry (IHC) staining. The staining results showed a similar MK3 expression pattern from the tumor tissues when compared with the para-tumor tissues ( Figures 1E, F).
Furthermore, we measured the protein expression levels of MK3 in five glioma cell lines (U118, U87, U251, T98G, and LN229), HEB, and NHA cell lines by Western blotting analysis. We observed that the protein levels of MK3 were also higher in the majority of the glioma cell lines ( Figure 1H). Collectively, our results suggest that MK3 is aberrantly overexpressed in glioma.

Correlation of MK3 Expression With Clinicopathologic Features in Glioma
The significant positive correlation between the MK3 level and WHO grade in glioma tissues (Figures 1E, G, Table 1) suggested a potential role of MK3 in the clinicopathologic features of LGG, low-grade glioma; HEB, human normal brain glial cells; NHA, normal human astrocyte. *p < 0.05, **p < 0.01, and ***p < 0.001.  Table 2, MK3 overexpression was significantly correlated with the WHO grade, IDH status, 1p/19q codeletion, and age. We also analyzed the expression of MK3 with OS, progression-free interval (PFI), and disease-specific survival (DSS) event of GBM and LGG, and the results showed that the expression of MK3 was lower from alive patients (Figures 2A-C). Moreover, we determined the clinical implication of MK3 expression from human GBM and LGG samples in TCGA database, and the significant differences were found in 1p/19q codeletion ( Figure 2D), WHO grade ( Figure 2E), and IDH status ( Figure 2F). The results were also validated with the CGGA database ( Figures 2G-L). Altogether, we demonstrate that MK3 expression is closely correlated with the poor clinicopathologic features of glioma patients.

The Prognostic Value of MK3 in Glioma
To investigate the relationship between the MK3 expression and the prognosis of glioma patients, we divided the cancer cases into high-risk and low-risk groups according to the cutoff value of the median risk score. The poorer prognosis and higher death rate were observed in the high-risk group ( Figure 3A). Glioma patients were also separated into the MK3-high and MK3-low groups based on the expression levels with the median value to generate a KM survival curve. The KM survival curve demonstrated that the OS rate of glioma patients in the highexpression group was significantly poorer than that in the lowexpression group in both the CGGA and TCGA datasets ( Figures 3B-D). The PFI and DSS rate of glioma patients were also significantly and negatively correlated with MK3 expression in TCGA dataset ( Figures 3E, F). In addition, we performed a subgroup survival analysis of OS with TCGA dataset, and we found that high expression of MK3 was associated with poor prognosis in WHO grade G3, non-codel 1p/19q codeletion, age less than or equal to 60 years, age greater than 60 years, and histological type astrocytoma subgroup of glioma We also conducted the ROC curve analysis to evaluate the diagnostic value of MK3, and the area under the curve (AUC) was 0.966, which indicated a high diagnostic value of MK3 in glioma ( Figure 3R).
To further evaluate the prognostic value of MK3, we performed the univariate and multivariate Cox regression analyses, which are the widely used approaches for identifying predictive biomarkers. Univariate Cox analysis unearthed that MK3 was a high-risk factor (hazard ratio (95% CI) = 3.037 (2.349-3.297), p < 0.001). Moreover, WHO grade, IDH status, 1p/19q codeletion, primary therapy outcome, age, and MK3 expression were all significantly correlated with poor OS (Table 3), which indicated their association with the OS in glioma. On the other hand, the multivariate Cox analysis demonstrated that only WHO grade, IDH status, primary therapy outcome, gender, and age were independent prognostic factors for glioma ( Table 3). In conclusion, these results demonstrate that MK3 might serve as a valuable prognostic biomarker for glioma patients.

Function and Pathway Enrichment Analyses of MK3 in Glioma
To identify the functional and physical interaction partners of MK3, we used STRING to analyze the protein-protein interactions. As shown in the interaction network ( Figure 4A), a total of 38 MK3 interaction proteins were found by experimental evidence. Furthermore, we recovered 3,067 positively correlated genes (Pearson's correlation coefficient r > 0.4) and 1,229 negatively correlated genes (Pearson's correlation coefficient r < −0.4) of MK3 from TCGA transcriptome dataset. The top 15 positively correlated genes and negatively correlated genes of MK3 were plotted in the heatmap ( Figure 4B). To further uncover the potential role of MK3 in glioma, we performed GO and KEGG analyses by TCGA dataset on MK3 correlated genes. Interestingly, the results of GO classification demonstrated that MK3 correlated genes were significantly enriched in immune-related functions among the top-ranked BPs, including the neutrophil activation, neutrophilmediated immunity, neutrophil activation involved immune  response, T-cell activation, regulation of innate immune response, leukocyte cell-cell adhesion, regulation of T-cell activation, leukocyte proliferation, antigen processing and presentation, response to interferon gamma, antigen processing and presentation of peptide antigen, cellular response to interferon gamma, T-cell proliferation, and interferon gammamediated signaling pathway ( Figures 4C, D). KEGG pathway analysis showed that MK3 correlated genes were also closely  We also performed GSEA by using the MK3-low and MK3high datasets from TCGA transcriptome dataset to identify signaling pathways that are affected by MK3 overexpression in glioma. Gene sets related to chemokine signaling, Toll-like receptor signaling, cytokine-cytokine receptor interaction, natural killer cell-mediated cytotoxicity, interferon signaling, interferon gamma-signaling, cell surface interactions at the vascular wall, neutrophil degranulation, extracellular matrix organization, signaling by interleukins, signaling by the B-cell receptor, and Toll-like receptor cascades pathways were significantly enriched in the high MK3 expression group ( Figures 5A-L). Together, our results suggest that MK3 might participate in the tumor immune microenvironment and immune regulation.

Correlation Between MK3 Expression and Immune Infiltration
Considering the tumor heterogeneity, we further analyzed the distribution of MK3 expression in glioma samples by assessing single-cell RNA-seq data. We observed that the extensive expression of MK3 was in immune cells (macrophage and T cell) (Figures 6A-C), which further support our IHC staining data in glioma tissues ( Figure 1E). These results suggest the potential role of MK3 in regulating immune response. Therefore, we determined the correlation between MK3 expression and immune infiltration in glioma by ssGSEA with Spearman's R. We found that MK3 expression was negatively correlated with infiltrating levels of T follicular helper (TFH), CD8 T cells, T gamma delta (Tgd), T central memory (Tcm), NK CD56bright cells, and plasmacytoid DC (pDC) cells and was positively correlated with infiltrating levels of macrophages, neutrophils, eosinophils, activated dendritic cells (aDC), immature DC (iDC), T, Th17, cytotoxic cells, NK CD56dim, Th2, NK, T helper, Th1, B cells, and DC (all p < 0.05, Figures 6D, E). These data indicate that MK3 might participate in immune infiltration in glioma.
Taking advantage of TISIDB website, we also explored the correlation between the expression of MK3 and immune-related genes in glioma. Our results suggest that the expression of MK3 was positively correlated to the majority of immunoinhibitors ( Figures 7A, B), pro-tumor chemokines ( Figures 7C, D), and chemokine receptors ( Figures 7E, F) in both GBM and LGG patients.
Programmed cell death protein 1 (PD-1/PDCD1), programmed cell death 1 ligand 1 (PD-L1/CD274), and cytotoxic T lymphocyte antigen 4 (CTLA-4) are vital immune checkpoints that play important roles in tumor immune escape. And they also served as predictive markers for the therapeutic efficacy of immune checkpoint inhibitors (ICIs). Thus, we especially assessed the relationship of MK3 with PD-1, PD-L1, and CTLA-4 in glioma by using TIMER tool. The results from TIMER website illustrated that there was a significant positive correlation between MK3 expression and PD-1, PD-L1, and CTLA-4 in glioma, which was adjusted by tumor purity (Figures 8A-F). Similarly, we also found that the expression of MK3 was significantly positively correlated with PD-1 (r = 0.520, p < 0.001), PD-L1 (r = 0.570, p < 0.001), and CTLA-4 (r = 0.360, p < 0.001) in glioma (Figures 8G-I). Taken together, these results suggest that MK3 might participate in tumorigenesis and the development of glioma by regulating the tumor immune escape.

DISCUSSION
MK2 and MK3 are MAPK-activated proteins that are often being classified as isoenzymes. They have a high degree of structural similarity, and they also share activators, substrates, and physiological functions (28). Both are ubiquitously expressed in various tissues of mammals. However, MK3 is mainly expressed in muscle, liver, heart, and T and NK cells (11). MK2 and MK3 usually cooperate to regulate a variety of cellular processes, such as cytokine production, endocytosis, cytoskeleton organization, cell migration, cell cycle control, chromatin remodeling, and gene expression (28), while the potential functions and effects of MK3 for tumorigenesis are poorly understood. One study suggested that MK3 could act as a reliable prognostic indicator in colorectal cancer patients (25), but the expression, clinical significance, and underlying molecular mechanisms of MK3 in glioma have yet not been reported.
In this study, we found that the MK3 was lower in the normal brain tissues and differentially expressed in many types of cancers. Importantly, MK3 expression was markedly elevated in GBM and LGG, which was further validated by using the GEO and CGGA database and the Western blotting approach in glioma cell lines. In addition, IHC staining results showed a similar result in glioma tissues (Figure 1). We also demonstrated that the overexpression of MK3 was closely associated with WHO grade, 1p/19q codeletion, IDH status, and age. We further validate that the higher MK3 expression correlated with the poor clinicopathologic features of glioma ( Figure 2). Moreover, the KM survival curve analyses showed that a high  level of MK3 was correlated with poor prognosis in glioma from both TCGA and CGGA datasets ( Figure 3). Together, our findings reveal the novel roles of MK3 for glioma and define MK3 as a valuable biomarker. The survival analyses found that high expression of MK3 was negatively correlated with the OS rate of glioma patients ( Figures 3B-D), while survival analyses of histological type for glioma patients showed that MK3 did not affect survival in GBM ( Figure 3O), which might be limited by the small number of patients enrolled due to the short survival time and low survival rate of GBM (29). More GBM cases were needed to further evaluate the prognostic value of MK3 in GBM in the future.
Multivariate Cox analysis demonstrated that gender was an independent prognostic factor for glioma ( Table 3). Sex differences have been well identified in many brain tumors including glioma. The sex-specific effects for the incidence, phenotype, and outcome of glioma have been well described; however, few insights are available to distinguish male and female glioma patients at the molecular level or allow specific targeting of these biological differences (30,31). More studies are needed to focus on sex differences in GBM in terms of pathophysiology, hormones, metabolism, tumor location, treatment response, recurrence, and outcome (31).
Our function and pathway enrichment analyses suggested that MK3 co-related genes enriched in interferon-gamma signaling pathway (Figures 4, 5), which is consistent with previous studies that MK3 regulated the transcription activity of type I interferon-dependent genes (20) and that MK3 could suppress interferon-gamma expression to control NK cell cytotoxicity and Th1 CD4 T-cell development (23). Interferon has been widely used in the clinic for cancer treatment, as the interferon signaling pathway regulates the immune checkpoint blockade and tumor microenvironment (32). This study also found that correlated genes of MK3 in glioma were mainly associated with tumor immune regulations. Moreover, singlecell RNA-seq data revealed that MK3 was expressed in immune cells ( Figures 6A-C). All these results indicated the potential roles of MK3 in regulating tumor immunity in glioma.
Although MK3 has been reported to express in NK and T cells, as well as regulate the inflammatory response of macrophages (33) and neutrophil recruitment (13), the exact functions and possible mechanisms of MK3 for regulating tumor immunity are unclarified. We found that MK3 expression was significantly correlated with infiltration of various immune cells, such as macrophages, neutrophils, Th17, CD8 T cells, NK CD56bright cells, and pDC ( Figures 6C, D). Immune infiltrations have been shown associated with the prognosis of a patient (34), which suggested that MK3 may regulate immune infiltration and, in turn, affects glioma prognosis. More work is still needed to explore the special role of MK3 in tumor immunity of glioma.
Immunoinhibitors, such as PD-1, PD-L1, and CTLA-4, have gained widespread attention, as they can serve as immune checkpoint targets in multiple cancers to block immunoinhibitory signals and enable to produce effective antitumor responses (35). Our study demonstrated that the expression of MK3 was positively correlated with PD-1, PD-L1, and CTLA-4 (Figure 8), suggesting a prominent role of MK3 in regulating the expression of immune checkpoints and immunotherapy. Future investigations of MK3 should advance the therapeutic efficiency for glioma.
Chemokines are the largest subfamily of cytokines that could recruit different immune cell subsets into the microenvironment via interactions with chemokine receptors, which regulate tumor progression and therapeutic outcomes (36). Here, we reported that MK3 expression was positively correlated with the majority of chemokines and chemokine receptors, such as CCL (2,5,8,18,22), CXCL12, CCR (1, 2, 5), CXCR2, CXCR4, and CX3CR1 ( Figures 7C-E), which play pro-tumor roles. Thus, our findings highlight the essential roles of MK3 in immune infiltrations in glioma.
In summary, we found that MK3 was significantly aberrantly upregulated in glioma. Through a series of comprehensive approaches, we demonstrated that the increased MK3 expression is strongly associated with clinicopathologic features, poor prognosis, and immune cell infiltration in glioma. Altogether, our results suggest that MK3 might serve as a valuable prognostic factor and a promising novel immunotherapy target for glioma.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/supplementary material. Further inquiries can be directed to the corresponding authors.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Affiliated Hospital of Xuzhou Medical University. The patients/participants provided their written informed consent to participate in this study.