WAKMAR2, a Prognosis-related Enhancer RNA in Gastric Cancer

Purpose An increasing number of long non-coding RNAs (lncRNAs) are thought to be associated with gastric cancer (GC). A lncRNA subclass that promotes enhancer function is called enhancer RNA (eRNA). We aimed to identify an eRNA that can predict GC prognosis and response to immune checkpoint inhibitors (ICIs). Methods Kaplan–Meier survival analysis was utilized to screen eRNA which can predict the prognosis of GC (P <0.05). The method of Spearman correlation analysis was employed in the ltration of target genes related to eRNA (r> 0.4, P <0.001). According to the median of WAKMAR2 expression, the patients were subdivided into low expression group and high expression group. Subsequently, differences of immune checkpoint-related genes and immune cell inltration between the two groups were further explored. Furthermore, we analyzed the correlation of WAKMAR2 with tumor mutation burden (TMB) and microsatellite instability (MSI) in GC and other types of cancer. Results WAKMAR2 and its target gene TNFAIP3 entered the subsequent analysis. Patients with high-WAKMAR2 expression had a favorable prognosis compared to patients with low-WAKMAR2 expression (P = 0.048). Immune checkpoint-related genes (PD-L1, CTLA4, PDCD1, LAG3) in the WAKMAR2 high-expression group were also highly expressed, except for B7-H3. In addition, inltration levels of B cells naive, T cells CD8, T cells CD4 memory activated, as well as Macrophages M1 in high-WAKMAR2 group were greater than in low-WAKMAR2 group. Last, the expression of WAKMAR2 in GC was signicantly correlated with TMB and MSI. Conclusion WAKMAR2, a new eRNA, is a promising biomarker that can be used to predict the overall survival (OS) of GC patients, and WAKMAR2 expression can be utilized to identify ICB responders in GC, providing new insights for immunotherapy strategies.


Introduction
Gastric cancer (GC) is a kind of malignant tumor with poor therapeutic effect, whose incidence and mortality rate shows a gradual upward trend [1]. Over the last few years, there has been major developments in GC diagnosis and treatment, but its prognosis is still poor [2]. At present, the most commonly used methods to evaluate the prognosis of GC are still histological diagnosis and TNM staging [3]. It is well described that GC is a highly heterogeneous disease [4]. Clinically, it is di cult to predict the prognosis of GC patients accurately by using the existing evaluation methods. Therefore, biomarkers that can accurately predict the prognosis of GC would help improve the clinical outcome of these patients.
There is mounting evidence showing that aberrant expression of long non-coding RNA (lncRNA) is closely related to numerous diseases, including cancer [5,6]. ERNA is a non-coding RNA transcribed from enhancer. Some recent studies have found that enhancer dysfunction is considered a key process of tumorigenesis [7,8]. Up to the present time, there is no research on GC related eRNA. Checkpoint inhibitors targeting cytotoxic T lymphocyte antigen 4 (CTLA4) and programmed cell death 1 (PD1) pave the way for the era of cancer immunotherapy, thus changing the way of cancer treatment [9,10]. However, although anti-PD-1 mAb is a promising approach for the treatment of patients with advanced GC, its response rate is still limited, so it is necessary to develop new strategies to maximize the e cacy of ICIs. However, although advanced GC patients can choose anti-PD-1 mAb treatment, the response rate is disappointing.
There is evidence to indicate the importance of the tumor microenvironment (TME) in tumor progression [11]. Tumor-in ltrating immune cells (TIICs) are one of the components of the TME, and it was a potential indicator to evaluate the therapeutic effect [12]. Many people had explored the link between TMB and immunotherapy response [13,14]. Higher TMB produce more neoantigens, the greater the possibility of being recognized by immune cells, so TMB is a good marker for predicting the response of immunotherapy [15,16]. Microsatellites (MS) are composed of repeating sequences of 1-6 nucleotides [17]. Its loss or gain is called MSI, which is also a powerful biomarker for immunotherapy [18].
In this study, we found that the expression of WAKMAR2 was related to the prognosis of GC, and it may be an immune-related eRNA. Its expression was related to the expression of immune checkpoint genes and the abundance of immune in ltrating cells. In addition, we also studied the correlation between WAKMAR2 expression and TMB and MSI in pan-cancer. This may help to accurately immunotherapy GC patients.

