Long non-coding RNA TP73-AS1 is a potential immune related prognostic biomarker for glioma

The expression of TP73-AS1 is correlated with poor clinicopathological features in gliomas

To explore the difference in expression of TP73-AS1 between gliomas with different clinicopathological features, we examined the expression of TP73-AS1 in three datasets: TCGA, CGGA and GSE16011. We found that the expression of TP73-AS1 increased with the increase of WHO grades (Figure 1A–1C). Age is a risk factor for glioma patients. We found that in the three datasets, the expression of TP73-AS1 in the older age group was significantly higher than that in the younger age group (Figure 1D–1F).

Considering the important influence of IDH mutation and 1p/19q chromosomal co-deletion on the prognosis of glioma patients, we evaluated the correlation between the expression of TP73-AS1 and the statuses of IDH and 1p/19q chromosome. We found the same result that TP73-AS1 expression was significantly downregulated in the IDH mutation group (Figure 2A–2F), and the 1p/19q chromosome co-deletion group in TCGA, CGGA and GSE16011 (Figure 2G–2L).

TP73-AS1 is a risk factor for the prognosis of gliomas

The expression of TP73-AS1 is higher in glioma patients with a worse prognosis. We therefore explored the direct relationship between TP73-AS1 expression and the prognosis of glioma. We found that TP73-AS1 was a significant prognostic risk factor, including for primary gliomas in the TCGA, CGGA, and GSE16011, and recurrent gliomas in CGGA (Figure 3A–3D). We further analyzed different WHO grade of gliomas. We found that for WHO grade III gliomas, the high expression of TP73-AS1 indicated a poor prognosis for primary gliomas in the TCGA, CGGA, and GSE16011 and for recurrent gliomas of CGGA (Figure 3I–3L); however, the relation of high TP73-AS1 expression with poor prognosis was not obvious in WHO grade II (Figure 3E–3H) and WHO grade IV (Figure 3M–3P) gliomas.

The methylation level of TP73-AS1 was negatively correlated with the expression of TP73-AS1, a favorable prognostic factor

In general, promoter methylation regulates gene expression. We analyzed the methylation levels of TP73-AS1 in gliomas with different clinicopathological features. We found that the methylation level of TP73-AS1 decreased with an increase in WHO grade in the TCGA database (Figure 4A). IDH mutation (Figure 4B) and 1p/19q co-deletion (Figure 4C) led to high levels of TP73-AS1 methylation. Moreover, the methylation level of TP73-AS1 was decreased in the older age group (Figure 4D). We confirmed that the methylation level of the TP73-AS1 promoter was negatively correlated with the expression level of TP73-AS1, suggesting that the methylation of the TP73-AS1 promoter led to a decrease in TP73-AS1 expression (Figure 4E). Survival analysis showed that the methylation level of TP73-AS1 was a favorable prognostic factor, especially in WHO grade III gliomas (Figure 4F–4I).

Functional enrichment analysis of TP73-AS1

To further study the function of TP73-AS1, we used the TCGA based database Linkedomics to perform a correlation analysis of TP73-AS1. Two heatmaps were constructed to illustrate the genes whose expression was most positively and negatively correlated with that of TP73-AS1 (Figure 5A). GSEA showed that TP73-AS1 probably participates in KEGG pathways, including those involving Staphylococcus aureus infection, complement and coagulation cascades, and Reactome pathways, including integrin cell surface interactions and immunoregulatory interactions between a Lymphoid and a non-Lymphoid cell. The genes whose expression was negatively related with that of TP73-AS1 were enriched in KEGG pathways, including those associated with the ribosome and nicotine addiction and Reactome pathways, including eukaryotic translation termination and cap-dependent translation initiation pathways (Figure 5B–5C).

To further explore potential functional pathways associated with TP73-AS1, we conducted a functional pathway analysis based on the most positively correlated genes identified in the abovementioned analysis (Pearson r>0.5, P<0.01) using the “org.Hs.eg.db,” “clusterProfiler,” and “enrichplot” R packages. These genes were mainly associated with extracellular structure organization, extracellular matrix organization, collagen fibril organization, and positive regulation of receptor binding in GO (Gene Ontology) categories (Figure 5D). KEGG pathway analysis indicated that the pathways in which the most highly-correlated genes were enriched involved focal adhesion, proteoglycans in cancer, ECM−receptor interactions, and leukocyte transendothelial migration (Figure 5F).

Tumor immune microenvironment analysis

