LINC01268 and CTB-31O20.2 as Prognostic and Functional Protective Immune Related lncRNAs in Glioblastoma

Glioblastoma (GBM) is the most common and deadly tumor in the central nervous system. Recent studies illuminated that long noncoding RNAs (lncRNAs) serve as competitive endogenous RNAs (ceRNAs) and play an important role in GBM by regulating immune responses. Here, GBM datasets from TCGA database were analyzed to obtain 356 signicantly differentially expressed lncRNAs (DE-lncRNAs), 4951 DE-mRNAs, and 34 DE-miRNAs in GBM, respectively. For mRNAs, 369 DE-mRNAs were identied as immune-related genes in Immport database. For DE-lncRNAs, univariate analysis identied 39 DE-lncRNAs with prognostic signicance, and 9 DE-lncRNAs are included in ImmLnc database. Then, combined analysis was conducted, by integrating 9 immune related DE-lncRNAs, 369 immune related DE-mRNAs and 34 DE-miRNAs, and generated a ceRNA network composed of 2 upregulated lncRNAs (LINC01268 and CTB-31O20.2), 3 downregulated miRNAs and 5 upregulated mRNAs. Then we focused on LINC01268 and CTB-31O20.2, and Kaplan-Meier survival, univariate and multivariate Cox regression analysis showed that LINC01268 and CTB-31O20.2 serve as independent protective prognostic markers in GBM. Finally, LINC01268 and CTB-31O20.2 overexpression was conducted in GBM cell U251. CCK8, transwell and scratch healing assay indicated that LINC01268 and CTB-31O20.2 inhibits GBM cell line U251 proliferation, invasion and migration. To sum up, LINC01268 and CTB-31O20.2 are independent prognostic immune related markers, and reduces cancer cell proliferation and metastasis in GBM.


Introduction
As one of the most feared type disease with poor prognosis, GBM accounts for 82% of primary malignant brain tumor and classi ed as grade IV tumors by the World Health Organization (WHO) [1]. Ample evidence suggested that the risk factors of GBM arising are cumulative genetic alterations and environment such as aging problem, family history, ionizing radiation, and virus [2,3]. However, the occurring of most GBMs are accidental and the causes are still unclear. Though some novel therapeutic methods such as molecular targeted therapy, gene therapy, and antiangiogenic agent treatment have appeared [4,5], the treatment of GBM remains a serious challenge. With the advent of advanced sequencing technologies, an abundant of gene expression data, distinct molecular, and immune-related genetic alterations have been revealed [6]. Through producing immunosuppressive factors and modulation of cell surface receptors and immune cell subsets, GBM impairs both local central nervous system (CNS) and systemic immune system functions [7]. Besides, GBM alters major immunogenic signaling pathways and altering cellular immunity inside and outside the brain [8]. A deeper understanding of the molecular biology of GBM can promote the development of various biomarkers and new therapeutic strategies.
Competitive endogenous RNA (ceRNA) hypothesis was rst proposed by Salmena and revealed a new mechanism for RNA interaction [9]. CeRNA competitively combined to microRNA through microRNA response elements (MREs) to regulate gene expression. Long non-coding RNA(LncRNA) was one of ceRNA, which was an RNA molecule greater than 200 bases in length [10]. With the advancement of highthroughput sequencing technology and the progress of abundant researches on lncRNA, the potential function of lncRNA has been revealed in many human disease [11][12][13][14][15][16]. Remarkably, a mass of RNA-Seq data helped identify a variety of biomarkers currently used for prognosis and treatment in cancer. For instance, Liang et al. identi ed prognostic six-lncRNA signature may improve prognosis prediction of   GBM, including C20orf166-AS1, LINC00645, LBX2-AS1, LINC00565, LINC00641, and PRRT3-AS1 [17].
P73-AS1 demonstrated that it could promote temozolomide resistance in glioblastoma stem cells, and revealed that high TP73-AS1 as a biomarker had a poor prognosis in glioblastoma [18]. LncRNA with prognostic value in GBM deserves our attention Many works on lncRNA as key regulator in complex immune response were emerging, which illustrated that lncRNA plays an important role in immune system [19,20]. Recently, the analysis of transcriptome identi ed many differentially expressed lncRNAs (DE-lncRNAs), which were involved in regulating the immune response. Such as SNHG14/miR-5590-3p/ZEB1 was pointed out that the positive feedback promoted large B cell lymphoma progression and immune evasion by regulating PD-1/PD-L1 checkpoints [21]. Moreover, lncRNA Sros1 accelerated IFN-γ-mediated activation of immune responses by stabilizing Stat1 [22]. With increasingly updated research, immune-related lncRNAs have been revealed in many type cancers. For example, OSTN-AS1 was an immune-related molecule and had a potential function of immunotherapy in triple-negative breast cancer [23]. The effect of lncRNA in immune regulation was complicated, and until now plenty of immune-related lncRNAs with prognostic role have not reported in GBM. Therefore, we focus on the prognostic role of immune-related lncRNA in GBM.
Here, we identi ed DE-lncRNAs, DE-mRNAs, and DE-miRNAs based on cohort study of GBM patients from TCGA database. Combining ImmLnc and ImmPort databases, we focus on the immune-related lncRNA and constructed an immune-associated ceRNA network. Then we identi ed the prognostic value of immune-related lncRNA in the ceRNA network. Our research identi ed 2 immune-related lncRNAs for prognosis of GBM. Finally, the function of LINC01268 and CTB-31O20.2 in GBM cell line U251 was investigated.