Data extraction
The gene expression, mutation data, clinical and survival data of GC and other type cancers were downloaded from the UCSC-Xena browser. From previous literature [19], we obtained a list of enhancer RNA and its predicted target genes.
Screening and identifying eRNA Patients were sorted into two groups according to the median level of eRNA expression. The Kaplan-Meier method was used to generate two groups of high and low expression survival curves. The difference of survival curve was determined by Log-rank test(p<0.05). The spearman method was used to analyze the correlation between eRNA and predicted target genes (correlation coe cient r > 0.4, P <0.001). It was considered as a candidate eRNA only if the following conditions were met, which was related to the overall survival of patients with GC (p <0.05) and the correlation coe cient between eRNA and target genes meets the conditions (r> 0.4, P <0.001).

GO and KEGG enrichment analysis
In addition to the predicted targets, other transcripts that were signi cantly related to eRNA were obtained through correlation analysis. In order to study the possible molecular mechanism of eRNA related coding genes, GO and KEGG pathway enrichment analysis was performed using the R package ("DOSE", "clusterPro ler", "enrichplot"). GO terms with both p-and q-value of <0.05 were considered as signi cant. KEGG terms with p < 0.05 were considered to be signi cantly enriched.

Relationship between eRNA expression and immune checkpoint-related genes
We extracted expression data of immune checkpoint-related genes (PD-L1, CTLA4, PDCD1, LAG3, B7-H3) from the expression matrix of GC. The Wilcoxon rank sum test was utilized to compare the expression difference of immune checkpoint-related genes between high and low expression eRNA groups.

CIBERSORT analysis
CIBERSORT is an analysis tool that can infer the expression matrix of 22 human leukocyte subtypes based on RNA-seq data. In order to quantify TIICs in GC samples, we use CIBERSORT algorithm. The Wilcoxon rank sum test was performed to compare the different abundance of immune in ltration between the two groups (eRNA high expression group and low expression group).

Correlation analysis of eRNA expression and TMB and MSI
We download the mutation data of 33 tumors from the UCSC website, and then calculate the TMB value of each sample in GC through a Perl script calculation. The MSI value of all tumor sample was obtained from the previous literature. Spearman correlation coe cient was calculated to assess the strength of the correlation between eRNA expression of GC TMB and MSI. In addition, we also conducted the above correlation analysis in other types of tumors.

Veri cation in pan-cancer
In order to determine whether eRNA expression can predict the OS of other types of tumors, we performed pan-cancer survival veri cation. In addition, we also veri ed the correlation between eRNA and target gene expression in pan-cancer.

Statistical analysis
Statistical analyses were conducted in R statistical package Version 3.6.3. Kaplan-Meier analysis was completed using R package 'survival'. The Wilcoxon rank sum test was used for comparison between two groups of clinicopathological parameters, and Kruskal-Wallis (K-W) test was used for three or more groups. All correlation coe cients were calculated by Spearman correlation analysis. For all analyses, two-tailed P-value < 0.05 was thought signi cant.

