Skip to main content
Download PDF

Multi-transcriptomics reveals niche-specific expression programs and endothelial cells in glioblastoma

Journal of Translational Medicine volume 23, Article number: 444 (2025) Cite this article

This article has been updated

Abstract

Background

Glioblastoma (GBM) is a highly lethal malignant intracranial tumor, distinguished from low-grade glioma by histopathological hallmarks such as pseudopalisading cells around necrosis (PAN) and microvascular proliferation (MVP). To date the spatial organization of the molecular and cellular components of these specific histopathological features has not been fully elucidated.

Methods

Here, using bulk RNA sequencing, spatial transcriptomic and single cell RNA sequencing (scRNA-seq) data of GBM patients, we identified niche-specific transcriptional programs and characterized the differences in molecular expression and cellular organization between PAN and MVP.

Results

Notably, we discovered spatially distinct domains within the tumor core and identified niche-specific signatures: NDRG1 and EPAS1, specifically expressed in the PAN and MVP regions. The clustering results showed two distinct phenotypes of endothelial cells (ECs) were enriched in the MVP and PAN regions, respectively. PAN-associated endothelial cells exhibit copy number variations similar to those in GBM cells. Single cell trajectory analysis reveals a pseudotime trajectory, indicating the differentiation of glioblastoma stem cells (GSCs) toward ECs.

Conclusions

Necrosis cores which are surrounded by hypoxic and perivascular niches and microvascular proliferation area within the glioblastoma tumor microenvironment, have been considered as standardized morphological indicators of aggressive GBM. Our findings provide a cellular and molecular insights into GBM progression.

Introduction

Glioblastoma (GBM) is a highly malignant, rapidly progressive brain tumor that is pathologically characterized by pseudopalisading cells around necrosis (PAN) and microvascular proliferation (MVP) [1, 2]. These histological features distinguish GBM from lower-grade brain tumors and have long been recognized as ominous prognostic features.

PAN refers to a hypercellular zone surrounding necrotic areas in GBM, while MVP is a form of angiogenesis morphologically recognized as endothelial proliferation within newly sprouted vessels. PAN and MVP have been considered as indicators of the high malignancy and invasive capacity of GBM. However, the pathogenesis of this organizational structure remains unclear. One generally accepted theory is that GBM cells overgrowth leads to intratumorally vascular damage, resulting in endothelial injury that exacerbates tissue hypoxia. Malignant cells migrate away from the hypoxic environment, forming waves that move toward the periphery and appear as pseudoapoptotic cells under the microscope. Additionally, tumor cells in the hypoxic area secrete vascular endothelial growth factor (VEGF) and IL-8, promoting microvascular proliferation and accelerating the outward expansion of tumor cells alongside new blood vessel formation. Despite considerable advances in molecular mechanisms involved in GBM tumorigenesis in recent decades, the molecular and cellular heterogeneity of these hallmark regions and their spatial characteristics have not been fully elucidated [3, 4].

GBM contains a host of non-tumor components like mesenchymal cells, vasculature, cytokines, immune cells, and fibroblasts, constituting a complex tumor microenvironment (TME) with many anatomically distinct regions. Necrosis cores which are surrounded by hypoxic and perivascular niches within the GBM TME, have been considered as a standardized morphological indicator of aggressive GBM. However, the pathogenesis of these regions in malignant GBM is not fully understood. Therefore, clarifying the molecular characteristics of specific regions and targeting therapies based on these signatures may prove a promising treatment strategy. The Ivy Glioblastoma Atlas Project (Ivy GAP) has collected genomic and transcriptomic profiles of anatomical regions like PAN, MVP and cellular tumor (CT). Accordingly, the database can enable understanding of molecular diversity relating to histological heterogeneity. Nevertheless, few studies revealed the spatial heterogeneity of molecular signatures in diverse cellular states from PAN and MVP. The development of spatially resolved transcriptomics (SRT) technology enables us to elaborate the correlation between characteristic gene expression profiles while retaining information of the spatial tissue context and histological structures within tumors.

Here, we integrate RNA sequencing, SRT, and single-cell RNA sequencing (scRNA-seq) to characterize differences between GBM structures, uncovering cellular heterogeneity and molecular stratification of histologic structures within GBMs. We identified prognostic and predictive markers corresponding to GBM localizations.

Methods

Description of the datasets