Data collection and preprocessing
The transcriptome data (include lncRNA, mRNA, and miRNA) and clinical data of GBM patients were downloaded from the Genomic Data Commons (GDC) data portal provided by TCGA database. A total of 155 GBM and 5 control samples were obtained. Then, 155 GBM patients were divided into high expression group and low expression group according to the 9 lncRNAs with prognostic value. The differentially expressed lncRNAs, mRNAs, and miRNAs between GBM and normal samples were analyzed by edgeR package of R, drawing volcano plot and heatmap by R script. The fold change > 2 and FDR < 0.05 were used for cutoff value.

Identi cation of immune-related DE-lncRNAs
ImmLnc database is a database developed by Yongsheng Li et al., which chart the landscape of immunerelated lncRNA regulation in 33 cancer types and can be used to identify potential carcinogenic biomarkers [24]. A total of 934 immune-related lncRNAs of GBM from ImmLnc database were downloaded as shown in Table 1, which used to recognize the potential function of immune-related DE-lncRNAs in GBM.

Identi cation of immune-related DE-mRNAs
ImmPort database contains a large number of immune-related genes and includes 17 immune-related pathways according to different molecular function [25]. Such as antimicrobials, antigen processing, presentation, BCR signaling, chemokine receptors, natural killer cell cytotoxicity, TNF family members, TGFβ family members, TCR signaling pathway. 1,424 immune-related genes were downloaded from the ImmPort database, which was identi ed immune-related DE-mRNAs in GBM.

Construction a ceRNA network
Database of LncBase [26] and miRDB [27] were used to predict corresponding target lncRNAs and mRNAs of miRNAs, respectively. Only immune-related lncRNAs and mRNAs that co-target miRNAs were eligible for the construction of ceRNA networks. Finally, immune-related lncRNAs, miRNAs, immune related mRNAs, immune-related pathways, and immune cells were selected to construct the immunerelated ceRNA network using Cytoscape visualization software.

Cell Culture
The human glioblastoma cell line U-251 was purchased from the Institute of Biochemistry and Cell Biology of the Chinese Academy of Sciences, Shanghai, China. U251 cells were cultured in Dulbecco's modi ed Eagle's medium (DMEM) (Gibco, Thermo sher Scienti c, Waltham, MA, USA). All cells were incubated in a humidi ed atmosphere containing 5% CO2 at 37°C.

LncRNA LNC01268 and CTB-31O20.2 overexpression
Full length of LNC01268 and CTB-31O20.2 was reverse transcribed into cDNA, PCR-ampli ed and cloned into the pCDNA3.1 vector (Invitrogen, Shanghai, China). The empty pcDNA3.1 vector was used as the control. The cell lines were transfected with a Lipofectamine 2000 transfection kit (Invitrogen, the United States) in accordance with the kit instructions. To detect the e ciency of transfection, the expression of LNC01268 and CTB-31O20.2 was veri ed by qRT-PCR.

Cell Proliferation Assays
Cell proliferation was measured using CCK-8 assay. After 48 hours of transfection, cells were seeded in 96-well plates with a density of 5×10 3 per well and cultured in an incubator for 0, 24, 48 and 72h at 37°C.

Scratch Wound Assay
Cell migration ability was evaluated using a scratch wound assay. U251 cells were seeded and cultured in 6-well plates for 24 h. Then a pipette tip was used to create wounds in culture plate. The wound healing process were observed at the time of 0h and 48 h.

Transwell Assay
Transwell assay were used to investigate the invasion ability of U251 cells. U251 cells were plated in the upper side of a Transwell chamber (Costar, USA), and 20% FBS medium was added into the lower compartment. After cultured and washed twice with PBS, the chamber was xed with methanol for 30 min, stained with 0.1% crystal violet. The cells were counted using a microscope.

Statistical analysis
Kaplan-Meier survival analysis, nomogram analysis, univariate and multivariate Cox regression analysis were performed by SPSS software (version 22) and R software (version 3.2.3). All statistical analyses were performed using GraphPad Prism (version 8). All cell experiments are replicated three times and Pvalue < 0.05 was considered statistically signi cant.