Results
Identi cation of eRNA associated with prognosis of GC We used a Perl script to convert the eRNA transcript ID into a gene symbol for subsequent analysis. Then, the expression of eRNA was extracted from the expression matrix of GC and combine it with the survival time. As showed in Table 1, 23 eRNAs signi cantly related to the OS of GC were ltered by Kaplan-Meier method (p < 0.05). Different from other eRNA, HAGLR has ve predicted target genes. Surprisingly, levels of these 23 eRNA were signi cantly correlated with their predicted levels of target gene mRNAs (r > 0.4, p < 0.001; Table 1).
LncRNA WAKMAR2 is a Key eRNA in GC LncRNA WAKMAR2 was selected as eRNA for further study, and its expression level was positively correlated with its predicted target gene TNFAIP3 level. Compared with patients in the WAKMAR2 high expression group, patients in the WAKMAR2 low expression group had a shorter OS ( Figure 1A, p < 0.05). As showed in Figure 1B, WAKMAR2 and TNFAIP3 mRNA levels are moderately correlated (r = 0.55, p < 0.001). It is worth noting that we studied the prognostic effect of WAKMAR2 in other cancer types and its correlation with TNFAIP3 mRNA levels. The impact of WAKMAR2 on OS and TNFAIP3 was speci c for 8 types of cancer only, which were GC, Adrenocortical carcinoma, Breast invasive carcinoma, Brain Lower Grade Glioma, Mesothelioma, Pancreatic adenocarcinoma, Pheochromocytoma and Paraganglioma and Thymoma (Table 2).

Association between WAKMAR2 expression and clinicopathological features
In order to verify the potential clinical utility of WAKMAR2 expression, the clinicopathological features of 371 GC patients were included in this study (Table 3). Figure 2 summarizes the associations between WAKMAR2 expression and clinicopathological features. Compared with Stage I and Stage II, WAKMAR2 expression level was higher in Stage III (Stage I vs III, p < 0.05; Stage II vs III, p < 0.01). However, compared with stage II and stage III, WAKMAR2 expression was lower in stage IV (Stage II vs IV, p < 0.05; Stage III vs IV, p < 0.001). The difference in WAKMAR2 expression level between Grade2 and Grade3 was statistically signi cant (G2 vs G3, p < 0.001). Furthermore, higher WAKMAR2 expression levels correlated with advanced T stage (T2 vs T3, p < 0.01; T2 vs 4, p < 0.01). However, no difference was observed between age, gender, N stage, M stage, family history, neoplasm status, pylori infection, and radiation therapy.

Functional enrichment analysis
To further elucidate WAKMAR2 function, we used correlation analysis to identify signi cantly coexpressed genes in GC. Including TNFAIP3, a total of 2335 transcripts were signi cantly associated with WAKMAR2 (p<0.001). Co-expressed genes underwent GO enrichment analysis to identify the functions ( Figure 3A). In BP category, "T cell activation", "T cell differentiation" and "lymphocyte differentiation" has been enriched, which means that co-expressed genes affect the function of the immune system in tumor microenvironment. Enriched CC terms included "phagocytic cup", "mast cell granule", "immunological synapse", and the enriched MF terms included "cytidine deaminase activity", "G protein−coupled chemoattractant receptor activity", "chemokine receptor activity". In addition, to determine the coexpressed gene enrichment pathway, we conducted KEGG enrichment analysis, including " Purine metabolism ", " Primary bile acid biosynthesis " and " Taurine and hypotaurine metabolism " ( Figure 3B). We analyzed the correlation between genes which were enriched in the "immune response-activating cell surface receptor signaling pathway" and WAKMAR2 (p < 0.001). Immune genes with Spearman correlation coe cient > 0.40 were listed in Table 4.
To evaluate the clinical signi cance of immune checkpoint-related genes The median value of WAKMAR2 expression was used as a cut-off value, patients were divided into high and low-expression groups. In this study, high-expression group patients exhibited higher gene expression of PDL1 (p < 0.001, Figure 4A), CTLA4 (p < 0.001, Figure 4B), PDCD1 (p < 0.001, Figure 4C) and LAG3 (p < 0.001, Figure 4D). However, the expression level of B7-H3 was higher in the low-expression group (p < 0.001, Figure 4e). The results of the study indicated that patients in the high-expression group were expected to be candidates for ICIs. Inhibitors against B7-H3 seem to have better therapeutic effects on patients in the low-expression group.