The tumor microenvironment (TME), especially the tumor immune microenvironment, is important for the survival of tumor patients. To evaluate the populations of various immune cells, we employed single-sample gene set enrichment analysis (ssGSEA; https://www.genepattern.org/modules/docs/ssGSEAProjection/4) to determine the relative infiltration level of 28 immune cells (Figure 6A). We made two major observations: (1) The high TP73-AS1 group had more complex immune cell components than the low TP73-AS1 group, including cells exhibiting anti-tumor reactivity, such as activated CD4 T cells, activated dendritic cells, central memory CD4 T cells, central memory CD8 T cells, natural killer cells, and natural killer T cells, and cells involved in tumor suppression, such as immature dendritic cells, MDSCs, plasmacytoid dendritic cells, regulatory T cells, and type 2 T helper cells (Figure 6B). (2) The correlation analysis revealed that the infiltration level of these two categories of immune cells were positively correlated (Figure 6C). This result revealed a feedback mechanism that anti-tumor inflammation could promote the recruitment or differentiation of immunosuppressive cells. Correlation analysis showed that TP73-AS1 expression was positively correlated with the infiltration scores of most immune cells (Figure 6D). To further verify the results, we used the Timer database (https://cistrome.shinyapps.io/timer/) to analyze the infiltration scores of immune cells in glioma tissue. We came to the same conclusion, that the high TP73-AS1 group had more complex immune cell components (Figure 6E) and that TP73-AS1 was positively correlated with the infiltration scores of most immune cells (Figure 6F).

Data analysis from the hospital

We found that TP73-AS1 was differentially expressed in gliomas with different clinicopathological features in the TCGA, CGGA, and GEO glioma datasets. The overexpression of TP73-AS1 was associated with poor clinical features, including age, stage, IDH mutation status, 1p/19q co-deletion status, and overall survival. RT-PCR was performed to detect the expression of the TP73-AS1 protein in 67 glioma tissues of different grades from our hospital. The results showed higher expression of TP73-AS1 in glioma samples of WHO grade IV than in glioma samples of WHO grade IV (p=0.002) (Figure 7A), and high TP73-AS1 expression predicted a poorer prognosis than low TP73-AS1 expression (P=0.048) (Figure 7B); these findings agreed with our preliminary results. We also found that TP73-AS1 expression levels were significantly correlated with the WHO grade of gliomas (p <0.0001), IDH1 mutation (p=0.0035), and 1p/19q co-deletion state (p=0.0119) (Table 1).

Data collection

The expression data and corresponding clinical information of glioma patients were obtained from (1) the TCGA dataset (https://xenabrowser.net/), which included expression data and data on the clinicopathological features of 664 patients; (2) CGGA RNA-seq data (http://cgga.org.cn/), which included expression data and data on the clinicopathological features of 984 patients; and (3) GSE16011 (https://www.ncbi.nlm.nih.gov), which included expression data and data on the clinicopathological features of 258 patients. This study was approved by the ethics committee of the Zhejiang Cancer Hospital. Written informed consent was obtained from all participants in this study.

qRT-PCR

TRIzol reagent (Life Technologies) and RevertAid First Strand cDNA Synthesis Kit (Thermo Scientific) were used to extract Total RNA and synthesize cDNA according to the manufacturer’s protocol. qRT-PCR analysis was performed using the standard protocol of the SYBR Select Master Mix kit. GAPDH was included as an internal reference gene. The primers used were as follows:

TP73-AS1-F 5ʹ - CCGTGGTTAGCGTGGACATCAG – 3ʹ

TP73-AS1-R 5ʹ - TGGGTGACAGGGCAAGACTCC – 3ʹ

GAPDH-F 5ʹ - CATGAGAAGTATGACAACAGCCT – 3ʹ

GAPDH-R 5ʹ - AGTCCTTCCACGATACCAAAGT – 3ʹ

Differential expression analysis and survival analysis

To investigate the expression of TP73-AS1, we used the R packages ggplot2, ggpubr, and grid Extra in the R software (3.6.2) to study gene expression patterns in gliomas with different clinicopathological features. For survival analysis, the R packages survminer, survival, ggpubr, and gridExtra were employed. The age segmentation value of each dataset was the median age of patients in the data set.

Gene set enrichment analyses and enrichment analysis

Data from the TCGA dataset were used for gene set enrichment analyses and enrichment analysis. GO analysis was performed using the EnrichGO function in the clusterProfiler R package, with the following parameters: ont = “all” pvalue-Cutoff = 0.05 and qvalue-Cutoff = 0.05. KEGG analysis was performed using the EnrichKEGG function of the R package “clusterProfiler,” with the following parameters: keyType = “kegg” pvalue-Cutoff = 0.05 and qvalue-Cutoff = 0.05. Gene set enrichment analyses were performed using the Linkedomics online dataset (http://www.linkedomics.org/login.php).

Immune cell infiltration

Single-sample GSEA was used to evaluate the infiltration level of 28 immune cell types. Immune cell signatures were obtained from a published paper [24]. The infiltration levels of immune cell types were represented by an enrichment score in the ssGSEA analysis of the GSVA R package. The ssGSEA score was normalized to unity distribution, with zero as the minimum and one as the maximum score for each immune cell type.

Ethics approval and consent to participate

Additional informed consent was obtained from all participants whose identifying information was included in this article.

Data availability statement

The datasets generated and/or analyzed during the current study are available from the corresponding author upon reasonable request.