Identi cation of DE-mRNAs, DE-miRNA, DE-lncRNAs in GBM.
The impact of DE coding genes is often the most intuitive and basic change of molecular, because many coding genes have an important function in human diseases. In order to intuitively clarify the expression changes of mRNA, we used gene expression pro le of 5 normal samples and 155 GBM samples from TCGA database for identifying the DE-lncRNAs. As shown in Fig. 1A, clustering of 4951 DE-mRNAs of high or low expression level were exhibited and the Fig. 1B showed distinct cluster distribution of these DE-mRNAs in normal and GBM samples. To construct immune-related ceRNA network, 34 DE-miRNAs in normal and GBM groups were also identi ed (FC > 2, FDR < 0.05, Fig. 1C). The results of heatmap demonstrated that these DE-miRNAs were hierarchical clustering in each group (Fig. 1D). Furthermore, 356 DE-lncRNAs were selected (FC > 2, FDR < 0.05, Fig. 1E-F).

Identi cation of immune-related pathways in GBM.
To explore the immune-related pathways in GBM, we obtained 369 immune DE-mRNAs by overlapping the 4951 DE-mRNAs and immune-related mRNAs of ImmPort database. As shown in Fig. 2A, 239 upregulated DE-mRNAs and 130 downregulated DE-mRNAs were identi ed. Figure 2B showed the clustering of these genes in all samples. Afterwards, we built a regulatory network between genes and immune-related pathways, and 16 immune-related to immune pathways were enriched in all. They were cytokines, chemokines, chemokine receptors, cytokine receptors, and interleukins receptors, etc. as shown in Fig. 2C. 3.3 Identi cation of immune-related DE-lncRNAs with prognostic signi cance in GBM Next, we found that 39 DE-lncRNAs can be used as independent factors which related to prognosis by univariate analysis, and the general description was showed in Fig. 3A. Then, the cluster heatmap of these 39 DE-lncRNAs in the normal group and the tumor group was shown in Fig. 3B. Considering that the expression of lncRNA was related to immune cell in ltration in tumors, we took the intersection of common genes by used the lncRNA and 39 DE-lncRNAs related to glioblastoma in the ImmLnc database Fig. 3C, and constructed a network of lncRNA-immune cell as shown in Fig. 3D. In the immune cell-related lncRNA network, a total of 7 DE-lncRNAs was up-regulation, 2 DE-lncRNAs was down-regulation.

Construction of immune-related ceRNA networks
To construct immune-related ceRNA network, database of LncBase and miRDB were used to predict corresponding target lncRNAs and mRNAs of DE-miRNAs, respectively. Only immune-related lncRNAs and mRNAs that co-target miRNAs were eligible for the construction of ceRNA networks. Searching the LncBase database and miRDB database to predict the targeting relationship of lncRNA-miRNA and miRNA-mRNA, respectively. We integrated all the information of the LncBase database, miRDB database, ImmLnc database, and ImmPort database, and constructed an immune-related ceRNA network as shown in Fig. 4A. This network contained 5 immune-related pathways, which were TNF family members, cytokines, antimicrobials, TCR signaling pathway, and cytokine receptors. And included 5 DE-mRNAs that they were TNFSF14, SLC11A1, NCK1, NOX4, and NR5A2. Furthermore, LINC01268 and CTB-31O20.2 were obtained from network, has-miR-23b-5p and has-miR-139-3p were also acquired.

Identifying LNC01268 and CTB-31O20.2 as immunerelated prognostic markers for GBM
To further determine whether the CTB-31O20.2 and LINC01268 could be acted an independent predictor of GBM patients, we carried out Cox regression analysis and multivariate. the univariate Cox regression analysis demonstrated LINC01268 and CTB-31O20.2 expression level were associated with the overall survival (OS) in the Fig. 4B. The multivariate analysis showed Gender, CTB-31O20.2, and LINC01268 were signi cantly associated with the OS in the Fig. 4C. Subsequently, a nomogram was constructed to project the 1-year, 3-year, and 5-year OS of GBM patients as shown in Fig. 4D.