Differences in the abundance of immune cells between groups with high and low WAKMAR2 expression
Co-expressed genes have been shown to be involved in immune regulation, so we wanted to further analyze the differences in the immune fractions between the WAKMAR2 high-expression group and the WAKMAR2 low-expression group. The fractions of in ltrating immune cells were calculated by the  Figure 6B).

Discussion
ERNA is a special type of lncRNAs, which originate in the region of gene enhancers and can affect the transcription of corresponding genes through cis-acting. According to our screening criteria, WAKMAR2 was selected to enter the follow-up study. The original name of WAKMAR2 is LOC100130476, its genomic locus is very special, its transcription direction is antisense, and there is an overlapping region with the TNFAIP3 gene body and promoter part [20,21]. The expression of WAKMAR2 and TNFAIP3 has a strong correlation in various tumor types. A previous study con rmed that LOC100130476 was down-regulated in gastric cardiac adenocarcinoma tissue, suggesting that it had a tumor suppressive effect [22]. In addition, another study proved that LOC100130476 was signi cantly down-regulated in esophageal cancer cell lines and primary esophageal squamous cell carcinoma tissues, and it's up-regulation could inhibit the proliferation and invasion ability of cancer cells [23]. In other words, the survival rate of patients with cardia adenocarcinoma and esophageal cancer with low expression of LOC100130476 is poor, which is consistent with the results of our study.
As we all know, the higher the expression level of immune checkpoint-related genes in tumor tissues, the better the therapeutic effect of ICI. There is now accumulated evidence suggests that ICI have therapeutic potential in GC [24,25]. Nevertheless, only a few patients bene t from ICI. In addition, tumor immune microenvironment plays a crucial role in tumor progression. It is particularly important to identify precisely those GC patients who will bene t most from ICI. Very coincidentally, the expression level of WAKMAR2 is a potential marker to accurately distinguish these patients. We explored the association between WAKMAR2 and immune checkpoint-related genes for expression. The expression levels of PDL1, CTLA4, PDCD1 and LAG3 were signi cantly increased in WAKMAR2 high expression group. However, the level of B7-H3 increased signi cantly in the low expression group. This result provides strong evidence for our hypothesis that eRNA expression levels can be used to screen which patients are more likely to respond to ICI. Mast cells not only to promote the metastasis of GC, but also regulate the tumor microenvironment by releasing IL-17 [32]. Neutrophils promote migration and invasion of GC cells by activating the ERK pathway and inducing epithelial-mesenchymal transition, suggesting that neutrophils may play a crucial role in GC metastasis [33]. Therefore, this study provides a new understanding of the role of WAKMAR2 in the regulation of TIICs.

Immune-cell in ltration is
With the development of high-throughput sequencing technology, detecting TMB and MSI from genome sequencing is a new method to predict the e cacy of ICB, and it has been proven to be effective in a variety of tumors [34]. A potential problem in daily clinical practice is that TMB and MSI detection was expensive and time-consuming. Therefore, it is particularly important to explore alternative biomarkers for TMB and MSI. In the study of GC, the expression of WAKMAR2 has a certain correlation with TMB and MSI, and may be used as a potential marker to replace them. In the process of pan-cancer veri cation, it was surprisingly found that TMB of thymoma had a negative correlation with WAKMAR2, and the correlation coe cient was as high as 0.768. This result will provide new insights for immunotherapy of thymoma.

Conclusions
In conclusion, we have identi ed a new prognostic marker for GC, the enhancer RNA WAKMAR2. Patients in the WAKMAR2 high expression group may have a better prognosis when receiving ICI therapy.

Declarations
Ethics approval and consent to participate No permissions were required to use the repository data.

Consent for publication
Not applicable.

Competing interests
The authors declare that they have no competing interests.