The RNA-seq datasets (n = 41) from different anatomical regions of primary human glioblastoma (GBM) were obtained from the Ivy Glioblastoma Atlas Project (Ivy-GAP). Additionally, spatially resolved transcriptomics (SRT) datasets (n = 16) of human GBM samples lacking IDH mutation, generated using Visium technology (10X Genomics), were sourced from Ravi et al. (https://doiorg.publicaciones.saludcastillayleon.es/10.5061/dryad.h70rxwdmj). The Smart-seq2 glioblastoma scRNA-seq dataset (n = 28, GSE131928) was acquired from the GEO database. Lastly, the Cancer Genome Atlas (TCGA) GBM dataset (n = 155), containing RNA-sequencing and clinical data, was retrieved from the UCSC Xena browser (https://xenabrowser.net/).

Histological structure based differential gene expression analysis

Gene expression profiles from the CT, PAN, and MVP regions were extracted from the Ivy-GAP datasets. To normalize the gene expression values in Fragments Per Kilobase Million (FPKM) format, we applied a log2(n + 1) transformation of counts using the dplyr package (R 1.0.10). Data quality control and pairwise differential expression analysis among the three groups were conducted using the limma package (R 3.54.0) with a logFoldChange cutoff of 1. Based on the differentially expressed genes (DEGs) identified between pairs, intersection queries were performed using the UPsetR package (R 1.4.0) to visualize shared and unique DEGs. Gene IDs were then converted from SYMBOL to ENTREZID format using the clusterProfiler package (R 4.6.0), and we conducted Gene Ontology (GO) enrichment analysis in the Biological Process (BP), Molecular Function (MF), and Cellular Component (CC) categories. Additionally, KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway enrichment analysis was performed to further interpret the functional roles of DEGs across these regions.

Spatial transcriptomics analysis

For the analysis of spatially resolved transcriptomics, we followed the standard protocol recommended by Seurat (R 4.0). To ensure reproducibility, data from 16 libraries were merged and filtered using consistent criteria: cells with gene counts between 500 and 8000, total transcript counts below 30,000, and a mitochondrial gene percentage below 30%. We applied the SCTransform function to normalize gene expression uniformly across all locations, maintaining default settings for consistency. Additionally, we used both Seurat and the SpaGCN toolkit (https://github.com/jianhuupenn/SpaGCN) to identify spatial variable genes (SVGs), with default parameters applied throughout. This dual approach allowed for cross-validation of SVG identification, enhancing reproducibility by verifying gene expression patterns across spatial domains. These combined methods facilitated the accurate characterization of tissue heterogeneity while ensuring methodological consistency.

Deconvolution

Bhaduri et al. [5] GBM scRNA-seq and Neftel et al. [6] SmartSeq2-derived human GBM scRNA-seq datasets were used as reference datasets. The “harmony” method was used to integrate the dataset from different specimens. The integrated reference is consisting of 40,791 cells, 22,346 genes and 20 cell types. The celltypes with < 1000 cells were filtered. we used the merged datasets as input for the subsequent deconvolution analysis. We used RCTD (v1.1.0, https://github.com/dmcable/spacexr) to deconvolute the transcriptome of each spot (run. RCTD function with the entire mode parameter) into the likely constituent cell types. Deconvoluted ratios of individual cell types were visualized by Stdeconvolve and ggplot2 packages.

Single cell RNA-sequencing data analysis

The scRNA-seq datasets contains 6863 GBM cells and 1067 normal cells from 28 patients. Count tables were analyzed used Seurat 4.0 package following the standard workflow with default settings. Next, InferCNV package (R 1.3.3) was applied to infer CNV from scRNA-seq data, with cutoff = 1. We then combined the CNV classification with the UMAP-based and the gene-set based classifications to define the GSCs and GBM associated macrophage/macroglia. The marker genes for macrophages, T cells, and oligodendrocytes were obtained from Neftel et al. (GSE131928) [9]. To investigate the molecular interaction networks among each cell subtype, CellPhoneDB (R 2.1.7) algorithm was used to infer intercellular communication in the GSC niche.

Immunofluorescence

Clinical glioma paraffin slides were obtained from Wuhan Union Hospital. All samples were collected with signed informed consent according to the internal review and ethics boards of Wuhan Union Hospital. The procedure was approved by the ethics boards of Wuhan Union Hospital. After de-waxing and rehydrating, slides were boiled in Tris-EDTA antigen retrieval buffer (PH = 8.0) for 10 min. Next, sections were treated with 3% H2O2 for 10 min to bleach endogenous peroxidase. After blocking with donkey serum in PBS for 30 min, slides were incubated at 4 °C overnight with primary antibodies. Glioma slides were probed with primary antibodies against the following proteins: CD44 (eBioscience, San Diego, USA), and HIF2A/EPAS1 (Signalway Antibody, Greenbelt, USA), IGFBP7, THY1, NDRG1 and CD31 (Cell Signaling Technology, MA, USA). Sides were then washed with PBS and incubated with DyLight 488/ FluorTM Plus 647 conjugated donkey anti-rabbit, mouse or rat secondary antibody (Invitrogen) for 2 h away from light at room temperature and followed by DAPI staining for another 10 min. After washing away from light, the slides were sealed with an anti-fluorescence quenching sealer and the fluorescence was observed using an Olympus VS200 Slide Scanner.

RNA isolation and qPCR

Total RNA was extracted from HUVEC cell samples subjected to PAN and MVP region environment simulation using the RNAiso Plus reagent (TAKARA, 9109) according to the manufacturer’s instructions. The extracted RNA was dissolved in RNase-free water, and its concentration and purity were determined using a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, USA), ensuring an A260/A280 ratio of 1.8–2.0. Subsequently, the RNA was reverse-transcribed into cDNA using the HiScript III RT SuperMix (Vazyme, R323-01) following the kit protocol. Quantitative real-time PCR was performed using the Taq Pro Universal SYBR qPCR Master Mix (Vazyme, Q712-02) on a real-time PCR instrument. All experiments were performed with three biological replicates and three technical replicates. qPCR data were analyzed using the ΔΔCt method to calculate the relative expression levels of target genes normalized to the GAPDH gene.

Results

Anatomical transcriptome analysis of glioblastoma revealed the molecular differences between PAN and MVP

To dissect heterogeneity within MVP and PAN anatomical structures hidden in bulk analysis and explore the distinct molecular signatures associated with GBM malignancy, we performed integrative bioinformatics analysis. This analysis incorporated histological structure-specific transcriptional profiles, spatial transcriptome datasets and single-cell transcriptome datasets (as illustrated in Fig. 1A).

In the Ivy-GAP dataset, we compared the profiles of 270 laser-microdissected GBM tissues, including samples from 111 cortical (CT), 22 hyperplastic blood vessels (HBV), 40 PAN, 28 MVP, 26 perinecrotic zone (PNZ), 24 infiltrating tumor (IT), and 19 leading edge (LE) regions (Fig. 1B). We compared expression profiles between perivascular regions (HBV and MVP) and CT, as well as between perinecrotic regions (PNZ and PAN) and CT, using criteria of > 2-fold change and FDR < 0.05. PAN and PNZ showed similar expression profiles, as did HBV and MVP, reflecting their histologic similarities. Notably, PAN and PNZ had elevated expression of angiogenesis-promoting genes, including IL-8, ADM, VEGFA, and NDRG1, compared to other regions. HBV and MVP exhibited high expression of endothelial and stemness-associated genes such as ESM1 and CRIP1 [7, 8]. Additionally, MVP exhibited increased expression of tumorigenesis and migration-related molecules like HIGD1B and MMP9. Differentially expressed genes were more pronounced in PAN and MVP compared to PNZ and HBV. A Venn diagram of DEGs highlighted 312 genes specifically up-regulated in PAN and 1048 in MVP (Fig. S1A).

Interestingly, 151 down-regulated genes overlapped between PAN and MVP, comprising one-third of PAN downregulation. These were involved in nervous system development and cell adhesion (Fig. S1B and S1C). Gene Ontology and KEGG analysis of MVP DEGs showed association with ameboidal migration and pathways like PI3K-Akt and focal adhesion (Fig. S1D). PAN DEGs were enriched in hypoxia response and HIF-1 signaling (Fig. S1E). This suggests cells in PAN and MVP are highly malignant, with reduced intercellular adhesion facilitating invasion and aggressiveness.

Fig. 1
figure 1

Characterizing intra-tumoral expression heterogeneity in GBM by distinct molecular profiles. A Schematic representation of processing and integration of multi-level RNA sequencing data. RNA-seq datasets (n = 41), Spatially resolved transcriptomics datasets (SRT, n = 28) scRNA-seq dataset (n = 28, GSE131928). B A heatmap of differentially expressed genes (rows) in distinct anatomical regions (columns). Expression differences in perivascular areas (including hyperplastic blood vessels [HBV], and areas of microvascular proliferation [MVP]) with cellular tumor (CT) and perinecrotic areas (perinecrotic zones [PNZ], pseudopalisading cells around necrosis [PAN]) with CT were compared in pairs (> 2 fold change, FDR < 0.05). IT, infiltrating tumor; LE, leading edge. C, D Venn diagrams representing the overlap between significantly up-regulated genes in Ivy GAP data, Seurat and SpaGCN spatially variable genes (SVGs) for PAN and MVP. E, F Expression of PAN and MVP associated signature genes in GBM patients. Genes with domain numbers greater than 3 points were listed. The bars at the top of the tile graph indicate the number of spatially variable expressed signature genes in each patient, and the bars on the right indicate the number of patients with spatially variable expression of each gene. The color indicates the level of expression

Full size image

Identification of domain-specific marker genes of PAN and MVP

To further characterize the spatial organization of morphologically distinct regions within GBM tumors, we employed the Seurat and SpaGCN packages to identify robust spatially variable genes (SVGs) in 10X Visium datasets. We integrated spots from the same cluster in each sample into pseudobulks using Seurat (Fig. S2).

To reveal the spatial distribution of upregulated genes, SRT datasets from 16 GBM patients were integrated, identifying 1000 and 524 SVGs through Seurat and SpaGCN, respectively (Table S1). Overlapping with PAN and MVP upregulated genes yielded 63 PAN and 53 MVP signature gene candidates (Fig. 1C and D). These were ranked by frequency in patient SVGs. The top PAN genes were IGFBP5, MTIX, CHI3L1, NDRG1, and VEGF. Top MVP genes were IGFBP, SPARC, IFITM3, IGFBP3 and A2M (Fig. 1E and F). Most of these associate with migration, invasion, stemness [9,10,11] and angiogenesis [12,13,14], which is consistent with the histological features and functions of PAN and MVP.

TCGA RNA seq data showed most signature genes had higher expression in GBM versus LGG (Fig. S3). Furthermore, high expression associated with poor prognosis (Fig. S4), including IGFBP5, IGFBP3, IGFBP7 and CHI3L1. This suggests PAN and MVP spatial upregulation of these genes may promote GBM malignancy.

The spatial localization of PAN and MVP signature genes

To clarify whether the PAN and MVP signature genes have specific spatial expression patterns, we compared their spatial expression in GBM patients. We selected UFK243 as a typical sample with obvious necrosis and vascular proliferation (Fig. 2A). UFK243 spots were clustered into 6 domains by Seurat and SpaGCN. We observed elevated APOC1 and NDRG1 in domain 0, EPAS1 and VEGFA in domain 1, and IGFBP7, CHI3L1 and TAGLN in domain 5 (Fig. 2B). PAN genes like IGFBP5, MT1X, CHI3L1 and NDRG1 had similar spatial patterns, predominantly in domains 0 and 5. Meanwhile, IGFBP7, IFITM3, A2M and EPAS1 were sporadic in domains 1, 3 and 5 (Fig. 2C). PAN and MVP genes also showed distinct spatial patterns in UFK243 (Fig. S5).

Fig. 2
figure 2

The spatial localization of PAN and MVP signature genes in GBM. A Representative section showing different domains based on H&E-staining and spatially variable genes (SVGs). B Characteristically expressed genes in each domain. C Gene expression for eight genes on a fine grid set of spatial locations in patient UKF423

Full size image

In situ hybridization from the Ivy GAP further verified elevated NDRG1 and EPAS1 in PAN (Fig. 3A). Immunofluorescence also validated the differences, showing high levels of NDRG1 and HIF2A in PAN (Fig. 3B). Consistently, immunofluorescence (IF) results from two additional clinical patient samples showed that IGFBP7 and THY1 were enriched around CD31-positive blood vessels (Fig. S6). These results confirm the accuracy of the obtained region-specific signature genes.

Fig. 3
figure 3

Expression of NDRG1 and SPAS1 in GBM. A In situ hybridization (ISH) images of Ivy GAP database show PAN expression patterns of NDRG1 (left) and MVP expression patterns of EPAS1 (right) in GBM tissues. B Representative multiplex immunofluorescence images of NDRG1 (green), HIF2A (green) and CD31 (pink), CD44 (red) in GBM slides

Full size image

Identification of cell type composition in spatially resolved transcriptomic GBMs

The intricate glioma microenvironment consists of non-neoplastic cell types including endothelial cells, pericytes, microglia, macrophages, and fibroblasts. To comprehensively characterize these, integrated data was visualized by UMAP dimension reduction and unsupervised clustering (Fig. 4A and B). Using SingleR annotation and manual curation, we identified 9 major cell types: astrocytes, macrophage/microglia, neurons, oligodendrocytes, OPCs, fibroblasts, endothelial cells (ECs), T cells, and adipocytes (Fig. 4C). Additionally, spatial deconvolution was performed (Fig. S7). Notably, these methods were consistent, and paired SRT/scRNA analysis showed the identified proportions were comparable to scRNA (Fig. S8).

Fig. 4
figure 4

UMAP visualization of human GBM SRT clusters and Cell type annotations. A UMAP of GBM spatial transcriptome data from 16 patients. Cells are color-coded according to the defined cluster. B UMAP of 16 patients. Cells are color-coded according to patient origin. C Percentage of each cell subset across the 16 patients. D Heat maps showed the expression of PAN (top) and MVP (bottom) regional characteristic genes in different cell types

Full size image

We examined the gene expression of 63 PAN and 53 MVP signature genes across cell types (Fig. 4D). Differences were most prominent in fibroblasts and endothelial cells. MVP genes showed elevated expression in fibroblasts and some endothelial cells, like COL4A1, COL4A2, and FN1. Meanwhile, PAN genes were highly expressed in most other endothelial cells, with CHI3L1 and NAMPT highest in endothelial cells followed by macrophage/microglia. The distinct endothelial expression patterns suggest these cells may play different roles in PAN and MVP histological features.

Exploring diverse endothelial cells in human glioblastoma

PAN and MVP signature module scores also showed relatively independent distribution (Fig. 5A and B), suggest that endothelial cells in these regions have distinct expression patterns and functional states. To understand endothelial heterogeneity, unsupervised clustering was performed (Fig. 5C and D). Differential expression analysis revealed that PAN-associated ECs highly expressed CHI3L1, SNRPD3, CLU, LAMTOR4, and AQP1, while MVP-associated ECs expressed MGP, FABP5, collagens, and modifying enzymes such as COL1A2, COL3A1, COL4A1 and COL4A2 (Fig. 5E and Table S3). Functionally, MVP endothelial cells were particularly active in cell adhesion, angiogenesis, and cytoplasmic translation pathways (Fig. 5F). To validate candidate gene expression in these regions, we conducted an in vitro module to measure mRNA expression levels in PAN and MVP regions. The results of the qPCR analysis indicated that PAN-associated genes, including CHI3L1, SNRPD3, AOP1, GUCD1, APBA2, and PTGDS, were significantly upregulated, consistent with previous analyses (Fig. S9A). Similarly, MVP-associated genes, such as FABP5, COL3A1, COL4A1, COL4A2, and SEC61G, were significantly upregulated in endothelial cells (Fig. S9B). These insights shed light on the distinct roles endothelial cells play in GBM PAN and MVP regions, highlighting the complexity and diversity of the tumor microenvironment.

Fig. 5
figure 5

Different states of endothelial cells in the PAN and MVP regions. A, B UMAP of GBM spatial transcriptome data. Composition of different kinds of cells in PAN and MVP regions was shown. C, D UMAP plot of the endothelial cells in glioma. UMAP plot colored based on the expression of PAN and MVP signatures. E Heatmap showing top 10 differentially expressed genes in PAN and MVP associated endothelial cells. F Heatmap showing top 10 enriched GO terms in different endothelial cells based on differentially expressed genes between PAN (748) and MVP (746)

Full size image

Differentiation of GSCs into endothelial cells

The above results establish a link between PAN/MVP-related characterized genes and different functional states of endothelial cells in GBM. Since SRT lacks single cell resolution, we analyzed 5 GBM scRNA-seq datasets, including 2 with paired SRT. Cell types were identified using Seurat: T cells, microglia/macrophages, oligodendrocytes, GSCs, astrocytes, monocytes, endothelial cells, and B cells. Notably, many PAN/MVP signatures were highly or specifically expressed in endothelial cells, like AKAP12 (Fig. 6A and B). Several genes like MT1X, MT2A, SPARC, EPAS1 and THY1 were highly expressed in both GSCs and the endothelial cells (Fig. 6A and B).

To explore intercellular communications, we visualized ligand-receptor pairs for each cell type (Fig S10). Multiple pairs involved endothelial cells, with the PTN/PTPRZ1 ligand-receptor pair a key mediator between GSCs and endothelium (Fig. 6C). Spatial heatmaps showed PTN and PTPRZ1 have similar localization (Fig. 6D). Pseudotime analysis revealed a transition trajectory from GSCs to endothelial cells (Fig. 6E and F). Correlation analysis confirmed a robust GSC-endothelial correlation (Fig. 6G). Double immunofluorescence validated this, showing colocalization of the GSC marker CD44 and the endothelial marker CD31 (Fig. 6H). This aligns with evidence that GSCs differentiate into tumor vascular endothelium [15, 16]. Together, these techniques facilitated a comprehensive analysis, providing insights into the unique gene expression profiles of GBM’s perivascular and perinecrotic zones, highlighting angiogenesis, stemness, and tumorigenic signatures specific to the MVP and PAN regions.

Fig. 6
figure 6

Differentiation of GSCs into endothelial cells. A, B violin plots showing the normalized expression levels of signature genes across the identified cell types. C Dot plots showing intercellular ligand and receptor pairs between endothelial and other cell types. D Spatial heatmaps of PTN and PTPRZ in patient UKF423. E The t-SNE plot of 9 identified cell types in GBM. F Pseudotime analysis of single-cell RNA-seq data. G The correlation heatmap illustrates the relationship between pairs of cell types infiltrated in GBM. The color is scaled by the correlation value. H Representative multiplex immunofluorescence images of CD31 (pink) and CD44 (red) in GBM. Scale bar represents 20 μm

Full size image

.

Discussion

Receptor tyrosine kinases (RTKs), such as Epidermal Growth Factor Receptor (EGFR) and Platelet-Derived Growth Factor Receptor Alpha (PDGFRA), have been central molecular targets in the development of therapies for GBM [17,18,19]. However, limitations in drug delivery [20] imposed by the blood-brain barrier, along with complex factors in glioblastoma—such as an immunosuppressive tumor microenvironment, low immune infiltration, and high tumor heterogeneity—have significantly hindered the efficacy of immune checkpoint inhibitors in GBM [20,21,22]. Comprehensively characterizing histological features in a spatial context could improve GBM diagnosis and treatment. Single-cell transcriptome technology has revealed invaluable heterogeneity insights and potential therapeutic molecular targets [6, 23,24,25]. SRT can categorize adjacent spots into histological regions, enabling spatial heterogeneity identification between patients.

Due to differences in sampling and measurement approaches among single-cell omics, spatial transcriptomics, and bulk transcriptomics, the biological information captured by each method is not entirely consistent. Integrating and interpreting data from these three sources to extract unified biological conclusions is a significant challenge. To address this, we effectively combined the molecular features of bulk transcriptomics with single-cell clustering information to annotate cell types at individual spots in spatial transcriptomics, thereby predicting single-cell expression profiles in spatial locations. Here, we identified PAN and MVP domain-specific molecules using multi-resolution RNA sequencing, shedding light on subgrouped endothelial cells and GSC differentiation into endothelium in GBM.

Intraoperative necrosis, tumor imaging with PET (positron emission tomography) tracer fluoromisonidazole (FMISO) and histological vascular pathology highlighted hypoxia as a unique GBM pathological feature, contributing to aggressiveness [26, 27]. The hypoxia-regulated gene in GBM, N-myc-downstream regulated gene-1 (NDRG1), has been implicated in the regulating of migration, angiogenesis and drug resistance in GBM [12, 28, 29]. NDRG1 overexpression has been shown to reduce angiogenesis and tumor growth in an orthotopic glioma model [12]. Aligning with this, we found high spatial NDRG1 expression in hypoxic PAN but not MVP. We also identified other PAN signature genes like CHI3L1, encoding secreted glycoprotein YKL-40. YKL-40 modulates glioma stem cell states [11]. Recent work by Yanming Ren et al. also highlighted the enrichment of CH3CL1 in hypoxic versus invasive glioma niches [30], validating our results. However, the biological significance of these differentially expressed genes requires further functional validation.

Angiogenesis is a critical for tumor growth and progression. Gliomas, including GBM, are highly vascularized tumors. ECs comprising blood vessels are most critical, forming capillaries that supply glioma cells. Studies using bulk or scRNA-seq have shown EC heterogeneity in GBM patients and models from different tumor sites [31, 32]. Our study introduces a new spatial dimension revealing the distinct functional states of endothelial cells in PAN and MVP regions. MVP ECs exhibited angiogenesis and high collagen expression, while PAN ECs showed hypoxia pathway activation. The MVP EC profile closely resembled the tumor core ECs reported by Yuan Xie et al. [32]. This consistency may be due to lower EC yields from necrotic and peripheral regions in their scRNA-seq preparation (centrifugation and FACS).

Research shows GSCs reside in perivascular niches near ECs [33, 34]. Here, GSCs secrete proangiogenic VEGF, driving EC proliferation, survival, migration, and permeability [33, 35,36,37,38]. Receptor PTPRZ1 is preferentially expressed in GSCs. Inhibiting PTPRZ1 suppresses GBM growth and prolongs survival [39, 40]. TAM-secreted PTN has been shown to stimulate PTPRZ1 on GSCs, promoting malignancy [39]. Our results suggest that PTN/PTPRZ1 and MDK/PTPRZ1 can also enable ECs-GSC communication.

It is believed GSCs can differentiate into pericyte-like cells supporting vessels and tumor growth [16, 41]. Recent live-cell imaging shows GSC also transdifferentiate into endothelial progenitor cells and ECs [42, 43], potentially explaining unsatisfactory antiangiogenic glioma therapy outcomes. Our pseudotime analysis and EC-GSC marker colocalization support this.

Conclusions

This study demonstrates two distinct phenotypes of endothelial cells were enriched in the MVP and PAN regions, respectively. PAN associated ECs possess copy number variation similar to GBM cells. Single cell trajectory infers a pseudotime trajectory of glioblastoma stem cell differentiation toward ECs. Our findings provide a cellular and molecular insights into GBM progression.

Data availability

The Spatially resolved transcriptomics (SRT) datasets were obtained from https://doiorg.publicaciones.saludcastillayleon.es/10.5 061/dryad.h70rxwdmj. The smart-seq2 glioblastoma scRNA-seq dataset (n = 28, GSE131928) was downloaded from the GEO database. The Cancer Genome Atlas (TCGA) GBM dataset (n = 155), including RNA-sequencing, clinical, was obtained from the UCSC Xena browser (https://xenabrowser.net/). No specific code was generated for analysis of these data. If you have any questions regarding the code, feel free to contact the corresponding author of the article at any time.

Change history

  • 12 May 2025

    Affiliation has been corrected.

Abbreviations

GBM:

Glioblastoma

PAN:

Pseudopalisading cells around necrosis

MVP:

Microvascular proliferation

TME:

Tumor microenvironment

CT:

Cellular tumor

SRT:

Spatially resolved transcriptomics

HBV:

Hyperplastic blood vessels

PNZ:

Perinecrotic zone

IT:

Infiltrating tumor

LE:

Leading edge

References

  1. Wen PY, Packer RJ. The 2021 WHO classification of tumors of the Central Nervous System: clinical implications. Neurooncology. 2021;23(8):1215–7.

    Google Scholar 

  2. Horbinski C, et al. Clinical implications of the 2021 edition of the WHO classification of central nervous system tumours. Nat Reviews Neurol. 2022;18(9):515–29.

    Article  Google Scholar 

  3. Barajas RF Jr., et al. Glioblastoma multiforme regional genetic and cellular expression patterns: influence on anatomic and physiologic MR imaging. Radiology. 2010;254(2):564–76.

    Article  PubMed  PubMed Central  Google Scholar 

  4. McBain C, et al. Treatment options for progression or recurrence of glioblastoma: a network meta-analysis. Cochrane Database Syst Rev. 2021;5(1):CD013579.

    PubMed  Google Scholar 

  5. Bhaduri A, et al. Outer radial glia-like Cancer stem cells contribute to heterogeneity of Glioblastoma. Cell Stem Cell. 2020;26(1):48–e636.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  6. Neftel C, et al. An Integrative Model of Cellular States, plasticity, and Genetics for Glioblastoma. Cell. 2019;178(4):835–e84921.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Pan KF, et al. Direct interaction of beta-catenin with nuclear ESM1 supports stemness of metastatic prostate cancer. EMBO J. 2021;40(4):e105450.

    Article  CAS  PubMed  Google Scholar 

  8. Wang J, et al. CRIP1 suppresses BBOX1-mediated carnitine metabolism to promote stemness in hepatocellular carcinoma. EMBO J. 2022;41(15):e110218.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Lin W, et al. IGFBP5 is an ROR1 ligand promoting glioblastoma invasion via ROR1/HER2-CREB signaling axis. Nat Commun. 2023;14(1):1578.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Dixit D, et al. The RNA m6A reader YTHDF2 maintains Oncogene expression and is a targetable dependency in Glioblastoma Stem cells. Cancer Discov. 2021;11(2):480–99.

    Article  CAS  PubMed  Google Scholar 

  11. Guetta-Terrier C, et al. Chi3l1 is a modulator of Glioma Stem Cell States and a therapeutic target in Glioblastoma. Cancer Res. 2023;83(12):1984–99.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Broggini T, et al. NDRG1 overexpressing gliomas are characterized by reduced tumor vascularization and resistance to antiangiogenic treatment. Cancer Lett. 2016;380(2):568–76.

    Article  CAS  PubMed  Google Scholar 

  13. Si M, Lang J. The roles of metallothioneins in carcinogenesis. J Hematol Oncol. 2018;11(1):107.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Lin S, et al. Transcription factor myeloid zinc-finger 1 suppresses human gastric carcinogenesis by interacting with metallothionein 2A. Clin Cancer Res. 2019;25(3):1050–62.

    Article  CAS  PubMed  Google Scholar 

  15. Ricci-Vitiani L, et al. Tumour vascularization via endothelial differentiation of glioblastoma stem-like cells. Nature. 2010;468(7325):824–8.

    Article  CAS  PubMed  Google Scholar 

  16. Hu B, et al. Epigenetic activation of WNT5A drives glioblastoma stem cell differentiation and invasive growth. Cell. 2016;167(5):1281–e129518.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Marin BM, et al. Heterogeneous delivery across the blood-brain barrier limits the efficacy of an EGFR-targeting antibody drug conjugate in glioblastoma. Neuro Oncol. 2021;23(12):2042–53.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Dunn GP, et al. Emerging insights into the molecular and cellular basis of glioblastoma. Genes Dev. 2012;26(8):756–84.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Verhaak RG, et al. Integrated genomic analysis identifies clinically relevant subtypes of glioblastoma characterized by abnormalities in PDGFRA, IDH1, EGFR, and NF1. Cancer Cell. 2010;17(1):98–110.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Samantasinghar A, et al. A comprehensive review of key factors affecting the efficacy of antibody drug conjugate. Biomed Pharmacother. 2023;161:114408.

    Article  CAS  PubMed  Google Scholar 

  21. Wang Q, et al. Tumor evolution of Glioma-Intrinsic Gene Expression Subtypes Associates with immunological changes in the Microenvironment. Cancer Cell. 2018;33(1):152.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Ahmed F et al. Network-based drug repurposing identifies small molecule drugs as immune checkpoint inhibitors for endometrial cancer. Mol Divers, 2024.

  23. Darmanis S, 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.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Richards LM, et al. Gradient of Developmental and Injury Response transcriptional states defines functional vulnerabilities underpinning glioblastoma heterogeneity. Nat Cancer. 2021;2(2):157–73.

    Article  CAS  PubMed  Google Scholar 

  25. Tirosh I, et al. Single-cell RNA-seq supports a developmental hierarchy in human oligodendroglioma. Nature. 2016;539(7628):309–13.

    Article  PubMed  PubMed Central  Google Scholar 

  26. Gerstner ER, et al. ACRIN 6684: Assessment of Tumor Hypoxia in newly diagnosed Glioblastoma using 18F-FMISO PET and MRI. Clin Cancer Res. 2016;22(20):5079–86.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Yamamoto Y, et al. Hypoxia assessed by 18F-fluoromisonidazole positron emission tomography in newly diagnosed gliomas. Nucl Med Commun. 2012;33(6):621–5.

    Article  CAS  PubMed  Google Scholar 

  28. Zhang X et al. N-myc Downstream-Regulated Gene 1 (NDRG1) Regulates Vascular Endothelial Growth Factor A (VEGFA) and Malignancies in Glioblastoma Multiforme (GBM). Biomed Res Int, 2022. 2022: p. 3233004.

  29. Weiler M, et al. mTOR target NDRG1 confers MGMT-dependent resistance to alkylating chemotherapy. Proc Natl Acad Sci U S A. 2014;111(1):409–14.

    Article  CAS  PubMed  Google Scholar 

  30. Ren Y, et al. Spatial transcriptomics reveals niche-specific enrichment and vulnerabilities of radial glial stem-like cells in malignant gliomas. Nat Commun. 2023;14(1):1028.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Carlson JC, et al. Identification of diverse tumor endothelial cell populations in malignant glioma. Neuro Oncol. 2021;23(6):932–44.

    Article  CAS  PubMed  Google Scholar 

  32. Xie Y et al. Key molecular alterations in endothelial cells in human glioblastoma uncovered through single-cell RNA sequencing. JCI Insight, 2021. 6(15).

  33. Prager BC, et al. Glioblastoma Stem Cells: driving resilience through Chaos. Trends Cancer. 2020;6(3):223–35.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Kumar S, et al. Identification of vascular cues contributing to cancer cell stemness and function. Angiogenesis. 2022;25(3):355–71.

    Article  CAS  PubMed  Google Scholar 

  35. Jain RK, et al. Angiogenesis in brain tumours. Nat Rev Neurosci. 2007;8(8):610–22.

    Article  CAS  PubMed  Google Scholar 

  36. Li D, et al. Glioma-associated human endothelial cell-derived extracellular vesicles specifically promote the tumourigenicity of glioma stem cells via CD9. Oncogene. 2019;38(43):6898–912.

    Article  CAS  PubMed  Google Scholar 

  37. Sharma A, Shiras A. Cancer stem cell-vascular endothelial cell interactions in glioblastoma. Biochem Biophys Res Commun. 2016;473(3):688–92.

    Article  CAS  PubMed  Google Scholar 

  38. Adjei-Sowah EA, et al. Investigating the Interactions of Glioma Stem Cells in the Perivascular Niche at single-cell resolution using a microfluidic Tumor Microenvironment Model. Adv Sci (Weinh). 2022;9(21):e2201436.

    Article  PubMed  Google Scholar 

  39. Shi Y, et al. Tumour-associated macrophages secrete pleiotrophin to promote PTPRZ1 signalling in glioblastoma stem cells for tumour growth. Nat Commun. 2017;8:15080.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Yang M, et al. PTN-PTPRZ1 signaling axis blocking mediates tumor microenvironment remodeling for enhanced glioblastoma treatment. J Control Release. 2023;353:63–76.

    Article  CAS  PubMed  Google Scholar 

  41. Cheng L, et al. Glioblastoma stem cells generate vascular pericytes to support vessel function and tumor growth. Cell. 2013;153(1):139–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Han X, et al. P4HA1 regulates CD31 via COL6A1 in the transition of Glioblastoma Stem-Like cells to Tumor Endothelioid cells. Front Oncol. 2022;12:836511.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Mei X, et al. Glioblastoma stem cell differentiation into endothelial cells evidenced through live-cell imaging. Neuro Oncol. 2017;19(8):1109–18.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

This work would like to thank the clinicians and organizations that contributed to the samples used in this study.

Funding

This work was supported by CAS Project for Young Scientists in Basic Research (Grant YSBR-067), the Basic Research Pilot Program of Suzhou (SJC2022007), National Natural Science Foundation of China (82273336), Pre-Research and Construction Project of Major Scientific Research Facilities in Jiangsu Province (BM2022010), Suzhou Science and Technology Plan Project (Grant NO. SZS2022008) Science and Technology Service Network Initiative of Chinese Academy of Sciences Huangpu special funds (STS-HP-202304), Suzhou Medical Innovation Applied Research Project (SKYD2022090).

Author information

Author notes

  1. Jiukun Hu and Xiaohan Sa contributed equally to this work.

Authors and Affiliations

  1. School of Biomedical Engineering (Suzhou), Division of Life Sciences and Medicine, University of Science and Technology of China, Hefei, 230026, China

    Jiukun Hu, Xiaohan Sa, Yuwen Han, Minxuan Sun, Shaheryar Shafi & Yuanshuai Zhou

  2. Department of Biomaterials and Stem Cells, Suzhou Institute of Biomedical Engineering and Technology, Chinese Academy of Sciences, Keling Road No.88, Suzhou, 215163, China

    Jiukun Hu, Xiaohan Sa, Yue Yang, Yuwen Han, Minxuan Sun, Shaheryar Shafi & Yuanshuai Zhou

  3. Department of Chemistry, College of Sciences, Shanghai University, Shanghai, 200444, China

    Yue Yang

  4. Department of Neurosurgery, Suzhou Hospital, Affiliated Hospital of Medical School, Nanjing University, Lijiang Road No. 1, Suzhou, 215153, China

    Jie Wu

  5. Institute of Biomedical & Genetic Engineering, 24-Mauve Area G-9/1, Islamabad, 44000, Pakistan

    Nafees Ahmad

  6. Institute of Pharmaceutical Sciences, Khyber Medical University, F1 Phase-6 Rd, Phase 5 Hayatabad, Peshawar, Khyber Pakhtunkhwa, 25100, Pakistan

    Sami Siraj

  7. Institute of Clinical Medicine Research, Suzhou Hospital, Affiliated Hospital of Medical School, Nanjing University, Lijiang Road No. 1, Suzhou, 215153, China

    Jiao Yang

Authors
  1. Jiukun Hu
    View author publications

    You can also search for this author inPubMed Google Scholar

  2. Xiaohan Sa
    View author publications

    You can also search for this author inPubMed Google Scholar

  3. Yue Yang
    View author publications

    You can also search for this author inPubMed Google Scholar

  4. Yuwen Han
    View author publications

    You can also search for this author inPubMed Google Scholar

  5. Jie Wu
    View author publications

    You can also search for this author inPubMed Google Scholar

  6. Minxuan Sun
    View author publications

    You can also search for this author inPubMed Google Scholar

  7. Shaheryar Shafi
    View author publications

    You can also search for this author inPubMed Google Scholar

  8. Nafees Ahmad
    View author publications

    You can also search for this author inPubMed Google Scholar

  9. Sami Siraj
    View author publications

    You can also search for this author inPubMed Google Scholar

  10. Jiao Yang
    View author publications

    You can also search for this author inPubMed Google Scholar

  11. Yuanshuai Zhou
    View author publications

    You can also search for this author inPubMed Google Scholar

Contributions

Project concept and design: YZ and JY. Collection of RNA-seq data: JH, YH, XS and YZ. Computational analysis: JH and YZ. Data interpretation and critical design and revision: JH, YY, NA, MS and JY. Figures editing: JW, XS, SaS and JY. Original Draft: JH, XS, SS and JY. Supervision and project administration: JY and YZ. All authors have revised the manuscript, approved the submitted version, and agreed to be accountable for all aspects of the work.

Corresponding authors

Correspondence to Jiao Yang or Yuanshuai Zhou.

Ethics declarations

Ethics approval and consent to participate

All glioma samples were collected with signed informed consent according to the internal review and ethics boards of Wuhan Union Hospital. The procedure was approved by the ethics boards of Wuhan Union Hospital.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Electronic supplementary material

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Hu, J., Sa, X., Yang, Y. et al. Multi-transcriptomics reveals niche-specific expression programs and endothelial cells in glioblastoma. J Transl Med 23, 444 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12967-025-06185-z

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12967-025-06185-z

Keywords

Journal of Translational Medicine

ISSN: 1479-5876

Contact us