Discussion
Over the past decade, primary GBM is known for its aggressiveness and resistance to treatment, and because of its high recurrence, patients often die from the disease. Many novel therapies have emerged with the development of medical technology, such as surgery, chemotherapy, radiotherapy, molecular targeted agents, immunotherapy and nanotechnology [3,28]. As we known, GBM alters major immunogenic signaling pathways and altering cellular immunity inside and outside the brain [8].Tumorrelated immune responses and immunotherapy have been used in the treatment of malignant cancer [29].With the advent of targeted therapy and immunotherapy, the discovery of new biomarkers is urgent. In addition, the advancement of high-throughput sequencing technology made it possible for more and more function of non-coding RNA to be explored. Especially the study of some speci cally expressed genes as a biomarker for treatment, diagnostic and prognostic in all kinds of tumors [30][31][32]. However, there are still few biomarkers used for diagnosis and prognosis in GBM, and further needed to explore.
Many researches have reported that lncRNAs not only regulate gene expression in some tumors, but also can be used as biomarkers for diagnosis and prognosis. For example, Qi Sun etc. found that LOXL1-AS1 upregulate the expression of USF1 as a ceRNA via sponging miR-708-5p [33]. LncRNA RPPH1 was signi cantly upregulated and was associated with poor prognosis in colorectal cancer [34]. LncRNA HOTTIP could mediate HOXA9 to enhance the Wnt/β-catenin pathway by combining with WDR5 in pancreatic cancer stem cells, HOTTIP/WDR5/HOXA9/Wnt axis was expected to be potential therapeutic targets for pancreatic cancer [35]. These results have given us a newer understanding of lncRNA, but they are not su cient. In our study, we used the sequencing data of GBM and normal groups in the TCGA database to obtain DE-lncRNAs, DE-mRNAs, and DE-miRNAs through gene differential expression analysis. These DEGs serve as the basis for subsequent research.
Recently, with the increasing development of research on immune in ltration in the tumor microenvironment, immune-related genes have immediately become the focus of attention. Just as, 63 immune-related genes were founded associated with overall survival of melanoma and showed a powerful predictive ability in study of Rongzhi Huang etc. [36]. Wen Wang etc. revealed 9 immune-related lncRNAs have prognostic value for anaplastic gliomas [37]. Meng Zhou etc. veri ed that six-lncRNA signature could be an independent prognostic factor in GBM multiforme, with signi cantly different survival in high-risk and low-risk groups [38]. Therefore, it is necessary to explore the characteristics of immune-related molecules and evaluate the function of immune genes in GBM, which not only revealing the immune mechanism but also nding new therapeutic targets in GBM. In this study, immune-related mRNAs were obtained from ImmPort database and overlap with DE-mRNAs, 369 immune-related genes were selected. Furthermore, we took the intersection of 925 in the ImmLnc database and 32 prognosticvalue DE-lncRNAs, and obtained 9 immune-related lncRNAs related to the prognosis. Similarly, we identi ed 34 DE-miRNAs, all the above data were used in the following analysis. An article was illustrated AC064875.2 as potential prognostic biomarker may regulate neutrophil in ltration in glioma [39]. MIR210HG was reported up-regulation upon hypoxia exposure in glioma cells, and was considered a resistance in GBM therapy [40,41]. Transcriptome sequencing and qRT-PCR con rmed that lncRNA SMIM25 was high expression in cerebral cavernous malformations on the before study [42]. LINC01268 could in uence the emotional regulation, and was involved in gene regulation of potential immune response [43]. The alias of GS1-358P8.4 is PDK3, the study reported that high expressions of PDK3 were poor prognostic factors of acute myeloid leukemia [44]. Besides, RP11-429B14.4, CTB-31O20.2, RP11-436K8.1, and RP11-268J15.5 have not been reported, and further needs to explore in the future.
The mechanism of lncRNAs as ceRNAs in GBM have been attracted many researchers [45,46]. Therefore, we constructed an immune ceRNA network by predicting the targeting relationship between lncRNA-miRNA and miRNA-mRNA. In this ceRNA network, two lncRNAs were enriched. LINC01268 could regulate the expression of SLC11A1, NCK1, TNFSF14 and NOX4 by adsorbing hsa-miR-23b-5p and hsa-miR-139-3p, thereby affecting related immune responses. CTB-31O20.2 can regulate NR5A2 by adsorbing hsa-miR-139-5p to affect related immune responses. Otherwise, LINC01268 and CTB-31O20.2 could directly affect dendritic cell in the GBM. CTB-31O20.2 is a novel LncRNA that haven't been reported. LINC01268 and CTB-31O20.2 are novel LncRNA that have never been reported in GBM. Besides, the univariate, multivariate, and nomogram analysis indicated that these two lncRNAs are prognostic markers for GBM. The results above revealed the prognostic role of LINC01268 and CTB-31O20.2 in GBM via immune ceRNA network, which are helpful to have a better understanding of the underlying mechanism in immune responses and GBM.
In conclusion, our study screen out two innovative molecules, LINC01268 and CTB-31O20.2, which play a crucial role in GBM and have the potential to predict the prognosis of GBM. Additionally, LINC01268 and CTB-31O20.2 reduces cancer cell proliferation and metastasis in GBM. Our results also provide new insights into the discovery of biomarkers for the prognosis of GBM. However, in vivo experiments need to be further veri ed.

Declarations
Ethics approval and consent to participate: All participants are consent to this study.