Tumor microenvironment variation implicates immune alterations and correlates with prognosis in patients with glioma

Tumor microenvironment variation implicates immune alterations and correlates with prognosis in patients with glioma | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Tumor microenvironment variation implicates immune alterations and correlates with prognosis in patients with glioma Danlei Chen, Yi He, Zhiyuan Feng, Longsheng Dong, Junfeng Yang This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-3829624/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract To investigate the microenvironment composition of gliomas and the associated clinical significance, we analyzed single-cell RNA sequencing and bulk RNA-seq data from glioma samples. Cell trajectory analysis identified five trajectories with distinct cell states and corresponding trajectory-related genes (TRGs). TRG-based clustering segregated patients with glioma with different overall survival, clinicopathological features, immune infiltration status, and immune checkpoint gene (ICG) expression levels. Notably, a worse prognosis was seen in patients with a higher immune score, lower tumor purity, higher M0 macrophage and regulatory T (Treg) cell infiltration, and increased ICG expression. Further survival analysis and functional enrichment analysis revealed a close relationship between prognosis and ICG-associated immunosuppressive pathways. Candidate prognostic genes were obtained using WGCNA analysis and differential expression analysis. LASSO and multivariate regression analysis were used to establish a prognostic prediction model. The prognostic risk-scoring signature including 12 genes successfully predicted patient survival with acceptable AUC values. A nomogram was constructed to evaluate the contribution of the risk signature to patient prognosis. This study highlights the potential involvement of tumor microenvironment variation and immune alteration in glioma progression and establishes a TRG-based prognostic model to predict patient clinical outcomes. single-cell RNA sequencing immune response glioma prognosis in silico analysis Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Figure 8 Figure 9 Figure 10 Figure 11 Figure 12 1 Introduction As the most common primary central nervous system neoplasms, gliomas encompass a large group of heterogeneous, aggressive tumors. Even if clinicians conduct current standard treatments, including a combination of surgery, radiotherapy, and chemotherapy, most tumors will recur due to chemo- and radioresistance. Patients usually succumb quickly to the disease, with approximately 5% five-year-survival rates ( 1 – 3 ). Recent advances in genomics, immunology, and other disciplines have resulted in various experimental therapies, including targeted therapy and immunotherapy, with the promise of providing new therapeutic strategies for gliomas. However, the efficacy of these new therapies remains unsatisfactory ( 4 , 5 ). Therefore, it is important to further explore the internal mechanisms of gliomas to uncover new therapeutic targets and approaches. The brain is considered an immune-privileged site because it lacks a traditional lymphatic system and contains only a few antigen-presenting cells. Recent findings changed this concept, demonstrating that the brain can coordinate a robust immune response involving both the innate and adaptive immune systems ( 6 , 7 ). The glioma-specific immunosuppressive microenvironment contributes greatly to the poor prognosis of patients with glioma ( 8 ). Novel immunotherapies targeting the tumor microenvironment (TME) and its immune components have shown great promise for the clinical management of many extracranial solid tumors, including melanoma and renal cell carcinoma ( 9 , 10 ). However, preliminary clinical trials using checkpoint inhibitors and vaccine therapy for glioma have been largely disappointing. The failure was likely due to the highly immunosuppressive environment, systemic immunosuppression, local immune dysfunction, and high tumor heterogeneity ( 4 , 11 ). Thus, the TME and immune alterations are of great interest for immunotherapy and glioma management. Single-cell sequencing is capable of examining sequence information from individual cells, probing cellular differences with high resolution, and offering a better understanding of how individual cells function in their microenvironments ( 12 ). Single-cell resolution genomic information can help uncover the roles of genetic mosaicism and intratumoral genetic heterogeneity during tumor development or treatment responses ( 13 ). Given the heterogeneity of glioma, single-cell profiling of genomic alterations in glioma could help dissect TMEs and provide important insights into the biology of tumor-infiltrating immune cells, including cell heterogeneity, dynamics, and potential roles in disease progression and responses to immunotherapy ( 14 ). Here, we used single-cell RNA-sequencing (scRNA-seq) data from glioma samples to identify five glioma cell trajectories with distinct identities. Then, we used the trajectory-related genes (TRGs) for immune and survival analysis. Figure 1 shows the study workflow, where glioma cell clustering, cell state trajectories, and TRGs were first acquired from the glioma scRNA-seq dataset in the GEO database. TRGs were then integrated into the Chinese Glioma Genome Atlas (CGGA) dataset to perform patient clustering, survival analysis, and immune-related analysis. Additionally, clinically-relevant TRGs were acquired through weighted correlation network analysis (WGCNA) and utilized to establish a prognostic model. A prognostic risk scoring signature was obtained through Least absolute shrinkage and selection operator (LASSO) and multivariate Cox regression analyses. Then, the scores were evaluated using ROC and Kaplan–Meier (K–M) curve analyses and quantified using a nomogram. 2 Materials and Methods 2.1 Acquisition and processing of scRNA-seq and bulk RNA-seq data scRNA-seq data from 3589 glioma cells were obtained from the Gene Expression Omnibus (GEO; accession number GSE84465; https://www.ncbi.nlm.nih. gov/geo/) database. The R package ‘Seurat’ was used to process the scRNA-seq data and perform quality control. Cells with a gene number < 100, read number 5% were excluded from the analysis. The LogNormalize method was then used to normalize the data, and the top 1500 genes with highly variable features were identified by variance analysis. Bulk RNA-seq data from 702 and 693 glioma samples were downloaded from the TCGA (http://cancergenome.nih.gov/) and CGGA mRNAseq-693 datasets (http://www.cgga.org.cn/index.jsp). For differential expression analysis of glioma compared to normal tissue, mRNA expression files from 207 normal, 518 low-grade glioma, and 163 glioblastoma samples were downloaded from the GEPIA2 database (http://gepia2.cancer-pku.cn/#index). Differentially expressed genes were identified using |log 2 (fold change)| > 1 and adjusted P-value < 0.05. 2.2 Dimensionality reduction and cell annotation A false discovery rate (FDR) < 0.05 cutoff was used to identify the top principal components (PCs), which were then used for dimension reduction with the tSNE algorithm to obtain the major single-cell clusters (15, 16). Based on |log 2 (fold change)| > 0.5 and FDR < 0.05, marker genes in each cluster were identified, and the top marker genes from the clusters were visualized using a heatmap. The clusters were annotated using the ‘SingleR’ package. 2.3 Pseudotime and trajectory analyses Pseudotime and trajectory analyses of single glioma cells were performed in R with the ‘Monocle’ package. Differentially expressed genes between trajectories of different states were obtained by comparing each trajectory with all other trajectories as a whole (using |log2(fold change)| > 0.5 and FDR < 0.05 thresholds) and designated as TRGs. The ‘org.Hs.eg.db,’ ‘clusterProfiler,’ ‘enrichplot,’ and ‘ggplot2’ packages were used for TRG gene ontology (GO) annotation and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis. 2.4 TRG-based molecular clustering of patients with glioma from the CGGA dataset The ‘ConsensusClusterPlus’ package is an R package capable of providing quantitative and visual evidence of stability, which is used to measure the number of unsupervised subtypes (17). Based on the TRGs, ‘ConsensusClusterPlus’ divided the CGGA glioma cohort into different clusters. To determine the optimal cluster number, the K-means algorithm and cumulative distribution function (CDF) curve were utilized with 50 iterations and a maximum tentative clustering number of 9. The survival difference was examined with K–M analysis using the ‘survival’ and ‘survminer’ packages. The clinicopathological features in each molecular cluster were visualized using the ‘ggplot2’ package. 2.5 TME, immune cell infiltration, and immune checkpoint analysis For each sample, TME scores, including immune score, stromal score, and tumor purity, were calculated using the ‘ESTIMATE’ package. The infiltration of 22 immune cells in each sample was identified using CIBERSORT software (18). Meanwhile, the immune cell infiltration density in different glioma clusters was compared using the ‘limma’ package and is shown as histograms and rainbow plots. Then, 38 immune checkpoint genes (ICGs) (acquired from a previous study (19)) across glioma clusters were examined for differential expression in the CGGA dataset. K–M analysis was used to evaluate the prognostic value of immune cells and ICGs. The prognosis-related ICGs were then analyzed using STRING enrichment to explore their underlying interactions and mechanisms. The core genes of the interaction network were identified using the MCODE application in Cytoscape software. The most significantly enriched pathways are shown in a bubble plot generated using the ggplot2 package. 2.6 WGCNA and establishing prognostic risk scoring signatures WGCNA is a systematic method for identifying modules of highly correlated genes, summarizing such modules using eigengene modules, and relating modules to one another and to external sample traits (20). In this study, WGCNA was used to analyze TRGs and identify key modules related to clinical features, especially patient survival. DEGs between glioma and normal tissues were obtained from the GEPIA2 online platform (21). The TRGs in the survival-related WGCNA modules were intersected with the DEGs from GEPIA2 to obtain candidate genes, which were further used in univariate, LASSO, and multivariate Cox regression analyses to generate prognostic risk-scoring signatures. LASSO regression identified 24 prognosis-related genes, among which 17 genes passed the collinearity test. These 17 genes were used for multivariate analysis to obtain the final 12 risk genes. The coefficients of the risk genes were determined with the training dataset, and the total risk score was calculated as the sum of the products of prognostic gene expression multiplied by their respective coefficients. The mean value of the total risk score in the training dataset was used as the threshold for patient risk stratification (high-/low-risk). K–M analysis and ROC curves were consecutively used to evaluate the accuracy and prediction efficacy of the prognostic risk-scoring signature. Multivariate Cox regression analysis was used to examine whether the risk-scoring signature was independent of other clinical features. 2.7 Nomogram construction Using the CGGA dataset, clinical features (e.g., age, sex, tumor grade, recurrent status, IDH mutation status, and 1p19q codeletion status) and risk scores were combined into a nomogram with the ‘rms’ package to predict the 1-, 3-, and 5-year overall survival (OS) for each patient with glioma. 3 Results 3.1 scRNA-seq quality control and normalization In this study, 3589 glioma cells were obtained from the GEO database (accession number GSE84465). After quality control and normalization, the samples showed good RNA features and count distribution. No correlation was detected between sequencing depth and mitochondrial gene sequences (Supplementary Figure 1A-C). We observed a significant positive correlation between sequencing depth and total intracellular sequences (R = 0.29, Supplementary Figure 1D). A total of 19312 genes were analyzed, of which the top 10 genes with the highest variation between cells were labeled (Figure 2A). 3.2 Dimension reduction and cell trajectory analysis identified 5 glioma cell states Given the limited computing power, principal component analysis (PCA) was first used for preliminary dimension reduction of scRNA-seq data and showed good segregation among glioma cells (Figure S2a). All top 15 principal components (PCs) showed significant differences and were enrolled for further analysis (Figure S2b). Top 20 variable genes for PC1-4 were displayed with scatter plots and heatmaps in Figure S2c, d. Then, utilizing t-distributed Stochastic Neighbor Embedding (tSNE) algorithm, glioma cells were aggregated into 9 clusters based on the top 15 PCs, and a total of 1061 marker genes were identified by differential analysis, and the top 10 marker genes in each cluster were demonstrated with heatmap (Figure 2b-c). Afterwards, pseudotime trajectory analysis revealed that clusters 3, 4, and 5 were distributed in state 1 and 5, which located in the relatively earlier stages of differentiation. In contrast, the clusters 0 and 2 were enriched in state 3, the latter stage of differentiation. (Figure 2d-f) 3.3 Molecular function enrichment of trajectory-related genes By comparing each pseudotime state trajectory with the remaining trajectories combined, we identified 601 glioma TRGs, with 313 TRGs for state 1 and 391 TRGs for state 3 (TRGs of different trajectories might overlap with each other). GO and KEGG enrichment analyses of TRGs of all states and showed significant enrichment in many immune-related pathways, including T cell activation, leukocyte-mediated immunity, and MHC protein complexes. Given that states 1 and 3 represented the early and late differentiation states respectively, we further ran enrichment analysis for TGRs of them, which showed results similar to total TRGs (Figure 3). 3.4 Genotyping based on TRGs segregated CGGA patients with glioma into three clusters with different prognoses TRG-based consensus clustering divided CGGA patients with glioma into three clusters (cluster C1–C3) with good differentiation and sufficient samples in each cluster (Figure 4A). K–M analysis was then used to determine the clinical significance of the consensus clustering. This analysis showed that cluster C3 had the worst OS (Figure 4B, P < 0.001). Clinical relevance studies further revealed that patients in cluster C3 had significantly older age, higher glioma grades, increased tumor recurrence rate, and a higher proportion of IDH wildtype and 1p19q non-co-deletion genotypes. Each of these parameters was associated with a worse prognosis (Figure 4C). 3.5 Clustering based on TRGs segregated glioma samples with different tumor microenvironment compositions and immune cell infiltration TME analysis indicated that samples in Cluster C3 had the highest stromal, immune, ESTIMATE scores, and the lowest tumor purity (Figure 5A). Further analysis of immune cell infiltration showed that cluster C3 had higher infiltration levels of many immune cells, including CD8+ T cells, regulatory T cells (Treg), and M0 macrophages (Figure 5B, C). Subsequent survival analysis of immune cell infiltration showed that 4 of 22 immune cells were prognosis-related, including monocytes, resting dendritic cells, memory B cells, and M0 macrophages. Notably, M0 macrophage infiltration was considerably increased in cluster C3, which displayed a worse patient prognosis. In addition, cluster C3 showed lower monocyte infiltration, which was correlated with better survival (Figure 5D–G). This result implies that the worse prognosis of patients with glioma in cluster C3 might be partly due to enriched M0 macrophages. 3.6 ICGs correlated with the prognosis of patients with glioma and might contribute to the survival difference between clusters segregated by TRGs Differential expression analysis showed that 36 of 38 ICGs were significantly differentially expressed between the TRG-based clusters (Figure 6A). In addition, survival analysis indicated that 25 of 38 ICGs were prognosis-related. ICGs that were correlated with longer survival were less expressed in cluster C3 (IL12A, LAMA3, LDHB, and VTCN), whereas ICGs that were associated with worse prognosis were upregulated in cluster C3 (B2M, CD8A, CD28, CD40, CD40LG, CD80, CD86, CD274, CTLA4, HAVCR2, ICOS, ICOSLG, LDHA, LGALS9, PDCD1, PDCD1LG2, PTPRC, PVR, SIGLEG15, TNFRSF9, and TNFSF4). K–M plots of ICGs with P < 0.01 (except VTCN, P = 0.012) are displayed in Figure 6B. These prognosis-related ICGs were then analyzed using the STRING online platform. Prognosis-related ICGs were significantly enriched in pathways associated with immune tolerance and immunodeficiency (Figure 7A, B). The protein-protein interaction networks obtained from the STRING analysis were further processed with MCODE to identify the core ICGs, which were used to construct a core network of 17 ICGs (Figure 7C). These results suggest possible underlying mechanisms for survival differences between patient clusters and directions for further immunotherapy in gliomas. 3.7 Identification of candidate prognostic genes WGCNA showed that TRGs were divided into nine modules, six of which (MEblue, MEturquoise, MEgreen, MEpink, MEbrown, and MEyellow; 461 genes in total) showed a significant correlation with both survival time and survival status (P < 0.05) (Figure 8A, B). Differential gene analysis using the TCGA glioma dataset identified 5745 and 7659 differentially expressed genes (DEGs) in low-grade glioma and glioblastoma, respectively, compared with normal tissue (Figure 8C, D). The union of those DEGs (8881 genes in total) were intersected with TRGs in the prognosis-related modules. Together, 381 intersecting genes were designated as candidate genes to construct the prognostic model (Figure 8E). 3.8 Establishment and evaluation of prognostic risk scoring signature The candidate genes were tested using univariate analysis, which showed that 360 of 381 candidate genes were associated with prognosis (Supplementary Figure 3). Given the numerous genes obtained by univariate analysis, we used LASSO regression analysis to filter the genes for downstream multivariate analysis. LASSO regression showed that 24 genes were correlated with survival. After further multivariate regression analysis, 12 genes were identified as independent risk genes (Figure 9) and were used to establish the prognostic risk scoring signature. Based on the score signature and corresponding gene expression in each sample, the TCGA-Low grade glioma (LGG) and glioblastoma (GBM) datasets were used to train the prognosis prediction model and obtain the coefficients for each risk gene. Then, the TCGA training cohort and CGGA testing cohort were divided into low- or high-risk groups according to the total risk score. Survival analysis showed that patients in the high-risk group had a significantly shorter OS (Figure 10A, B). ROC analysis indicated that the model had excellent sensitivity and specificity for predicting the 1-, 3-, and 5-year survival rate in both the TCGA and CGGA cohorts (Area Under Curves (AUCs) of 0.878, 0.925, and 0.896 in the TCGA cohort and 0.781, 0.789, and 0.768 in the CGGA cohort for 1-, 3-, and 5-year predicted survival, respectively; Figure 10C, D). Further univariate and multivariate analyses incorporated other clinical features and showed that the risk score is a prognostic risk factor for patients with glioma independent of other clinical characteristics, such as age, tumor grade, and recurrence status (Figure 10E–H). In addition, the expression levels of 12 risk genes in 9 single-cell clusters were plotted (Figure 11). ZCCHC24 and SULF2 were overexpressed mainly in cluster 3, SERPINE2 was chiefly expressed in cluster 6, and MT1X was expressed in cluster 1. KDM1A and IGFBP5 were present in cluster 4; MMP19, G0S2, CH25H, and C1orf162 were upregulated mostly in clusters 0 and 2; and TSPAN3 and ITGB8 were extensively expressed in multiple clusters, though these genes were mainly detected in clusters 1 and 4. 3.9 Nomogram to qualitatively predict OS in patients with glioma With the CGGA dataset, six prognostic factors, including age, sex, tumor grade, recurrence status, IDH mutation status, and 1p19q codeletion status were combined with the risk score to construct a nomogram for predicting the 1-, 3-, and 5-year survival probabilities (Figure 12). 4 Discussion As one of the most lethal brain tumors, gliomas exhibit an extremely dreadful prognosis. Even with complete standard therapy, patients often quickly succumb to the disease owing to therapy resistance and tumor recurrence ( 23 ). Tumor heterogeneity, including genetic, epigenetic, and TME heterogeneity, contributes greatly to this doomed prognosis. Glioma heterogeneity is difficult to study with bulk genomic approaches but is better dissected using scRNA-seq (24). An in-depth understanding of the TME and related immune components via scRNA-seq analysis may facilitate the development and improvement of immunotherapy ( 25 ). In this study, starting from a glioma single-cell sequencing dataset, we divided single glioma cells into five state trajectories through dimensionality reduction and pseudo-temporal analysis. Then, we identified the differential genes between trajectories. Based on the TRGs, the patients in the CGGA dataset were separated into three clusters. The clinical features, survival prognosis, TME, immune infiltration, and ICG expression levels were significantly different among the clusters. Survival analysis identified immune cells and ICGs with clinical significance. Protein-protein interaction analysis and functional enrichment revealed the involvement of immune suppression pathways. Subsequently, clinically relevant TRGs were obtained from WGCNA analysis and intersected with glioma DEGs to identify candidate genes for establishing a prognosis prediction model. LASSO regression, multivariate regression, survival analysis, and ROC analysis were used to identify risk genes. Then, a prognosis prediction model was constructed and evaluated. This model has acceptable predictive efficacy. Finally, a nomogram was built to quantify the contribution of each risk factor in specific patients. Given the widespread success of immunomodulatory therapy across diverse cancer types, characterizing the immune fraction of the glioma TME and developing corresponding immunotherapies is of significant interest ( 26 ). The brain represents a relatively immunosuppressive microenvironment for glioma cells because it contains unique immune cells and regulatory cell types such as microglia and astrocytes, lacks classic lymphoid components, and is protected from circulating immune cells by the blood–brain barrier ( 27 ). Tumor-associated macrophages (TAMs) are essential glioma TME immune cells that can arise from two distinct sources: peripheral bone marrow-derived monocytes and yolk sac-derived microglia ( 28 , 29 ). TAMs are enriched in gliomas and can adopt multiple pro-tumorigenic phenotypes that play key roles in promoting invasion, angiogenesis, metastasis, and immunosuppression ( 30 ). Survival analysis of the immune cell signature revealed that the M0 macrophage signature significantly correlated with a worse prognosis. M0 macrophage infiltration, but not that of M1 or M2 macrophages, was increased in cluster C3. This result implies that macrophage M0 infiltration might partly contribute to the poor survival of patients with glioma. Macrophages are characterized by their plasticity and are commonly divided into classic-activated macrophages (M1) that stimulate anti-tumor immune responses and alternatively activated macrophages (M2), which suppress immune responses ( 34 ). M0 macrophages are non-polarized naïve macrophages. The M1 and M2 phases of TAMs are widely recognized in numerous cancers, while M2 macrophages have been reported mainly in high-grade glioma. However, a recent study reported the enrichment of nonpolarized M0 macrophages rather than M1 or M2 macrophage assemblies in glioblastoma ( 35 , 36 ). Our results are consistent with these findings and suggest that inducing macrophage polarization from M0 to M1 cells might have a therapeutic effect on gliomas. Further clustering of the CGGA cohort demonstrated that TRGs distinguish different clusters of patients with different OS and various TME compositions. Notably, cluster C3 showed the worst prognosis and had the lowest tumor purity and the highest immune score, indicating that immune cell infiltration in gliomas could be an adverse prognostic factor. This is consistent with their general immunosuppression functions. Thus, evaluating the TME composition using TRGs might have prognostic value for patients with glioma. Many other immune cells also play essential roles in glioma. Increased migration of dendritic cells (DCs) to nearby draining lymph nodes can lead to prolonged progression-free and improved OS in patients with glioblastoma ( 37 ). Neutrophils contribute to glioblastoma progression by supporting glioma stem cell pool expansion ( 38 ). CD8 + T cells function as anti-tumor components, and their activation with immune checkpoint inhibitors leads to longer survival in mouse glioma models ( 39 ). CD4 + lymphocyte enrichment, which could enhance the function of CD8 + T cells, profoundly decreases in glioma. Within the reduced number of CD4 + lymphocytes, CD4 + helper T (Th) cell numbers are reduced, while immunosuppressive Tregs increase ( 40 ). This could possibly explain why cluster C3 showed the highest fraction of anti-tumor CD8 + T cells but still displayed the worst prognosis. The survival advantage conferred by CD8 + T cells was likely offset by the immunosuppressive effects of Treg cells, which were significantly enriched in cluster C3. Immune checkpoints, regulators of the immune system, are crucial for self-tolerance and prevent the immune system from indiscriminately attacking cells. However, this mechanism could also be employed by cancer cells to protect themselves from immune attack ( 41 ). Among all the studied ICGs, CTLA4 and PD-L1 (CD274) are the most well-known. CTLA4 functions by competitive binding to CD80/86 ligands in DCs to prevent CD4 + helper T cell activation while boosting Treg, thereby leading to a pro-tumor immunosuppressive phenotype ( 42 ). PD-L1 and PD-L2 are expressed by tumor cells, and their interaction with PD-1 in immune cells can reduce immune responses ( 43 , 44 ). Numerous studies have been performed on other ICGs in glioma. A higher LDHA/LDHB ratio is associated with the IDH wild-type genotype and poorer OS in patients with glioma ( 45 ). B2M expression is correlated with malignancy of glioma ( 46 ). Inhibiting the ICOSLG axis in GBM may provide promising immunotherapeutic effects ( 47 ). LGALS9 regulates the immunosuppressive features of glioma ( 48 ). Depleting endogenous PVR from human glioma cells can inhibit tumor migration ( 49 ). HAVCR2 (TIM-3) and SIGLEC15 are well-described immunosuppressive components ( 50 , 51 ). Consistently, CD274, CTLA4, HAVCR2, PDCD1, PDCD1LG2, SIGLEC15, LDHA, B2M, ICOSLG, LGALS9, and PVR were up-regulated in cluster C3 samples, correlating with a worse prognosis. Likewise, LDHB, IL12, and LAMA3, which are anti-tumor ICGs and are associated with better survival ( 52 , 53 ), were downregulated in cluster C3. However, the roles of some ICGs are still not well defined. High PTPRC expression in glioma could result in both pro- and anti-immune responses ( 54 ). TNFRSF9 is upregulated in human gliomas compared to that in normal brain tissue ( 55 ), but its effect on gliomas has not been extensively studied. Further, TNFSF4 has not yet been studied in glioma. Our study sheds light on their possible roles in gliomas and uncovers novel potential therapeutic directions. Interestingly, CD8A, CD28, CD40, CD40LG, and CD80 are essential immune checkpoint components associated with T cell activation. These proteins are designated as anti-tumor immune markers in many extracranial solid tumors ( 51 , 56 ). However, in the present study, higher expression of these markers was associated with worse prognosis. Whether these proteins function differently in the unique glioma microenvironment or whether they are upregulated due to negative feedback from the specific immunosuppressive state of gliomas requires further study. Different clusters based on TRGs in the CGGA cohort showed distinct survival profiles, suggesting that TRG-based patient classification can be used to predict patient OS. Various TRG expression patterns indicate the diverse composition of the TME and immune components. Successfully predicting patient prognosis using TRG clustering indirectly shows that alterations in the TME and immunological components have a significant impact on the prognosis of patients with glioma. Thus, TRGs were used to identify risk genes and construct a prognosis prediction model. Our model had high accuracy and efficiency in both the training and testing datasets. To the best of our knowledge, this is the first TRG-based signature used to predict prognosis in patients with glioma. Moreover, a nomogram combining the TRG-based risk score and prognostic clinicopathological variables was constructed to provide a visual method for predicting patient OS. The developed nomogram was more accurate and effective than using the risk score alone. This study has several limitations. First, the correlation between glioma prognosis and immune components (immune cells and ICGs) was determined using bioinformatic analysis. The verification of conclusion and the in-depth exploration of internal mechanism depend on further molecular biology experiments. Second, the cell state trajectories and TRGs constructed by the scRNA-seq analysis were derived from the same GEO dataset and therefore lack systematic validation of repeated multi-sample experiments. In addition, the GEO, CGGA, and TCGA data used to construct the prognostic model were retrospective data. Large-scale prospective studies are necessary to assess the effectiveness and practicality of our models. 5 Conclusions We identified five state trajectories of single glioma cells with distinct cell identities based on scRNA-seq and demonstrated that TRG-based genotyping can accurately predict patient OS, clinicopathological features, TME composition, immune infiltration status, and ICG expression. Immune cells and ICGs that are correlated with prognosis of patients with glioma were identified and could be promising therapeutic targets. We developed a nomogram that integrates the TRG-based risk score signature and clinicopathological variables to provide a preliminary quantitative method for predicting patient OS. In summary, this study highlighted alterations in the TME and immune components in glioma, identified a TRG-based risk gene signature for predicting patient clinical outcomes, and explored associated immune responses to uncover potential immunotherapy directions for patients with glioma. Declarations Conflict of Interest The authors declare no conflicts of interest. Funding This study was supported by the Applied Basic Research joint project of Science and Technology Department of Yunnan Province and Kunming Medical University (No. 202364). Author Contribution Conceptualization, D.C. and J.F.; methodology, Z.F.; software, D.C. and Y.H.; validation, D.C. and L.D.; formal analysis, Y.H.; investigation, D.C. and Y.H.; resources, Y.H.; data curation, Y.H.; writing—original draft preparation, D.C. and Y.H.; writing—review and editing, D.C. and Y.H.; visualization, D.C.; supervision, D.C. and Y.H.; project administration, J.F.; and funding acquisition, J.F. All authors have read and agreed to the published version of the manuscript. Data Availability Statement Links to publicly archived datasets analyzed: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE84465 ; http://gepia2.cancer-pku.cn/#index ; https://xenabrowser.net/hub/ ; and http://www.cgga.org.cn/ , accessed on Jan 12, 2022. References Ostrom QT, Cote DJ, Ascha M, Kruchko C, Barnholtz-Sloan JS. Adult Glioma Incidence and Survival by Race or Ethnicity in the United States From 2000 to 2014. JAMA Oncol. 2018;4(9):1254–62. Ostrom QT, Gittleman H, Truitt G, Boscia A, Kruchko C, Barnholtz-Sloan JS. CBTRUS Statistical Report: Primary Brain and Other Central Nervous System Tumors Diagnosed in the United States in 2011–2015. Neuro Oncol. 2018;20(suppl_4):iv1-iv86. Wesseling P, Capper D. WHO 2016 Classification of gliomas. Neuropathol Appl Neurobiol. 2018;44(2):139–50. Qi Y, Liu B, Sun Q, Xiong X, Chen Q. Immune Checkpoint Targeted Therapy in Glioma: Status and Hopes. Front Immunol. 2020;11:578877. Yang K, Wu Z, Zhang H, Zhang N, Wu W, Wang Z, et al. Glioma targeted therapy: insight into future of molecular approaches. Mol Cancer. 2022;21(1):39. Louveau A, Harris TH, Kipnis J. Revisiting the Mechanisms of CNS Immune Privilege. Trends Immunol. 2015;36(10):569–77. Bailey SL, Carpentier PA, McMahon EJ, Begolka WS, Miller SD. Innate and adaptive immune responses of the central nervous system. Crit Rev Immunol. 2006;26(2):149–88. DeCordova S, Shastri A, Tsolaki AG, Yasmin H, Klein L, Singh SK, et al. Molecular Heterogeneity and Immunosuppressive Microenvironment in Glioblastoma. Front Immunol. 2020;11:1402. Weber JS, D'Angelo SP, Minor D, Hodi FS, Gutzmer R, Neyns B, et al. Nivolumab versus chemotherapy in patients with advanced melanoma who progressed after anti-CTLA-4 treatment (CheckMate 037): a randomised, controlled, open-label, phase 3 trial. Lancet Oncol. 2015;16(4):375–84. Motzer RJ, Escudier B, McDermott DF, George S, Hammers HJ, Srinivas S, et al. Nivolumab versus Everolimus in Advanced Renal-Cell Carcinoma. N Engl J Med. 2015;373(19):1803–13. Medikonda R, Dunn G, Rahman M, Fecci P, Lim M. A review of glioblastoma immunotherapy. J Neurooncol. 2021;151(1):41–53. Eberwine J, Sul JY, Bartfai T, Kim J. The promise of single-cell sequencing. Nat Methods. 2014;11(1):25–7. Gawad C, Koh W, Quake SR. Single-cell genome sequencing: current state of the science. Nat Rev Genet. 2016;17(3):175–88. Papalexi E, Satija R. Single-cell RNA sequencing to explore immune cell heterogeneity. Nat Rev Immunol. 2018;18(1):35–45. Lall S, Sinha D, Bandyopadhyay S, Sengupta D. Structure-Aware Principal Component Analysis for Single-Cell RNA-seq Data. J Comput Biol. 2018. Satija R, Farrell JA, Gennert D, Schier AF, Regev A. Spatial reconstruction of single-cell gene expression data. Nat Biotechnol. 2015;33(5):495–502. Wilkerson MD, Hayes DN. ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics. 2010;26(12):1572–3. Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells with CIBERSORT. Methods Mol Biol. 2018;1711:243–59. Xiang R, Rong Y, Ge Y, Song W, Ren J, Fu T. Cell differentiation trajectory predicts patient potential immunotherapy response and prognosis in gastric cancer. Aging (Albany NY). 2021;13(4):5928–45. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559. Tang Z, Kang B, Li C, Chen T, Zhang Z. GEPIA2: an enhanced web server for large-scale expression profiling and interactive analysis. Nucleic Acids Res. 2019;47(W1):W556-W60. Darmanis S, Sloan SA, Croote D, Mignardi M, Chernikova S, Samghababi P, et al. Single-Cell RNA-Seq Analysis of Infiltrating Neoplastic Cells at the Migrating Front of Human Glioblastoma. Cell Rep. 2017;21(5):1399–410. Ghotme KA, Barreto GE, Echeverria V, Gonzalez J, Bustos RH, Sanchez M, et al. Gliomas: New Perspectives in Diagnosis, Treatment and Prognosis. Curr Top Med Chem. 2017;17(12):1438–47. Nicholson JG, Fine HA. Diffuse Glioma Heterogeneity and Its Therapeutic Implications. Cancer Discov. 2021;11(3):575–90. Hernandez Martinez A, Madurga R, Garcia-Romero N, Ayuso-Sacido A. Unravelling glioblastoma heterogeneity by means of single-cell RNA sequencing. Cancer Lett. 2022;527:66–79. Waldman AD, Fritz JM, Lenardo MJ. A guide to cancer immunotherapy: from T cell basic science to clinical practice. Nat Rev Immunol. 2020;20(11):651–68. Quail DF, Joyce JA. The Microenvironmental Landscape of Brain Tumors. Cancer Cell. 2017;31(3):326–41. Chen Z, Feng X, Herting CJ, Garcia VA, Nie K, Pong WW, et al. Cellular and Molecular Identity of Tumor-Associated Macrophages in Glioblastoma. Cancer Res. 2017;77(9):2266–78. Bowman RL, Klemm F, Akkari L, Pyonteck SM, Sevenich L, Quail DF, et al. Macrophage Ontogeny Underlies Differences in Tumor-Specific Education in Brain Malignancies. Cell Rep. 2016;17(9):2445–59. Mantovani A, Marchesi F, Malesci A, Laghi L, Allavena P. Tumour-associated macrophages as treatment targets in oncology. Nat Rev Clin Oncol. 2017;14(7):399–416. Graeber MB, Scheithauer BW, Kreutzberg GW. Microglia in brain tumors. Glia. 2002;40(2):252–9. Shi C, Pamer EG. Monocyte recruitment during infection and inflammation. Nat Rev Immunol. 2011;11(11):762–74. Ginhoux F, Schultze JL, Murray PJ, Ochando J, Biswas SK. New insights into the multidimensional concept of macrophage ontogeny, activation and function. Nat Immunol. 2016;17(1):34–40. Murray PJ, Allen JE, Biswas SK, Fisher EA, Gilroy DW, Goerdt S, et al. Macrophage activation and polarization: nomenclature and experimental guidelines. Immunity. 2014;41(1):14–20. Gabrusiewicz K, Rodriguez B, Wei J, Hashimoto Y, Healy LM, Maiti SN, et al. Glioblastoma-infiltrated innate immune cells resemble M0 macrophage phenotype. JCI Insight. 2016;1(2). Huang L, Wang Z, Chang Y, Wang K, Kang X, Huang R, et al. EFEMP2 indicates assembly of M0 macrophage and more malignant phenotypes of glioma. Aging (Albany NY). 2020;12(9):8397–412. Mitchell DA, Batich KA, Gunn MD, Huang MN, Sanchez-Perez L, Nair SK, et al. Tetanus toxoid and CCL3 improve dendritic cell vaccines in mice and glioblastoma patients. Nature. 2015;519(7543):366–9. Liang J, Piao Y, Holmes L, Fuller GN, Henry V, Tiao N, et al. Neutrophils promote the malignant glioma phenotype through S100A4. Clin Cancer Res. 2014;20(1):187–98. Cohen JV, Kluger HM. Systemic Immunotherapy for the Treatment of Brain Metastases. Front Oncol. 2016;6:49. Fecci PE, Mitchell DA, Whitesides JF, Xie W, Friedman AH, Archer GE, et al. Increased regulatory T-cell fraction amidst a diminished CD4 compartment explains cellular immune defects in patients with malignant glioma. Cancer Res. 2006;66(6):3294–302. Pardoll DM. The blockade of immune checkpoints in cancer immunotherapy. Nat Rev Cancer. 2012;12(4):252–64. Wing K, Onishi Y, Prieto-Martin P, Yamaguchi T, Miyara M, Fehervari Z, et al. CTLA-4 control over Foxp3 + regulatory T cell function. Science. 2008;322(5899):271–5. Freeman GJ, Long AJ, Iwai Y, Bourque K, Chernova T, Nishimura H, et al. Engagement of the PD-1 immunoinhibitory receptor by a novel B7 family member leads to negative regulation of lymphocyte activation. J Exp Med. 2000;192(7):1027–34. Latchman Y, Wood CR, Chernova T, Chaudhary D, Borde M, Chernova I, et al. PD-L2 is a second ligand for PD-1 and inhibits T cell activation. Nat Immunol. 2001;2(3):261–8. Givechian KB, Garner C, Benz S, Rabizadeh S, Soon-Shiong P. Glycolytic expression in lower-grade glioma reveals an epigenetic association between IDH mutation status and PDL1/2 expression. Neurooncol Adv. 2021;3(1):vdaa162. Zhang H, Cui B, Zhou Y, Wang X, Wu W, Wang Z, et al. B2M overexpression correlates with malignancy and immune signatures in human gliomas. Sci Rep. 2021;11(1):5045. Iwata R, Hyoung Lee J, Hayashi M, Dianzani U, Ofune K, Maruyama M, et al. ICOSLG-mediated regulatory T-cell expansion and IL-10 production promote progression of glioblastoma. Neuro Oncol. 2020;22(3):333–44. Liang T, Wang X, Wang F, Feng E, You G. Galectin-9: A Predictive Biomarker Negatively Regulating Immune Response in Glioma Patients. World Neurosurg. 2019;132:e455-e62. Sloan KE, Stewart JK, Treloar AF, Matthews RT, Jay DG. CD155/PVR enhances glioma cell dispersal by regulating adhesion signaling and focal adhesion dynamics. Cancer Res. 2005;65(23):10930–7. Huang YH, Zhu C, Kondo Y, Anderson AC, Gandhi A, Russell A, et al. CEACAM1 regulates TIM-3-mediated tolerance and exhaustion. Nature. 2015;517(7534):386–90. Morad G, Helmink BA, Sharma P, Wargo JA. Hallmarks of response, resistance, and toxicity to immune checkpoint blockade. Cell. 2021;184(21):5309–37. Yan J, Smyth MJ, Teng MWL. Interleukin (IL)-12 and IL-23 and Their Conflicting Roles in Cancer. Cold Spring Harb Perspect Biol. 2018;10(7). Tang L, Wang P, Wang Q, Zhong L. Correlation of LAMA3 with onset and prognosis of ovarian cancer. Oncol Lett. 2019;18(3):2813–8. Brandenburg S, Turkowski K, Mueller A, Radev YT, Seidlitz S, Vajkoczy P. Myeloid cells expressing high level of CD45 are associated with a distinct activated phenotype in glioma. Immunol Res. 2017;65(3):757–68. Blank AE, Baumgarten P, Zeiner P, Zachskorn C, Loffler C, Schittenhelm J, et al. Tumour necrosis factor receptor superfamily member 9 (TNFRSF9) is up-regulated in reactive astrocytes in human gliomas. Neuropathol Appl Neurobiol. 2015;41(2):e56-67. Zhang Y, Huang T, Hu Y, Wang Y. Activation of CD40 by soluble recombinant human CD40 ligand inhibits human glioma cells proliferation via nuclear factor-kappaB signaling pathway. J Huazhong Univ Sci Technolog Med Sci. 2012;32(5):691–6. Additional Declarations No competing interests reported. Supplementary Files FigS1.tif FigS2.tif FigS3.tif Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-3829624","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":265133058,"identity":"051203ba-f261-44de-98c4-fa3837f8306e","order_by":0,"name":"Danlei Chen","email":"","orcid":"","institution":"First People’s Hospital of Yunnan Province","correspondingAuthor":false,"prefix":"","firstName":"Danlei","middleName":"","lastName":"Chen","suffix":""},{"id":265133059,"identity":"417be398-ba26-40c0-a662-ccdc2afd10ce","order_by":1,"name":"Yi He","email":"","orcid":"","institution":"First Affiliated Hospital of Fujian Medical University","correspondingAuthor":false,"prefix":"","firstName":"Yi","middleName":"","lastName":"He","suffix":""},{"id":265133060,"identity":"1feda91f-b20b-45c3-9617-673609a0d70d","order_by":2,"name":"Zhiyuan Feng","email":"","orcid":"","institution":"First People’s Hospital of Yunnan Province","correspondingAuthor":false,"prefix":"","firstName":"Zhiyuan","middleName":"","lastName":"Feng","suffix":""},{"id":265133061,"identity":"e9a3a2df-f1b3-49ec-890a-426777f84e78","order_by":3,"name":"Longsheng Dong","email":"","orcid":"","institution":"First People’s Hospital of Yunnan Province","correspondingAuthor":false,"prefix":"","firstName":"Longsheng","middleName":"","lastName":"Dong","suffix":""},{"id":265133062,"identity":"a86ba249-cd37-4636-9c91-7128345933b9","order_by":4,"name":"Junfeng Yang","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAtElEQVRIiWNgGAWjYFCCg49BpBwbe/sBYrUcNgaRxnw8ZxKI1cIM1pI4T8LBgDgNBgcPMxvdzLFLb5NgSGD4UbGNCC0HDjMn525Lzm2TbjzA2HPmNjFazh8+nLvtQG6bzIEEZsY2orQcZgZpSWeTSDAgXgvQYQcSiNciCdRiDPSLYRswkA8S5Re+G4eZpXO32cnLt7cffPCjgggtCjcOIDgHcKlCAfL9DUSpGwWjYBSMgpEMAJwKQJ4tkDTGAAAAAElFTkSuQmCC","orcid":"","institution":"First People’s Hospital of Yunnan Province","correspondingAuthor":true,"prefix":"","firstName":"Junfeng","middleName":"","lastName":"Yang","suffix":""}],"badges":[],"createdAt":"2024-01-02 13:44:12","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-3829624/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-3829624/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":49237620,"identity":"08a8c296-4237-4804-a29e-c8ae16d21e7a","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":477658,"visible":true,"origin":"","legend":"\u003cp\u003eThe study workflow.\u003c/p\u003e","description":"","filename":"Fig1.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/adf884c7907ee99547c6d4c2.png"},{"id":49237621,"identity":"182a027d-39e0-48a9-af05-93aac64859c2","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":665962,"visible":true,"origin":"","legend":"\u003cp\u003eDimensionality reduction and cell trajectory analysis identified cell clusters, types, states, and pesudotime correlation.\u003cstrong\u003e (A) \u003c/strong\u003eA scatter plot displaying the expression level and intercellular variation of each gene, with the top 10 variable genes labeled in the plot.\u003cstrong\u003e (B) \u003c/strong\u003eThe top 10 marker genes in each cluster are displayed on the heat map. \u003cstrong\u003e(C) \u003c/strong\u003eNine clusters of glioma single cells were segregated by the top 15 PC. \u003cstrong\u003e(D-F)\u003c/strong\u003ePesudotime and trajectory analyses identified cell clusters, states, and pesudotime correlation of glioma single cells.\u003c/p\u003e","description":"","filename":"Fig2.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/cd97e07dfff7138b97e6e61f.png"},{"id":49237956,"identity":"74a8282f-f9cf-4223-8b80-437b6475c334","added_by":"auto","created_at":"2024-01-05 18:09:36","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":463395,"visible":true,"origin":"","legend":"\u003cp\u003eGO and KEGG enrichment for TRGs of state 1, 3, or all combined states. TRGs: trajectory-related genes, GO: gene ontology, BP: biological process, MF: molecular function, CC: cellular component, KEGG: Kyoto Encyclopedia of Genes and Genomes.\u003c/p\u003e","description":"","filename":"Fig3.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/8cfe2b714c8c1a8bf0c1f64e.png"},{"id":49237624,"identity":"cfab5369-19d7-4854-935c-ce5d4a4dc8ce","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":223207,"visible":true,"origin":"","legend":"\u003cp\u003eTRG-based clustering of the CGGA cohort segregates patients with glioma with different prognosis and clinical features. \u003cstrong\u003e(A)\u003c/strong\u003e The CGGA cohort was divided into 3 clusters (C1–3) according to TRG expression. \u003cstrong\u003e(B)\u003c/strong\u003e Survival analysis showing patients in cluster C3 have significantly worse prognosis compared with patients in the other two clusters. \u003cstrong\u003e(C)\u003c/strong\u003e Clinical correlation analysis showing that patients in cluster C3 have a higher proportion of clinical features associated with worse prognosis (namely older age, higher tumor grade, increased recurrence rate, IDH wild type, and 1p19q non-codeletion).\u003c/p\u003e","description":"","filename":"Fig4.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/53b5756c9fca245b3d1c30f3.png"},{"id":49237622,"identity":"5363e8e1-9b18-46bc-93ce-1a0e65870cb0","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":654303,"visible":true,"origin":"","legend":"\u003cp\u003eTRG-based clustering isolates subgroups of glioma samples with different microenvironment composition and altered immune infiltration status. \u003cstrong\u003e(A)\u003c/strong\u003e Glioma samples in different clusters show various stromal scores, immune scores, and tumor purity. \u003cstrong\u003e(B)\u003c/strong\u003eClusters segregated by TRGs exhibit diverse immune cell infiltration. \u003cstrong\u003e(C)\u003c/strong\u003eThe proportions of immune infiltration in three clusters are shown in the bar chart. \u003cstrong\u003e(D–G)\u003c/strong\u003e Kaplan–Meier plots of survival related immune cells. *P \u0026lt; 0.05; **P \u0026lt; 0.01; ***P \u0026lt; 0.001.\u003c/p\u003e","description":"","filename":"Fig5.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/00a7868958c477f7d40c6d36.png"},{"id":49239825,"identity":"fd9e0ab3-17f7-4240-b6d9-8233f0e536d2","added_by":"auto","created_at":"2024-01-05 18:17:36","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":595666,"visible":true,"origin":"","legend":"\u003cp\u003eImmune checkpoint genes are differently expressed among TRG-based glioma clusters and are correlated with patient survival status. \u003cstrong\u003e(A)\u003c/strong\u003e Differential analysis showing glioma samples in three clusters have different ICG expression. \u003cstrong\u003e(B)\u003c/strong\u003eKaplan–Meier plots of prognosis-related ICGs. *P \u0026lt; 0.05; **P \u0026lt; 0.01; ***P \u0026lt; 0.001.\u003c/p\u003e","description":"","filename":"Fig6.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/acb1abe05474655f65a35a5a.png"},{"id":49237961,"identity":"7af1b89b-ca35-44fa-b7c7-5fd128062558","added_by":"auto","created_at":"2024-01-05 18:09:36","extension":"png","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":812462,"visible":true,"origin":"","legend":"\u003cp\u003eProtein-protein interaction analysis identifies core ICG elements and reveals the involvement of immune suppression pathways. \u003cstrong\u003e(A)\u003c/strong\u003e STRING enrichment results of prognosis-related ICGs. \u003cstrong\u003e(B)\u003c/strong\u003eFunctional enrichment analysis with STRING indicating that prognosis-related ICGs are closely associated with immune suppression pathways. \u003cstrong\u003e(C)\u003c/strong\u003ePrognosis-related ICGs were further processed with MCODE module in Cytoscape to identify core components.\u003c/p\u003e","description":"","filename":"Fig7.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/1f5c037e33dee43d2654c8be.png"},{"id":49237958,"identity":"212fa668-196f-4dd6-8d3e-7b56bf4ed22c","added_by":"auto","created_at":"2024-01-05 18:09:36","extension":"png","order_by":8,"title":"Figure 8","display":"","copyAsset":false,"role":"figure","size":588305,"visible":true,"origin":"","legend":"\u003cp\u003eWGCNA and differential expression analysis identifies candidate genes for building the prognosis prediction model. \u003cstrong\u003e(A)\u003c/strong\u003e Dendrogram displaying the genes clustered in different modules. \u003cstrong\u003e(B)\u003c/strong\u003e Module traits showing the correlations between modules and clinical features. The upper numbers indicate the correlation coefficient and the lower number indicates the statistical significance (P \u0026lt; 0.05 is considered significant). \u003cstrong\u003e(C–D) \u003c/strong\u003eThe chromosomal distribution of differentially expressed genes in bulk glioma samples compared with normal tissue. \u003cstrong\u003e(C)\u003c/strong\u003e Low grade glioma, LGG; \u003cstrong\u003e(D)\u003c/strong\u003e Glioblastoma, GBM. \u003cstrong\u003e(E)\u003c/strong\u003eVenn plot showing the intersection of clinically-related module genes from WGCNA and DEGs from GEPIA2.\u003c/p\u003e","description":"","filename":"Fig8.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/e692ebdefbeba1e8afd03e14.png"},{"id":49239826,"identity":"f732c183-895a-4551-8bb6-70950b569b2c","added_by":"auto","created_at":"2024-01-05 18:17:36","extension":"png","order_by":9,"title":"Figure 9","display":"","copyAsset":false,"role":"figure","size":111924,"visible":true,"origin":"","legend":"\u003cp\u003eMultivariate analysis of candidate genes which passed LASSO regression and collinearity tests. *P \u0026lt; 0.05; **P \u0026lt; 0.01; ***P \u0026lt; 0.001.\u003c/p\u003e","description":"","filename":"Fig9.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/322712ce122fb2d6d3bc972b.png"},{"id":49237632,"identity":"47187c11-e9b2-410e-a293-51571ffcda37","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":10,"title":"Figure 10","display":"","copyAsset":false,"role":"figure","size":137208,"visible":true,"origin":"","legend":"\u003cp\u003eThe risk scoring signature has good prognostic value. Patients in the TGCA training cohort \u003cstrong\u003e(A, C, E, G)\u003c/strong\u003e and CGGA testing cohort \u003cstrong\u003e(B, D, F, H)\u003c/strong\u003e were divided into high- and low- risk groups based on the risk score. Kaplan–Meier plots were used to analyze the survival difference between the high- and low- risk groups \u003cstrong\u003e(A, B)\u003c/strong\u003e. ROC analysis was used to examine the sensitivity and specificity of the established model \u003cstrong\u003e(C, D)\u003c/strong\u003e. Uni- \u003cstrong\u003e(E, F)\u003c/strong\u003e and multi- \u003cstrong\u003e(G, H)\u003c/strong\u003e variate analyses were conducted to confirm the risk signature on patient prognosis and the independence of the model.\u003c/p\u003e","description":"","filename":"Fig10.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/040849879a7853469dd6a0f7.png"},{"id":49237631,"identity":"98268529-22cf-4370-9bcc-4607f96060c6","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":11,"title":"Figure 11","display":"","copyAsset":false,"role":"figure","size":453273,"visible":true,"origin":"","legend":"\u003cp\u003eExpression pattern of risk genes in single cell clusters. \u003cstrong\u003e(A)\u003c/strong\u003e Scatter plots showing the expression pattern of 12 risk genes in the 9 single glioma cell clusters. \u003cstrong\u003e(B)\u003c/strong\u003e Bubble plot showing the expression levels of 12 risk genes in each single cell cluster.\u003c/p\u003e","description":"","filename":"Fig11.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/a2f0408c16d435bd86582d64.png"},{"id":49237627,"identity":"9b63f0a7-0759-42c6-a48b-4111263b9601","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"png","order_by":12,"title":"Figure 12","display":"","copyAsset":false,"role":"figure","size":67138,"visible":true,"origin":"","legend":"\u003cp\u003eA nomogram for predicting the 1-, 3-, and 5-year survival possibility of patients with glioma. The points for each item are calculated according to the corresponding status and the matching location in the first row of points. The sum of points is then used to predict the 1-, 3-, and 5-year survival possibility.\u003c/p\u003e","description":"","filename":"Fig12.png","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/d5e820612c7c41cd683b7f42.png"},{"id":49779240,"identity":"26fb0cf0-a962-4167-94a0-ca02b32c879f","added_by":"auto","created_at":"2024-01-17 23:07:18","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":5293620,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/30890959-5524-4340-84c9-b633011aff30.pdf"},{"id":49237959,"identity":"dd31b678-3a13-495f-a394-74981c6fbc9b","added_by":"auto","created_at":"2024-01-05 18:09:36","extension":"tif","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":1512660,"visible":true,"origin":"","legend":"","description":"","filename":"FigS1.tif","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/c19e14ce6e6a87b9ef380c92.tif"},{"id":49237635,"identity":"ac31d5ad-1bf3-4f6d-9614-45f1be83accc","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"tif","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":1930776,"visible":true,"origin":"","legend":"","description":"","filename":"FigS2.tif","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/90cabb89b8ebe85ca11e0c9c.tif"},{"id":49237633,"identity":"99ddaa51-ee0a-41a1-8e04-c62bd9fe35db","added_by":"auto","created_at":"2024-01-05 18:01:36","extension":"tif","order_by":3,"title":"","display":"","copyAsset":false,"role":"supplement","size":6969712,"visible":true,"origin":"","legend":"","description":"","filename":"FigS3.tif","url":"https://assets-eu.researchsquare.com/files/rs-3829624/v1/91fac3b7044580c424868f55.tif"}],"financialInterests":"No competing interests reported.","formattedTitle":"Tumor microenvironment variation implicates immune alterations and correlates with prognosis in patients with glioma","fulltext":[{"header":"1 Introduction","content":"\u003cp\u003eAs the most common primary central nervous system neoplasms, gliomas encompass a large group of heterogeneous, aggressive tumors. Even if clinicians conduct current standard treatments, including a combination of surgery, radiotherapy, and chemotherapy, most tumors will recur due to chemo- and radioresistance. Patients usually succumb quickly to the disease, with approximately 5% five-year-survival rates (\u003cspan additionalcitationids=\"CR2\" citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e). Recent advances in genomics, immunology, and other disciplines have resulted in various experimental therapies, including targeted therapy and immunotherapy, with the promise of providing new therapeutic strategies for gliomas. However, the efficacy of these new therapies remains unsatisfactory (\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e, \u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e). Therefore, it is important to further explore the internal mechanisms of gliomas to uncover new therapeutic targets and approaches.\u003c/p\u003e \u003cp\u003eThe brain is considered an immune-privileged site because it lacks a traditional lymphatic system and contains only a few antigen-presenting cells. Recent findings changed this concept, demonstrating that the brain can coordinate a robust immune response involving both the innate and adaptive immune systems (\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e, \u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e). The glioma-specific immunosuppressive microenvironment contributes greatly to the poor prognosis of patients with glioma (\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e). Novel immunotherapies targeting the tumor microenvironment (TME) and its immune components have shown great promise for the clinical management of many extracranial solid tumors, including melanoma and renal cell carcinoma (\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e, \u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e). However, preliminary clinical trials using checkpoint inhibitors and vaccine therapy for glioma have been largely disappointing. The failure was likely due to the highly immunosuppressive environment, systemic immunosuppression, local immune dysfunction, and high tumor heterogeneity (\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e, \u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e). Thus, the TME and immune alterations are of great interest for immunotherapy and glioma management.\u003c/p\u003e \u003cp\u003eSingle-cell sequencing is capable of examining sequence information from individual cells, probing cellular differences with high resolution, and offering a better understanding of how individual cells function in their microenvironments (\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e). Single-cell resolution genomic information can help uncover the roles of genetic mosaicism and intratumoral genetic heterogeneity during tumor development or treatment responses (\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e). Given the heterogeneity of glioma, single-cell profiling of genomic alterations in glioma could help dissect TMEs and provide important insights into the biology of tumor-infiltrating immune cells, including cell heterogeneity, dynamics, and potential roles in disease progression and responses to immunotherapy (\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eHere, we used single-cell RNA-sequencing (scRNA-seq) data from glioma samples to identify five glioma cell trajectories with distinct identities. Then, we used the trajectory-related genes (TRGs) for immune and survival analysis. Figure\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e shows the study workflow, where glioma cell clustering, cell state trajectories, and TRGs were first acquired from the glioma scRNA-seq dataset in the GEO database. TRGs were then integrated into the Chinese Glioma Genome Atlas (CGGA) dataset to perform patient clustering, survival analysis, and immune-related analysis. Additionally, clinically-relevant TRGs were acquired through weighted correlation network analysis (WGCNA) and utilized to establish a prognostic model. A prognostic risk scoring signature was obtained through Least absolute shrinkage and selection operator (LASSO) and multivariate Cox regression analyses. Then, the scores were evaluated using ROC and Kaplan\u0026ndash;Meier (K\u0026ndash;M) curve analyses and quantified using a nomogram.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e"},{"header":"2 Materials and Methods","content":"\u003ch2\u003e2.1\u0026nbsp; \u0026nbsp; \u0026nbsp;Acquisition and processing of scRNA-seq and bulk RNA-seq data\u003c/h2\u003e\n\u003cp\u003escRNA-seq data from 3589 glioma cells were obtained from the Gene Expression Omnibus (GEO; accession number GSE84465; https://www.ncbi.nlm.nih. gov/geo/) database. The R package ‘Seurat’ was used to process the scRNA-seq data and perform quality control. Cells with a gene number \u0026lt; 100, read number \u0026lt; 50, and mitochondrial gene content \u0026gt; 5% were excluded from the analysis. The LogNormalize method was then used to normalize the data, and the top 1500 genes with highly variable features were identified by variance analysis. Bulk RNA-seq data from 702 and 693 glioma samples were downloaded from the TCGA (http://cancergenome.nih.gov/) and CGGA mRNAseq-693 datasets (http://www.cgga.org.cn/index.jsp). For differential expression analysis of glioma compared to normal tissue, mRNA expression files from 207 normal, 518 low-grade glioma, and 163 glioblastoma samples were downloaded from the GEPIA2 database (http://gepia2.cancer-pku.cn/#index). Differentially expressed genes were identified using |log\u003csub\u003e2\u003c/sub\u003e (fold change)| \u0026gt; 1 and adjusted P-value \u0026lt; 0.05.\u003c/p\u003e\n\u003ch2\u003e2.2\u0026nbsp; \u0026nbsp; \u0026nbsp;Dimensionality reduction and cell annotation\u003c/h2\u003e\n\u003cp\u003eA false discovery rate (FDR) \u0026lt; 0.05 cutoff was used to identify the top principal components (PCs), which were then used for dimension reduction with the tSNE algorithm to obtain the major single-cell clusters\u0026nbsp;(15, 16). Based on |log\u003csub\u003e2\u003c/sub\u003e (fold change)| \u0026gt; 0.5 and FDR \u0026lt; 0.05, marker genes in each cluster were identified, and the top marker genes from the clusters were visualized using a heatmap. The clusters were annotated using the ‘SingleR’ package.\u0026nbsp;\u003c/p\u003e\n\u003ch2\u003e2.3\u0026nbsp; \u0026nbsp; \u0026nbsp;Pseudotime and trajectory analyses\u003c/h2\u003e\n\u003cp\u003ePseudotime and trajectory analyses of single glioma cells were performed in R with the ‘Monocle’ package. Differentially expressed genes between trajectories of different states were obtained by comparing each trajectory with all other trajectories as a whole (using |log2(fold change)| \u0026gt; 0.5 and FDR \u0026lt; 0.05 thresholds) and designated as TRGs. The ‘org.Hs.eg.db,’ ‘clusterProfiler,’ ‘enrichplot,’ and ‘ggplot2’ packages were used for TRG gene ontology (GO) annotation and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis.\u003c/p\u003e\n\u003ch2\u003e2.4\u0026nbsp; \u0026nbsp; \u0026nbsp;TRG-based molecular clustering of patients with glioma from the CGGA dataset\u003c/h2\u003e\n\u003cp\u003eThe ‘ConsensusClusterPlus’ package is an R package capable of providing quantitative and visual evidence of stability, which is used to measure the number of unsupervised subtypes\u0026nbsp;(17). Based on the TRGs, ‘ConsensusClusterPlus’ divided the CGGA glioma cohort into different clusters. To determine the optimal cluster number, the K-means algorithm and cumulative distribution function (CDF) curve were utilized with 50 iterations and a maximum tentative clustering number of 9. The survival difference was examined with K–M analysis using the ‘survival’ and ‘survminer’ packages. The clinicopathological features in each molecular cluster were visualized using the ‘ggplot2’ package.\u0026nbsp;\u003c/p\u003e\n\u003ch2\u003e2.5\u0026nbsp; \u0026nbsp; \u0026nbsp;TME, immune cell infiltration, and immune checkpoint analysis\u003c/h2\u003e\n\u003cp\u003eFor each sample, TME scores, including immune score, stromal score, and tumor purity, were calculated using the ‘ESTIMATE’ package. The infiltration of 22 immune cells in each sample was identified using CIBERSORT software\u0026nbsp;(18). Meanwhile, the immune cell infiltration density in different glioma clusters was compared using the ‘limma’ package and is shown as histograms and rainbow plots. Then, 38 immune checkpoint genes (ICGs) (acquired from a previous study\u0026nbsp;(19)) across glioma clusters were examined for differential expression in the CGGA dataset. K–M analysis was used to evaluate the prognostic value of immune cells and ICGs. The prognosis-related ICGs were then analyzed using STRING enrichment to explore their underlying interactions and mechanisms. The core genes of the interaction network were identified using the MCODE application in Cytoscape software. The most significantly enriched pathways are shown in a bubble plot generated using the ggplot2 package.\u003c/p\u003e\n\u003ch2\u003e2.6\u0026nbsp; \u0026nbsp; \u0026nbsp;WGCNA and establishing prognostic risk scoring signatures\u003c/h2\u003e\n\u003cp\u003eWGCNA is a systematic method for identifying modules of highly correlated genes, summarizing such modules using eigengene modules, and relating modules to one another and to external sample traits\u0026nbsp;(20). In this study, WGCNA was used to analyze TRGs and identify key modules related to clinical features, especially patient survival. DEGs between glioma and normal tissues were obtained from the GEPIA2 online platform\u0026nbsp;(21). The TRGs in the survival-related WGCNA modules were intersected with the DEGs from GEPIA2 to obtain candidate genes, which were further used in univariate, LASSO, and multivariate Cox regression analyses to generate prognostic risk-scoring signatures. LASSO regression identified 24 prognosis-related genes, among which 17 genes passed the collinearity test. These 17 genes were used for multivariate analysis to obtain the final 12 risk genes. The coefficients of the risk genes were determined with the training dataset, and the total risk score was calculated as the sum of the products of prognostic gene expression multiplied by their respective coefficients. The mean value of the total risk score in the training dataset was used as the threshold for patient risk stratification (high-/low-risk). K–M analysis and ROC curves were consecutively used to evaluate the accuracy and prediction efficacy of the prognostic risk-scoring signature. Multivariate Cox regression analysis was used to examine whether the risk-scoring signature was independent of other clinical features.\u003c/p\u003e\n\u003ch2\u003e2.7\u0026nbsp; \u0026nbsp; \u0026nbsp;Nomogram construction\u003c/h2\u003e\n\u003cp\u003eUsing the CGGA dataset, clinical features (e.g., age, sex, tumor grade, recurrent status, IDH mutation status, and 1p19q codeletion status) and risk scores were combined into a nomogram with the ‘rms’ package to predict the 1-, 3-, and 5-year overall survival (OS) for each patient with glioma.\u003c/p\u003e"},{"header":"3 Results","content":"\u003ch2\u003e3.1\u0026nbsp; \u0026nbsp; \u0026nbsp;scRNA-seq quality control and normalization\u003c/h2\u003e\n\u003cp\u003eIn this study, 3589 glioma cells were obtained from the GEO database (accession number GSE84465). After quality control and normalization, the samples showed good RNA features and count distribution. No correlation was detected between sequencing depth and mitochondrial gene sequences (Supplementary Figure 1A-C). We observed a significant positive correlation between sequencing depth and total intracellular sequences (R = 0.29, Supplementary Figure 1D). A total of 19312 genes were analyzed, of which the top 10 genes with the highest variation between cells were labeled (Figure 2A).\u003c/p\u003e\n\u003ch2\u003e3.2\u0026nbsp; \u0026nbsp; \u0026nbsp;Dimension reduction and cell trajectory analysis identified 5 glioma cell states\u003c/h2\u003e\n\u003cp\u003eGiven the limited computing power, principal component analysis (PCA) was first used for preliminary dimension reduction of scRNA-seq data and showed good segregation among glioma cells (Figure S2a). All top 15 principal components (PCs) showed significant differences and were enrolled for further analysis (Figure S2b). Top 20 variable genes for PC1-4 were displayed with scatter plots and heatmaps in Figure S2c, d. Then, utilizing t-distributed Stochastic Neighbor Embedding (tSNE) algorithm, glioma cells were aggregated into 9 clusters based on the top 15 PCs, and a total of 1061 marker genes were identified by differential analysis, and the top 10 marker genes in each cluster were demonstrated with heatmap (Figure 2b-c).\u003c/p\u003e\n\u003cp\u003eAfterwards, pseudotime trajectory analysis revealed that clusters 3, 4, and 5 were distributed in state 1 and 5, which located in the relatively earlier stages of differentiation. In contrast, the clusters 0 and 2 were enriched in state 3, the latter stage of differentiation. (Figure 2d-f)\u003c/p\u003e\n\u003ch2\u003e3.3\u0026nbsp; \u0026nbsp; \u0026nbsp;Molecular function enrichment of trajectory-related genes\u003c/h2\u003e\n\u003cp\u003eBy comparing each pseudotime state trajectory with the remaining trajectories combined, we identified 601 glioma TRGs, with 313 TRGs for state 1 and 391 TRGs for state 3 (TRGs of different trajectories might overlap with each other). GO and KEGG enrichment analyses of TRGs of all states and showed significant enrichment in many immune-related pathways, including T cell activation, leukocyte-mediated immunity, and MHC protein complexes. Given that states 1 and 3 represented the early and late differentiation states respectively, we further ran enrichment analysis for TGRs of them, which showed results similar to total TRGs (Figure 3).\u003c/p\u003e\n\u003ch2\u003e3.4\u0026nbsp; \u0026nbsp; \u0026nbsp;Genotyping based on TRGs segregated CGGA patients with glioma into three clusters with different prognoses\u003c/h2\u003e\n\u003cp\u003eTRG-based consensus clustering divided CGGA patients with glioma into three clusters (cluster C1–C3) with good differentiation and sufficient samples in each cluster (Figure 4A). K–M analysis was then used to determine the clinical significance of the consensus clustering. This analysis showed that cluster C3 had the worst OS (Figure 4B, P \u0026lt; 0.001). Clinical relevance studies further revealed that patients in cluster C3 had significantly older age, higher glioma grades, increased tumor recurrence rate, and a higher proportion of IDH wildtype and 1p19q non-co-deletion genotypes. Each of these parameters was associated with a worse prognosis (Figure 4C).\u003c/p\u003e\n\u003ch2\u003e3.5\u0026nbsp; \u0026nbsp; \u0026nbsp;Clustering based on TRGs segregated glioma samples with different tumor microenvironment compositions and immune cell infiltration\u003c/h2\u003e\n\u003cp\u003eTME analysis indicated that samples in Cluster C3 had the highest stromal, immune, ESTIMATE scores, and the lowest tumor purity (Figure 5A). Further analysis of immune cell infiltration showed that cluster C3 had higher infiltration levels of many immune cells, including CD8+ T cells, regulatory T cells (Treg), and M0 macrophages (Figure 5B, C). Subsequent survival analysis of immune cell infiltration showed that 4 of 22 immune cells were prognosis-related, including monocytes, resting dendritic cells, memory B cells, and M0 macrophages. Notably, M0 macrophage infiltration was considerably increased in cluster C3, which displayed a worse patient prognosis. In addition, cluster C3 showed lower monocyte infiltration, which was correlated with better survival (Figure 5D–G). This result implies that the worse prognosis of patients with glioma in cluster C3 might be partly due to enriched M0 macrophages.\u003c/p\u003e\n\u003ch2\u003e3.6\u0026nbsp; \u0026nbsp; \u0026nbsp;ICGs correlated with the prognosis of patients with glioma and might contribute to the survival difference between clusters segregated by TRGs\u003c/h2\u003e\n\u003cp\u003eDifferential expression analysis showed that 36 of 38 ICGs were significantly differentially expressed between the TRG-based clusters (Figure 6A). In addition, survival analysis indicated that 25 of 38 ICGs were prognosis-related. ICGs that were correlated with longer survival were less expressed in cluster C3 (IL12A, LAMA3, LDHB, and VTCN), whereas ICGs that were associated with worse prognosis were upregulated in cluster C3 (B2M, CD8A, CD28, CD40, CD40LG, CD80, CD86, CD274, CTLA4, HAVCR2, ICOS, ICOSLG, LDHA, LGALS9, PDCD1, PDCD1LG2, PTPRC, PVR, SIGLEG15, TNFRSF9, and TNFSF4). K–M plots of ICGs with P \u0026lt; 0.01 (except VTCN, P = 0.012) are displayed in Figure 6B. These prognosis-related ICGs were then analyzed using the STRING online platform. Prognosis-related ICGs were significantly enriched in pathways associated with immune tolerance and immunodeficiency (Figure 7A, B). The protein-protein interaction networks obtained from the STRING analysis were further processed with MCODE to identify the core ICGs, which were used to construct a core network of 17 ICGs (Figure 7C). These results suggest possible underlying mechanisms for survival differences between patient clusters and directions for further immunotherapy in gliomas.\u003c/p\u003e\n\u003ch2\u003e3.7\u0026nbsp; \u0026nbsp; \u0026nbsp;Identification of candidate prognostic genes\u003c/h2\u003e\n\u003cp\u003eWGCNA showed that TRGs were divided into nine modules, six of which (MEblue, MEturquoise, MEgreen, MEpink, MEbrown, and MEyellow; 461 genes in total) showed a significant correlation with both survival time and survival status (P \u0026lt; 0.05) (Figure 8A, B). Differential gene analysis using the TCGA glioma dataset identified 5745 and 7659 differentially expressed genes (DEGs) in low-grade glioma and glioblastoma, respectively, compared with normal tissue (Figure 8C, D). The union of those DEGs (8881 genes in total) were intersected with TRGs in the prognosis-related modules. Together, 381 intersecting genes were designated as candidate genes to construct the prognostic model (Figure 8E).\u003c/p\u003e\n\u003ch2\u003e3.8\u0026nbsp; \u0026nbsp; \u0026nbsp;Establishment and evaluation of prognostic risk scoring signature\u003c/h2\u003e\n\u003cp\u003eThe candidate genes were tested using univariate analysis, which showed that 360 of 381 candidate genes were associated with prognosis (Supplementary Figure 3). Given the numerous genes obtained by univariate analysis, we used LASSO regression analysis to filter the genes for downstream multivariate analysis. LASSO regression showed that 24 genes were correlated with survival. After further multivariate regression analysis, 12 genes were identified as independent risk genes (Figure 9) and were used to establish the prognostic risk scoring signature.\u003c/p\u003e\n\u003cp\u003eBased on the score signature and corresponding gene expression in each sample, the TCGA-Low grade glioma (LGG) and glioblastoma (GBM) datasets were used to train the prognosis prediction model and obtain the coefficients for each risk gene. Then, the TCGA training cohort and CGGA testing cohort were divided into low- or high-risk groups according to the total risk score. Survival analysis showed that patients in the high-risk group had a significantly shorter OS (Figure 10A, B). ROC analysis indicated that the model had excellent sensitivity and specificity for predicting the 1-, 3-, and 5-year survival rate in both the TCGA and CGGA cohorts (Area Under Curves (AUCs) of 0.878, 0.925, and 0.896 in the TCGA cohort and 0.781, 0.789, and 0.768 in the CGGA cohort for 1-, 3-, and 5-year predicted survival, respectively; Figure 10C, D). Further univariate and multivariate analyses incorporated other clinical features and showed that the risk score is a prognostic risk factor for patients with glioma independent of other clinical characteristics, such as age, tumor grade, and recurrence status (Figure 10E–H).\u003c/p\u003e\n\u003cp\u003eIn addition, the expression levels of 12 risk genes in 9 single-cell clusters were plotted (Figure 11). ZCCHC24 and SULF2 were overexpressed mainly in cluster 3, SERPINE2 was chiefly expressed in cluster 6, and MT1X was expressed in cluster 1. KDM1A and IGFBP5 were present in cluster 4; MMP19, G0S2, CH25H, and C1orf162 were upregulated mostly in clusters 0 and 2; and TSPAN3 and ITGB8 were extensively expressed in multiple clusters, though these genes were mainly detected in clusters 1 and 4.\u003c/p\u003e\n\u003ch2\u003e3.9\u0026nbsp; \u0026nbsp; \u0026nbsp;Nomogram to qualitatively predict OS in patients with glioma\u003c/h2\u003e\n\u003cp\u003eWith the CGGA dataset, six prognostic factors, including age, sex, tumor grade, recurrence status, IDH mutation status, and 1p19q codeletion status were combined with the risk score to construct a nomogram for predicting the 1-, 3-, and 5-year survival probabilities (Figure 12).\u0026nbsp;\u003c/p\u003e"},{"header":"4 Discussion","content":"\u003cp\u003eAs one of the most lethal brain tumors, gliomas exhibit an extremely dreadful prognosis. Even with complete standard therapy, patients often quickly succumb to the disease owing to therapy resistance and tumor recurrence (\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e). Tumor heterogeneity, including genetic, epigenetic, and TME heterogeneity, contributes greatly to this doomed prognosis. Glioma heterogeneity is difficult to study with bulk genomic approaches but is better dissected using scRNA-seq (24). An in-depth understanding of the TME and related immune components via scRNA-seq analysis may facilitate the development and improvement of immunotherapy (\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eIn this study, starting from a glioma single-cell sequencing dataset, we divided single glioma cells into five state trajectories through dimensionality reduction and pseudo-temporal analysis. Then, we identified the differential genes between trajectories. Based on the TRGs, the patients in the CGGA dataset were separated into three clusters. The clinical features, survival prognosis, TME, immune infiltration, and ICG expression levels were significantly different among the clusters. Survival analysis identified immune cells and ICGs with clinical significance. Protein-protein interaction analysis and functional enrichment revealed the involvement of immune suppression pathways. Subsequently, clinically relevant TRGs were obtained from WGCNA analysis and intersected with glioma DEGs to identify candidate genes for establishing a prognosis prediction model. LASSO regression, multivariate regression, survival analysis, and ROC analysis were used to identify risk genes. Then, a prognosis prediction model was constructed and evaluated. This model has acceptable predictive efficacy. Finally, a nomogram was built to quantify the contribution of each risk factor in specific patients.\u003c/p\u003e \u003cp\u003eGiven the widespread success of immunomodulatory therapy across diverse cancer types, characterizing the immune fraction of the glioma TME and developing corresponding immunotherapies is of significant interest (\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e). The brain represents a relatively immunosuppressive microenvironment for glioma cells because it contains unique immune cells and regulatory cell types such as microglia and astrocytes, lacks classic lymphoid components, and is protected from circulating immune cells by the blood\u0026ndash;brain barrier (\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e). Tumor-associated macrophages (TAMs) are essential glioma TME immune cells that can arise from two distinct sources: peripheral bone marrow-derived monocytes and yolk sac-derived microglia (\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e, \u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e). TAMs are enriched in gliomas and can adopt multiple pro-tumorigenic phenotypes that play key roles in promoting invasion, angiogenesis, metastasis, and immunosuppression (\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e). Survival analysis of the immune cell signature revealed that the M0 macrophage signature significantly correlated with a worse prognosis. M0 macrophage infiltration, but not that of M1 or M2 macrophages, was increased in cluster C3. This result implies that macrophage M0 infiltration might partly contribute to the poor survival of patients with glioma. Macrophages are characterized by their plasticity and are commonly divided into classic-activated macrophages (M1) that stimulate anti-tumor immune responses and alternatively activated macrophages (M2), which suppress immune responses (\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e). M0 macrophages are non-polarized na\u0026iuml;ve macrophages. The M1 and M2 phases of TAMs are widely recognized in numerous cancers, while M2 macrophages have been reported mainly in high-grade glioma. However, a recent study reported the enrichment of nonpolarized M0 macrophages rather than M1 or M2 macrophage assemblies in glioblastoma (\u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e, \u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e). Our results are consistent with these findings and suggest that inducing macrophage polarization from M0 to M1 cells might have a therapeutic effect on gliomas.\u003c/p\u003e \u003cp\u003eFurther clustering of the CGGA cohort demonstrated that TRGs distinguish different clusters of patients with different OS and various TME compositions. Notably, cluster C3 showed the worst prognosis and had the lowest tumor purity and the highest immune score, indicating that immune cell infiltration in gliomas could be an adverse prognostic factor. This is consistent with their general immunosuppression functions. Thus, evaluating the TME composition using TRGs might have prognostic value for patients with glioma.\u003c/p\u003e \u003cp\u003eMany other immune cells also play essential roles in glioma. Increased migration of dendritic cells (DCs) to nearby draining lymph nodes can lead to prolonged progression-free and improved OS in patients with glioblastoma (\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e). Neutrophils contribute to glioblastoma progression by supporting glioma stem cell pool expansion (\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e). CD8\u0026thinsp;+\u0026thinsp;T cells function as anti-tumor components, and their activation with immune checkpoint inhibitors leads to longer survival in mouse glioma models (\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e). CD4\u0026thinsp;+\u0026thinsp;lymphocyte enrichment, which could enhance the function of CD8\u0026thinsp;+\u0026thinsp;T cells, profoundly decreases in glioma. Within the reduced number of CD4\u0026thinsp;+\u0026thinsp;lymphocytes, CD4\u0026thinsp;+\u0026thinsp;helper T (Th) cell numbers are reduced, while immunosuppressive Tregs increase (\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e). This could possibly explain why cluster C3 showed the highest fraction of anti-tumor CD8\u0026thinsp;+\u0026thinsp;T cells but still displayed the worst prognosis. The survival advantage conferred by CD8\u0026thinsp;+\u0026thinsp;T cells was likely offset by the immunosuppressive effects of Treg cells, which were significantly enriched in cluster C3.\u003c/p\u003e \u003cp\u003eImmune checkpoints, regulators of the immune system, are crucial for self-tolerance and prevent the immune system from indiscriminately attacking cells. However, this mechanism could also be employed by cancer cells to protect themselves from immune attack (\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e). Among all the studied ICGs, CTLA4 and PD-L1 (CD274) are the most well-known. CTLA4 functions by competitive binding to CD80/86 ligands in DCs to prevent CD4\u0026thinsp;+\u0026thinsp;helper T cell activation while boosting Treg, thereby leading to a pro-tumor immunosuppressive phenotype (\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e). PD-L1 and PD-L2 are expressed by tumor cells, and their interaction with PD-1 in immune cells can reduce immune responses (\u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e, \u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e44\u003c/span\u003e). Numerous studies have been performed on other ICGs in glioma. A higher LDHA/LDHB ratio is associated with the IDH wild-type genotype and poorer OS in patients with glioma (\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e). B2M expression is correlated with malignancy of glioma (\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e). Inhibiting the ICOSLG axis in GBM may provide promising immunotherapeutic effects (\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e). LGALS9 regulates the immunosuppressive features of glioma (\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e). Depleting endogenous PVR from human glioma cells can inhibit tumor migration (\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e). HAVCR2 (TIM-3) and SIGLEC15 are well-described immunosuppressive components (\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e, \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e). Consistently, CD274, CTLA4, HAVCR2, PDCD1, PDCD1LG2, SIGLEC15, LDHA, B2M, ICOSLG, LGALS9, and PVR were up-regulated in cluster C3 samples, correlating with a worse prognosis. Likewise, LDHB, IL12, and LAMA3, which are anti-tumor ICGs and are associated with better survival (\u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e, \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e53\u003c/span\u003e), were downregulated in cluster C3.\u003c/p\u003e \u003cp\u003eHowever, the roles of some ICGs are still not well defined. High PTPRC expression in glioma could result in both pro- and anti-immune responses (\u003cspan citationid=\"CR54\" class=\"CitationRef\"\u003e54\u003c/span\u003e). TNFRSF9 is upregulated in human gliomas compared to that in normal brain tissue (\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e), but its effect on gliomas has not been extensively studied. Further, TNFSF4 has not yet been studied in glioma. Our study sheds light on their possible roles in gliomas and uncovers novel potential therapeutic directions.\u003c/p\u003e \u003cp\u003eInterestingly, CD8A, CD28, CD40, CD40LG, and CD80 are essential immune checkpoint components associated with T cell activation. These proteins are designated as anti-tumor immune markers in many extracranial solid tumors (\u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e, \u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e56\u003c/span\u003e). However, in the present study, higher expression of these markers was associated with worse prognosis. Whether these proteins function differently in the unique glioma microenvironment or whether they are upregulated due to negative feedback from the specific immunosuppressive state of gliomas requires further study.\u003c/p\u003e \u003cp\u003eDifferent clusters based on TRGs in the CGGA cohort showed distinct survival profiles, suggesting that TRG-based patient classification can be used to predict patient OS. Various TRG expression patterns indicate the diverse composition of the TME and immune components. Successfully predicting patient prognosis using TRG clustering indirectly shows that alterations in the TME and immunological components have a significant impact on the prognosis of patients with glioma. Thus, TRGs were used to identify risk genes and construct a prognosis prediction model. Our model had high accuracy and efficiency in both the training and testing datasets. To the best of our knowledge, this is the first TRG-based signature used to predict prognosis in patients with glioma. Moreover, a nomogram combining the TRG-based risk score and prognostic clinicopathological variables was constructed to provide a visual method for predicting patient OS. The developed nomogram was more accurate and effective than using the risk score alone.\u003c/p\u003e \u003cp\u003eThis study has several limitations. First, the correlation between glioma prognosis and immune components (immune cells and ICGs) was determined using bioinformatic analysis. The verification of conclusion and the in-depth exploration of internal mechanism depend on further molecular biology experiments. Second, the cell state trajectories and TRGs constructed by the scRNA-seq analysis were derived from the same GEO dataset and therefore lack systematic validation of repeated multi-sample experiments. In addition, the GEO, CGGA, and TCGA data used to construct the prognostic model were retrospective data. Large-scale prospective studies are necessary to assess the effectiveness and practicality of our models.\u003c/p\u003e"},{"header":"5 Conclusions","content":"\u003cp\u003eWe identified five state trajectories of single glioma cells with distinct cell identities based on scRNA-seq and demonstrated that TRG-based genotyping can accurately predict patient OS, clinicopathological features, TME composition, immune infiltration status, and ICG expression. Immune cells and ICGs that are correlated with prognosis of patients with glioma were identified and could be promising therapeutic targets. We developed a nomogram that integrates the TRG-based risk score signature and clinicopathological variables to provide a preliminary quantitative method for predicting patient OS. In summary, this study highlighted alterations in the TME and immune components in glioma, identified a TRG-based risk gene signature for predicting patient clinical outcomes, and explored associated immune responses to uncover potential immunotherapy directions for patients with glioma.\u003c/p\u003e"},{"header":"Declarations","content":"\u003ch2\u003eConflict of Interest\u003c/h2\u003e \u003cp\u003eThe authors declare no conflicts of interest.\u003c/p\u003e \u003ch2\u003eFunding\u003c/h2\u003e \u003cp\u003eThis study was supported by the Applied Basic Research joint project of Science and Technology Department of Yunnan Province and Kunming Medical University (No. 202364).\u003c/p\u003e\u003ch2\u003eAuthor Contribution\u003c/h2\u003e\u003cp\u003eConceptualization, D.C. and J.F.; methodology, Z.F.; software, D.C. and Y.H.; validation, D.C. and L.D.; formal analysis, Y.H.; investigation, D.C. and Y.H.; resources, Y.H.; data curation, Y.H.; writing\u0026mdash;original draft preparation, D.C. and Y.H.; writing\u0026mdash;review and editing, D.C. and Y.H.; visualization, D.C.; supervision, D.C. and Y.H.; project administration, J.F.; and funding acquisition, J.F. All authors have read and agreed to the published version of the manuscript.\u003c/p\u003e\u003ch2\u003eData Availability Statement\u003c/h2\u003e \u003cp\u003eLinks to publicly archived datasets analyzed: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE84465\u003c/span\u003e\u003cspan address=\"https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE84465\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e; \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://gepia2.cancer-pku.cn/#index\u003c/span\u003e\u003cspan address=\"http://gepia2.cancer-pku.cn/#index\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e; \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://xenabrowser.net/hub/\u003c/span\u003e\u003cspan address=\"https://xenabrowser.net/hub/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e; and \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.cgga.org.cn/\u003c/span\u003e\u003cspan address=\"http://www.cgga.org.cn/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e, accessed on Jan 12, 2022.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eOstrom QT, Cote DJ, Ascha M, Kruchko C, Barnholtz-Sloan JS. Adult Glioma Incidence and Survival by Race or Ethnicity in the United States From 2000 to 2014. JAMA Oncol. 2018;4(9):1254\u0026ndash;62.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOstrom QT, Gittleman H, Truitt G, Boscia A, Kruchko C, Barnholtz-Sloan JS. CBTRUS Statistical Report: Primary Brain and Other Central Nervous System Tumors Diagnosed in the United States in 2011\u0026ndash;2015. Neuro Oncol. 2018;20(suppl_4):iv1-iv86.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWesseling P, Capper D. WHO 2016 Classification of gliomas. Neuropathol Appl Neurobiol. 2018;44(2):139\u0026ndash;50.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eQi Y, Liu B, Sun Q, Xiong X, Chen Q. Immune Checkpoint Targeted Therapy in Glioma: Status and Hopes. Front Immunol. 2020;11:578877.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYang K, Wu Z, Zhang H, Zhang N, Wu W, Wang Z, et al. Glioma targeted therapy: insight into future of molecular approaches. Mol Cancer. 2022;21(1):39.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLouveau A, Harris TH, Kipnis J. Revisiting the Mechanisms of CNS Immune Privilege. Trends Immunol. 2015;36(10):569\u0026ndash;77.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBailey SL, Carpentier PA, McMahon EJ, Begolka WS, Miller SD. Innate and adaptive immune responses of the central nervous system. Crit Rev Immunol. 2006;26(2):149\u0026ndash;88.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDeCordova S, Shastri A, Tsolaki AG, Yasmin H, Klein L, Singh SK, et al. Molecular Heterogeneity and Immunosuppressive Microenvironment in Glioblastoma. Front Immunol. 2020;11:1402.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWeber JS, D'Angelo SP, Minor D, Hodi FS, Gutzmer R, Neyns B, et al. Nivolumab versus chemotherapy in patients with advanced melanoma who progressed after anti-CTLA-4 treatment (CheckMate 037): a randomised, controlled, open-label, phase 3 trial. Lancet Oncol. 2015;16(4):375\u0026ndash;84.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMotzer RJ, Escudier B, McDermott DF, George S, Hammers HJ, Srinivas S, et al. Nivolumab versus Everolimus in Advanced Renal-Cell Carcinoma. N Engl J Med. 2015;373(19):1803\u0026ndash;13.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMedikonda R, Dunn G, Rahman M, Fecci P, Lim M. A review of glioblastoma immunotherapy. J Neurooncol. 2021;151(1):41\u0026ndash;53.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eEberwine J, Sul JY, Bartfai T, Kim J. The promise of single-cell sequencing. Nat Methods. 2014;11(1):25\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGawad C, Koh W, Quake SR. Single-cell genome sequencing: current state of the science. Nat Rev Genet. 2016;17(3):175\u0026ndash;88.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePapalexi E, Satija R. Single-cell RNA sequencing to explore immune cell heterogeneity. Nat Rev Immunol. 2018;18(1):35\u0026ndash;45.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLall S, Sinha D, Bandyopadhyay S, Sengupta D. Structure-Aware Principal Component Analysis for Single-Cell RNA-seq Data. J Comput Biol. 2018.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSatija R, Farrell JA, Gennert D, Schier AF, Regev A. Spatial reconstruction of single-cell gene expression data. Nat Biotechnol. 2015;33(5):495\u0026ndash;502.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWilkerson MD, Hayes DN. ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics. 2010;26(12):1572\u0026ndash;3.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells with CIBERSORT. Methods Mol Biol. 2018;1711:243\u0026ndash;59.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXiang R, Rong Y, Ge Y, Song W, Ren J, Fu T. Cell differentiation trajectory predicts patient potential immunotherapy response and prognosis in gastric cancer. Aging (Albany NY). 2021;13(4):5928\u0026ndash;45.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLangfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTang Z, Kang B, Li C, Chen T, Zhang Z. GEPIA2: an enhanced web server for large-scale expression profiling and interactive analysis. Nucleic Acids Res. 2019;47(W1):W556-W60.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDarmanis S, Sloan SA, Croote D, Mignardi M, Chernikova S, Samghababi P, et al. Single-Cell RNA-Seq Analysis of Infiltrating Neoplastic Cells at the Migrating Front of Human Glioblastoma. Cell Rep. 2017;21(5):1399\u0026ndash;410.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGhotme KA, Barreto GE, Echeverria V, Gonzalez J, Bustos RH, Sanchez M, et al. Gliomas: New Perspectives in Diagnosis, Treatment and Prognosis. Curr Top Med Chem. 2017;17(12):1438\u0026ndash;47.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNicholson JG, Fine HA. Diffuse Glioma Heterogeneity and Its Therapeutic Implications. Cancer Discov. 2021;11(3):575\u0026ndash;90.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHernandez Martinez A, Madurga R, Garcia-Romero N, Ayuso-Sacido A. Unravelling glioblastoma heterogeneity by means of single-cell RNA sequencing. Cancer Lett. 2022;527:66\u0026ndash;79.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWaldman AD, Fritz JM, Lenardo MJ. A guide to cancer immunotherapy: from T cell basic science to clinical practice. Nat Rev Immunol. 2020;20(11):651\u0026ndash;68.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eQuail DF, Joyce JA. The Microenvironmental Landscape of Brain Tumors. Cancer Cell. 2017;31(3):326\u0026ndash;41.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen Z, Feng X, Herting CJ, Garcia VA, Nie K, Pong WW, et al. Cellular and Molecular Identity of Tumor-Associated Macrophages in Glioblastoma. Cancer Res. 2017;77(9):2266\u0026ndash;78.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBowman RL, Klemm F, Akkari L, Pyonteck SM, Sevenich L, Quail DF, et al. Macrophage Ontogeny Underlies Differences in Tumor-Specific Education in Brain Malignancies. Cell Rep. 2016;17(9):2445\u0026ndash;59.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMantovani A, Marchesi F, Malesci A, Laghi L, Allavena P. Tumour-associated macrophages as treatment targets in oncology. Nat Rev Clin Oncol. 2017;14(7):399\u0026ndash;416.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGraeber MB, Scheithauer BW, Kreutzberg GW. Microglia in brain tumors. Glia. 2002;40(2):252\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShi C, Pamer EG. Monocyte recruitment during infection and inflammation. Nat Rev Immunol. 2011;11(11):762\u0026ndash;74.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGinhoux F, Schultze JL, Murray PJ, Ochando J, Biswas SK. New insights into the multidimensional concept of macrophage ontogeny, activation and function. Nat Immunol. 2016;17(1):34\u0026ndash;40.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMurray PJ, Allen JE, Biswas SK, Fisher EA, Gilroy DW, Goerdt S, et al. Macrophage activation and polarization: nomenclature and experimental guidelines. Immunity. 2014;41(1):14\u0026ndash;20.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGabrusiewicz K, Rodriguez B, Wei J, Hashimoto Y, Healy LM, Maiti SN, et al. Glioblastoma-infiltrated innate immune cells resemble M0 macrophage phenotype. JCI Insight. 2016;1(2).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHuang L, Wang Z, Chang Y, Wang K, Kang X, Huang R, et al. EFEMP2 indicates assembly of M0 macrophage and more malignant phenotypes of glioma. Aging (Albany NY). 2020;12(9):8397\u0026ndash;412.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMitchell DA, Batich KA, Gunn MD, Huang MN, Sanchez-Perez L, Nair SK, et al. Tetanus toxoid and CCL3 improve dendritic cell vaccines in mice and glioblastoma patients. Nature. 2015;519(7543):366\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiang J, Piao Y, Holmes L, Fuller GN, Henry V, Tiao N, et al. Neutrophils promote the malignant glioma phenotype through S100A4. Clin Cancer Res. 2014;20(1):187\u0026ndash;98.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCohen JV, Kluger HM. Systemic Immunotherapy for the Treatment of Brain Metastases. Front Oncol. 2016;6:49.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFecci PE, Mitchell DA, Whitesides JF, Xie W, Friedman AH, Archer GE, et al. Increased regulatory T-cell fraction amidst a diminished CD4 compartment explains cellular immune defects in patients with malignant glioma. Cancer Res. 2006;66(6):3294\u0026ndash;302.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePardoll DM. The blockade of immune checkpoints in cancer immunotherapy. Nat Rev Cancer. 2012;12(4):252\u0026ndash;64.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWing K, Onishi Y, Prieto-Martin P, Yamaguchi T, Miyara M, Fehervari Z, et al. CTLA-4 control over Foxp3\u0026thinsp;+\u0026thinsp;regulatory T cell function. Science. 2008;322(5899):271\u0026ndash;5.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFreeman GJ, Long AJ, Iwai Y, Bourque K, Chernova T, Nishimura H, et al. Engagement of the PD-1 immunoinhibitory receptor by a novel B7 family member leads to negative regulation of lymphocyte activation. J Exp Med. 2000;192(7):1027\u0026ndash;34.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLatchman Y, Wood CR, Chernova T, Chaudhary D, Borde M, Chernova I, et al. PD-L2 is a second ligand for PD-1 and inhibits T cell activation. Nat Immunol. 2001;2(3):261\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGivechian KB, Garner C, Benz S, Rabizadeh S, Soon-Shiong P. Glycolytic expression in lower-grade glioma reveals an epigenetic association between IDH mutation status and PDL1/2 expression. Neurooncol Adv. 2021;3(1):vdaa162.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang H, Cui B, Zhou Y, Wang X, Wu W, Wang Z, et al. B2M overexpression correlates with malignancy and immune signatures in human gliomas. Sci Rep. 2021;11(1):5045.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eIwata R, Hyoung Lee J, Hayashi M, Dianzani U, Ofune K, Maruyama M, et al. ICOSLG-mediated regulatory T-cell expansion and IL-10 production promote progression of glioblastoma. Neuro Oncol. 2020;22(3):333\u0026ndash;44.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiang T, Wang X, Wang F, Feng E, You G. Galectin-9: A Predictive Biomarker Negatively Regulating Immune Response in Glioma Patients. World Neurosurg. 2019;132:e455-e62.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSloan KE, Stewart JK, Treloar AF, Matthews RT, Jay DG. CD155/PVR enhances glioma cell dispersal by regulating adhesion signaling and focal adhesion dynamics. Cancer Res. 2005;65(23):10930\u0026ndash;7.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHuang YH, Zhu C, Kondo Y, Anderson AC, Gandhi A, Russell A, et al. CEACAM1 regulates TIM-3-mediated tolerance and exhaustion. Nature. 2015;517(7534):386\u0026ndash;90.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMorad G, Helmink BA, Sharma P, Wargo JA. Hallmarks of response, resistance, and toxicity to immune checkpoint blockade. Cell. 2021;184(21):5309\u0026ndash;37.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYan J, Smyth MJ, Teng MWL. Interleukin (IL)-12 and IL-23 and Their Conflicting Roles in Cancer. Cold Spring Harb Perspect Biol. 2018;10(7).\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTang L, Wang P, Wang Q, Zhong L. Correlation of LAMA3 with onset and prognosis of ovarian cancer. Oncol Lett. 2019;18(3):2813\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBrandenburg S, Turkowski K, Mueller A, Radev YT, Seidlitz S, Vajkoczy P. Myeloid cells expressing high level of CD45 are associated with a distinct activated phenotype in glioma. Immunol Res. 2017;65(3):757\u0026ndash;68.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBlank AE, Baumgarten P, Zeiner P, Zachskorn C, Loffler C, Schittenhelm J, et al. Tumour necrosis factor receptor superfamily member 9 (TNFRSF9) is up-regulated in reactive astrocytes in human gliomas. Neuropathol Appl Neurobiol. 2015;41(2):e56-67.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang Y, Huang T, Hu Y, Wang Y. Activation of CD40 by soluble recombinant human CD40 ligand inhibits human glioma cells proliferation via nuclear factor-kappaB signaling pathway. J Huazhong Univ Sci Technolog Med Sci. 2012;32(5):691\u0026ndash;6.\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"single-cell RNA sequencing, immune response, glioma, prognosis, in silico analysis","lastPublishedDoi":"10.21203/rs.3.rs-3829624/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-3829624/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eTo investigate the microenvironment composition of gliomas and the associated clinical significance, we analyzed single-cell RNA sequencing and bulk RNA-seq data from glioma samples. Cell trajectory analysis identified five trajectories with distinct cell states and corresponding trajectory-related genes (TRGs). TRG-based clustering segregated patients with glioma with different overall survival, clinicopathological features, immune infiltration status, and immune checkpoint gene (ICG) expression levels. Notably, a worse prognosis was seen in patients with a higher immune score, lower tumor purity, higher M0 macrophage and regulatory T (Treg) cell infiltration, and increased ICG expression. Further survival analysis and functional enrichment analysis revealed a close relationship between prognosis and ICG-associated immunosuppressive pathways. Candidate prognostic genes were obtained using WGCNA analysis and differential expression analysis. LASSO and multivariate regression analysis were used to establish a prognostic prediction model. The prognostic risk-scoring signature including 12 genes successfully predicted patient survival with acceptable AUC values. A nomogram was constructed to evaluate the contribution of the risk signature to patient prognosis. This study highlights the potential involvement of tumor microenvironment variation and immune alteration in glioma progression and establishes a TRG-based prognostic model to predict patient clinical outcomes.\u003c/p\u003e","manuscriptTitle":"Tumor microenvironment variation implicates immune alterations and correlates with prognosis in patients with glioma","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2024-01-05 18:01:31","doi":"10.21203/rs.3.rs-3829624/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"c0fd318d-25df-4d53-8e78-8051dec82cb8","owner":[],"postedDate":"January 5th, 2024","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[],"tags":[],"updatedAt":"2024-01-23T06:14:38+00:00","versionOfRecord":[],"versionCreatedAt":"2024-01-05 18:01:31","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-3829624","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-3829624","identity":"rs-3829624","version":["v1"]},"buildId":"_2-kVJe1T_tPrBINL-cwx","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}