Bulk RNA-seq data were obtained from the Chinese Glioma Genome Atlas (CGGA, http://www.cgga.org.cn/), specifically dataset mRNAseq_325, which includes transcriptomic and clinical data from 325 glioma samples. For specific analyses, subsets of samples were used according to World Health Organization (WHO) grades: 103 Grade Ⅱ, 79 Grade Ⅲ, and 139 Grade Ⅳ cases [25]. Single-cell RNA-seq data were sourced from the Gene Expression Omnibus (GEO, https://www.ncbi.nlm.nih.gov/geo/) dataset GSE182109, from this cohort, we selected newly diagnosed glioblastoma (ndGBM) and low-grade glioma (LGG) samples for downstream analysis. It is noteworthy that a key commonality between the CGGA bulk dataset and this single-cell dataset is the inclusion of glioma samples across different grades, enabling multimodal exploration of our research questions [25, 26, 27]. Spatial transcriptomics data from Datadryad (https://doi.org/10.5061/dryad.h70rxwdmj) were accessed. Gene sets related to EMT and hypoxia were retrieved from the HALLMARK database (https://www.gsea-msigdb.org/gsea/msigdb/index.jsp), while those associated with PERK signaling were obtained from the REACTOME database (https://curator.reactome.org/) [28, 29, 30].

PERK, EMT, and hypoxia Gene sets from the HALLMARK and REACTOME databases were employed to transform the expression matrix of glioma samples from the CGGA dataset into enrichment score (ES) matrices using the R package in RStudio R 4.1.2 (R Core Team, Vienna, Austria). This process generated GSVA enrichment scores for each gene set across all samples, which were visualized using the ComplexHeatmap 2.13.1 (Guangchuang Yu, Hong Kong, China).

Clinical data and PERK scores from the CGGA database were used to perform survival analysis. Statistical analyses and survival data visualization were conducted using the R package survival in RStudio.

Differential gene expression analysis on transcriptome count data was conducted between two groups, with the reference group designated as “Low”. The DESeq2 1.43.0 (Michael Love, Berkeley, CA, USA) was used to identify differentially expressed genes on the raw count matrix, following standard procedures. Differential expression was visualized through volcano plots using ggplot2 3.4.2 (Hadley Wickham, New Zealand). Next, PPI network analysis was performed using the STRING database (https://string-db.org) on the top 100 upregulated genes identified in the high PERK score group from the CGGA database [31]. Finally, we utilized igraph for data processing, including data cleaning and organization, followed by visualization of the network graph by using ggraph.

We conducted GSEA on the differentially expressed genes, including Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), and Hallmark pathways, using gene sets from the MSigDB Collections [29]. Initially, molecular identifiers in the input data were converted to ensure compatibility with gene set databases. Subsequently, we employed the clusterProfiler 4.2.2 (Guangchuang Yu, Hong Kong, China) in R to perform GSEA. GSEA results were visualized to identify enriched pathways. In addition, we computed immune infiltration levels in the CGGA dataset using the ssGSEA Algorithm from the GSVA 1.42.0 (Robert K. Gideons, Barcelona, Catalonia, Spain), based on markers for 24 immune cell types [32].

We analyzed single-cell RNA sequencing data using the Seurat package in R 4.1.2 (R Core Team, Vienna, Austria) by RStudio 2021.09.1 Build 372 (Posit, Boston, MA, USA). Initially, genes expressed in fewer than three cells and cells expressing fewer than 100 genes were filtered out. Next, the RemoveDoublets function was applied to eliminate 4% of doublets, followed by exclusion of samples with mitochondrial gene expression exceeding 10%. Gene expression variability was determined using Seurat’s FindVariableFeatures function. Dimensionality reduction was performed with Principal Component Analysis (PCA), and batch effects were corrected using the RunFastMNN algorithm. UMAP clustering identified 26 distinct clusters; cell types were annotated by integrating reported and highly variable marker genes. Tumor cell subclusters were then extracted for further dimensionality reduction and clustering using t-Distributed Stochastic Neighbor Embedding (tSNE). The annotation of glioma cell subtypes was based on gene sets and signature genes provided by previously published literature [33].

We analyzed the single-cell transcriptome data by employing the AddMododuleScore method to assess the activity of specific gene expression modules at the cellular level. This method transforms gene expression data into numerical scores that reflect the involvement of each cell in the target biological processes. Results were visualized using ggplot2.

We utilized the to conduct Pseudo-temporal analysis of tumor cell populations using machine learning methods were performed using the monocle2 2.22.0 (Cole Trapnell, Seattle, WA, USA) [34]. Changes in the expression of key marker genes associated with the EMT- and hypoxia-related marker genes were observed along their transitions in tumor cell subpopulations.

Cell-cell interactions between PERK-activated glioma cells and their immune microenvironment were explored using the CellChat 1.6.1 (Suoqin Jin, Hong Kong, China) [35]. Intercellular communication was simulated interactions involving ligands, receptors, and their respective auxiliary factors were evaluated. Signaling pathways were visualized using the “netVisual_aggregate” function. Ligands and receptors were categorized as outgoing and incoming signals, respectively.

We utilized the Seurat package in RStudio (version 4.1.2), to read and assess spatial transcriptomics data. Data were normalized using SCTransform, ensuring robustness in downstream analyses. We assigned pathway scores to each spatial spot through AddMoudleScore, allowing the evaluation of pathway activities across spatial locations. Gene and pathway expression were visually represented using SpatialFeaturePlot in Seurat 4.3.0 (Satija Lab, NY, USA).

Single-cell RNA sequencing data were deconvoluted into spatial transcriptomics data using the SPOTlight method. Gene expression matrices for each cell were extracted from our single-cell RNA sequencing data. We integrated spatial information and gene expression patterns employing SPOTlight to map single-cell data onto a spatial grid and visualize it using spatial feature plots. Furthermore, the relationship between the proportion of glioma cells with high activation of PERK and the EMT score in each SPOT was investigated through correlation analysis.

The glioma cell lines (LN229) obtained from Procell Life Science & Technology Co., Ltd. (Wuhan, China). These cell lines were cultured in a 5% O₂ normal oxygen cell incubator and a 1% O₂ hypoxia cell incubator according to the experimental protocol. Next, the shPERK lentivirus was transduced into glioma cells, which were expanded to the required number, and a stably transfected PERK-low-expression glioma cell line was selected. The cells were randomly categorized into the negative control and the experimental groups. All cell lines were validated by STR profiling and tested negative for mycoplasma.

PERK and negative control (NC) shRNAs were first constructed according to the effective PERK shRNA (5^′-GCAGGUCAUUAGUAAUUAU-3^′) and NC sequences (5^′-UUCUCCGAACGUGUCACGU-3^′), which were validated in this study. For transient transfection experiments, PERK and negative control (NC) shRNAs were directly used in this study. For stable transfection, PERK and NC shRNAs were cloned into the pENTR/U6-GFP vector containing GFP, followed by recombination between the pENTR/U6-shRNA-GFP and the pLenti6/Block-it-DEST vector. The resultant vector and package plasmid were co-transfected into HEK 293T cells, and viral particles were collected after 48 h. LN229 cells were then infected with lentiPERK shRNA or NC vector, and cells expressing GFP were chosen for further culture.

Antibodies used for Western blotting were as follows: anti-PERK (Santa Cruz Biotechnology, Inc., cat. no. sc-13073, 1:100, Dallas, TX, USA), anti-p-PERK (Thr981) (Santa Cruz Biotechnology, Inc., cat. no. sc-32577, 1:100, Dallas, TX, USA), anti-$\beta$-ACTIN (Santa Cruz Biotechnology, Inc., cat. no. sc-8432, 1:100, Dallas, TX, USA), anti-E-cadherin (Abcam, cat. no. ab231303, 1:1000, Shanghai, China), anti-N-cadherin (Santa Cruz Biotechnology, Inc., cat. no. SC-53488, 1:200, Dallas, TX, USA), and anti-Vimentin (Abcam, cat. no. ab92547, 1:1000, Shanghai, China). The tissue homogenate or the harvested cells were lysed in protein lysis buffer containing PMSF on an ice bath, and centrifuged. The supernatant was used for protein quantification by the BCA method. A 10–100 µg of protein in each sample was separated by sodium dodecyl sulfate-polyacrylamide gel electrophoresis at 300 mA and transferred to the PVDF membrane. The membrane was placed in the blocking solution for 1 h at 37 ℃ and incubated overnight with the primary antibody at 4 ℃ in a refrigerator. The membrane was blocked for 1 h at 37 °C and then incubated overnight with the primary antibody at 4 °C. After washing with TBST, the membrane was subsequently incubated with a goat anti-mouse IgG horseradish peroxidase-conjugated secondary antibody (Beyotime Institute of Biotechnology, cat. no. A0216, 1:1000, Shanghai, China) for 1 h at room temperature. The membrane was rinsed with TBST solution, and an enhanced chemiluminescence reagent was added to detect signals using a versatile iBright CL1500 Imaging System (Thermo Fisher Scientific, Waltham, MA, USA), and semi-quantitative statistical analysis was performed using ImageLab 6.1 (Bio-Rad Laboratories, Hercules, CA, USA).

The cells were cultured in serum-free medium for 24 h, trypsinized, centrifuged at 1500 rpm for 3 min, and the supernatant was aspirated. The cells were resuspended in phosphate-buffered saline (PBS), the supernatant was discarded following centrifugation, resuspended again in serum-free medium containing 0.1% BSA to a density of 1 $\times$ 10⁵ cells/mL. To a 24-well plate, 500 µL of complete medium and 200 µL of cell suspension was added to the prefabricated Matrigel transwell chamber placed in the 24-well plate containing complete medium. The cells were cultured in a 24-well plate incubator for 24 h, and the cells in the upper chamber were removed. The membrane was then fixed in 10% methanol for 30 s and stained with crystal violet for 20 min. The cells were analysed by a microscope (Leica, DMi1, Shanghai, China) and counted after subtracting the background.

The cells were seeded on 6-well plates in triplicate. After overnight culture, the cells were transversely scratched in the middle of the well with a pre-sterilized pipette tip. The floating cells were washed with PBS and fresh medium was added. The cells were visualized under a microscope (Leica, DMi1, Shanghai, China) and images were taken every 24 h. Image-Pro Plus 6.0 (Media Cybernetics, Rockville, MD, USA) was used to calculate the scratch area before and after healing was calculated using the IPP software, and the wound closure rate between the different groups was compared. Wound closure rate using the wound contraction ratioformula, defined as (initial area – current area) / initial area $\times$ 100%.

Animal experiments were reviewed and approved by the Institutional Animal Care and Use Committee (IACUC) of Shanghai General Hospital (ID:2021SQ054). All experimental procedures involving animals were conducted adhering strictly with with the Guide for the Care and Use of Laboratory Animals (NIH Publication No. 86-23, revised 1985, U.S. National Academy of Sciences), with predefined humane endpoints (tumor volume $\le$1500 mm³, ulceration, or body weight loss $\ge$20%). Twelve female BALB/c nude mice (6-week-old, weighing 18–22 g) were obtained from Beijing Weitong Lihua Laboratory Animal Technology Co., Ltd., China, and were used for the tumor formation assay. These mice were housed under specific pathogen-free conditions in individually ventilated cages with controlled temperature (22 $\pm$ 1 ℃), humidity (50 $\pm$ 10%), a 12 h light/dark cycle, and ad libitum access to autoclaved food and water. After 1 week of acclimatization, 12 mice were randomly assigned into two groups (n = 6) and were injected subcutaneously with 2 $\times$ 10⁶ lenti-NC or lenti-PERK shRNA LN229 stably transfected glioma cells. Implantation in the hip region was preferred due to low vascularity and ease of monitoring. Tumour dimensions were measured using callipers every 3 days, and their volumes were calculated as 0.5 $\times$ length $\times$ width². Mice reaching endpoints were euthanized by cervical dislocation under isoflurane anesthesia; the remaining mice were sacrificed on day 30. The specific parameters were as follows: induction at 1.5–2% isoflurane in 100% oxygen at a flow rate of 2 L/min, followed by maintenance at 1–1.5% isoflurane. The tumours were weighed, photographed on a blue background, and then divided into two groups. One group was fixed using 4% paraformaldehyde, embedded in paraffin, and paraffin sections were prepared for immunohistochemistry. The second group of samples was frozen at –80 ℃ and used subsequently to extract total protein for western blot analysis.

Tissue samples were fixed with 4% paraformaldehyde, dehydrated in an alcohol gradient, embedded in xylene-treated paraffin, sliced, and heated at 60 ℃ in an oven for 2 h. The slices were dewaxed and hydrated. The antigens were retrieved by placing the slices in the antigen retrieval solution for 10 min under high pressure. After soaking in PBS, the endogenous peroxidase was blocked with 3% H₂O₂ for 15 min at room temperature. The samples were rinsed with PBS at room temperature and incubated with goat serum for 15 min at room temperature, followed by primary antibodies overnight at 4 °C. The specific antibodies and dilutions were as follows: anti-Vimentin (Abcam, cat. no. ab92547, 1:200, Shanghai, China), anti-E-cadherin (Abcam, cat. no. ab231303, 1:200, Shanghai, China), and anti-N-cadherin (Santa Cruz Biotechnology, Inc., cat. no. SC-53488, 1:50, Dallas, TX, USA). The prepared wash buffer was added and stored in the refrigerator at 4 ℃ overnight. The samples were then soaked in PBS, and the prepared horseradish-labelled sheep anti-mouse/rabbit IgG polymer was added and incubated at 37 ℃ for 30 min. Then, the DAB reagent and haematoxylin stain were applied for 3 min, and the slides were rinsed with 1% hydrochloric acid ethanol. The samples were separated, regressed to blue, dehydrated using gradient alcohol, and cleared with xylene. The stained sections were observed under a neutral resin seal microscope and analysed statistically.

Data are expressed as the mean $\pm$ standard error of the mean (SEM) from three different experiments. The significance was evaluated by Student’s t-test for comparisons between two groups, or one-way analysis of variance (ANOVA) followed by Tukey’s post-hoc testfor comparisons among multiple groups using GraphPad Prism 10.5.0 (GraphPad Software, San Diego, CA, USA). A p-value of less than 0.05 was considered statistically significant.

The activation status of hypoxia across different grades of gliomas was initially investigated by GSVA to score glioma samples from the CGGA dataset using hypoxia gene sets from HALLMARKER, and visualizing the results using the scatter plot (Fig. 1A). The hypoxia score in gliomas increases with tumor grade. Next, we divided the glioma samples into high- and low-hypoxia groups. Enrichment analysis of differentially expressed genes showed that the high hypoxia group was enriched in signaling pathways related to ERS and EMT, such as endoplasmic reticulum lumen, extracellular matrix structural constituent, and leukocyte-mediated immunity pathways (Fig. 1B). Further correlation analysis revealed that, compared with low-grade gliomas, the association between hypoxia and the PERK pathway was stronger in high-grade gliomas, suggesting a possible differential activation of the PERK pathway in high-grade gliomas than in low-grade gliomas under hypoxic conditions (Fig. 1C). We further stratified the cohort into high and low PERK groups based on GSVA scores, and observed lower overall survival (OS) in the high PERK group (Fig. 1D).

Fig. 1.

Bulk transcriptomic analysis showed a correlation between PKR-like endoplasmic reticulum kinase (PERK) pathway activation and epithelial-mesenchymal transition (EMT) in hypoxic gliomas. (A) Scatter plots of hypoxia Gene Set Variation Analysis (GSVA) scores for gliomas of different grades. (B) Enrichment analysis of differentially expressed genes in high- and low-hypoxia groups in glioma. (C) Scatter plot presenting the correlation between hypoxia and PERK score in gliomas of different grades. (D) Survival curve of high- and low-PERK groups in glioma. (E) Protein-Protein Interaction network analysis of the top 100 upregulated genes in the high-PERK activation group. (F,G) Gene Set Enrichment Analysis of differentially expressed genes between high- and low-PERK activation groups. (H) GSVA analysis on gliomas from the Chinese Glioma Genome Atlas (CGGA) dataset and a heatmap to represent the findings. ** p 0.01, *** p 0.001.

We conducted differential analysis between high and low PERK score groups to explore whether high activation of PERK in gliomas correlates with EMT (Supplementary Fig. 1A). Subsequently, we performed PPI network analysis on significantly upregulated genes in the high PERK activation group, revealing elevated expression of MMP family genes associated with EMT (Fig. 1E). Further GSEA indicated significant enrichment of EMT and hypoxia-related pathways in the high PERK group (Fig. 1F,G). ssGSEA analysis further revealed higher levels of macrophage infiltration in the high PERK group, indicating that PERK activation alters the immune microenvironment of gliomas (Supplementary Fig. 1B). Thus, these findings revealed that high-grade gliomas exhibit significantly higher activation in hypoxia, PERK, and EMT pathways than low-grade gliomas (Fig. 1H).

The cellular characteristics and interactions within glioma cells and their tumor microenvironment following PERK activation were examined employing publicly available single-cell sequencing data. We selected two cases each of low-grade and high-grade gliomas for subsequent dimensional reduction clustering to distinguish 11 cell subtypes (Fig. 2A,B, Supplementary Fig. 1C). GSVA revealed that compared with low-grade gliomas, high-grade gliomas exhibited significant activation of EMT, hypoxia, and UPR pathways (Fig. 2C). Additionally, AddmodleScore confirmed significantly higher activation of these pathways in high-grade gliomas than in low-grade gliomas, with a cluster of tumor cell subpopulations only in high-grade gliomas exhibiting a state of hypoxia, PERK, and EMT activation (Fig. 2D–G, Supplementary Fig. 1D). This finding is consistent with our previous observations.

Fig. 2.

PKR-like endoplasmic reticulum kinase (PERK), hypoxia and epithelial-mesenchymal transition (EMT) activation in different grades of glioma and cell subpopulations. (A,B) Uniform Manifold Approximation and Projection (UMAP) plots showing annotation of high- and low-grade glioma cells. (C) Gene Set Variation Analysis scoring of high- and low-grade gliomas based on the Hallmark database. (D) Violin plots of AddModuleScore for PERK, hypoxia, and EMT pathways in different grades of gliomas. (E–G) UMAP featureplot of PERK, hypoxia, and EMT pathways scores in gliomas. Black circles highlight tumor cell populations exhibiting high expression of the three specified pathways. **** p 0.0001.

Specific tumor subclusters were explored after isolating tumor cells and further subdividing them into AC_like, MES_like, NPC_like, OPC_like, and Prolifeing_Tumor subtypes by their markers. Then we computed the GSVA scores for each tumor subset based on gene sets derived from published literature. The subgroups exhibited the highest enrichment scores corresponding to their respective gene signatures. This consistency supports the validity and reliability of our cell type annotations (Supplementary Fig. 1E,F and Supplementary Table 1), mesenchymal-like (MES-like) subtype were predominantly present in high-grade gliomas (Fig. 3A,B). The MES-like subtype of tumor cells revealed the highest PERK scoring (Fig. 3C), based on which, we categorized glioma cells into high and low PERK activation groups. The proportion of high PERK activation tumor cells was significantly higher in high-grade gliomas (Fig. 3D). Subsequent GSEA analyses exhibited enhanced activation of EMT, hypoxia, and UPR pathways in the high activation group, further establishing the close association among hypoxia, PERK, and EMT (Fig. 3E,F). Next, a pseudo-temporal analysis was conducted to further investigate the relationship among PERK activation and EMT and hypoxia in glioma cells. The quasi-temporal heatmap simulates the developmental trajectory of the cell (Fig. 3G). At the cellular subtype level, the MES subtype predominantly appeared in the later stages of the pseudo-temporal trajectory, corresponding closely to the developmental branch associated with high PERK activation group (Fig. 3H, Supplementary Fig. 2A,B). We further examined the distribution of key EMT-related genes across high- and low-PERK activation groups and observed that the transcriptional expression levels of most hypoxia- and EMT-related genes increased with the transition to the high-PERK activation subgroup (Fig. 3I). These results suggest that the activation of the PERK pathway is accompanied by upregulation of EMT-related genes.

Fig. 3.

Single-cell analysis of tumor cell subpopulations shows a strong association between high-PKR-like endoplasmic reticulum kinase (PERK) activation and epithelial-mesenchymal transition (EMT). (A,B) Uniform Manifold Approximation and Projection plots of glioma cell subpopulations in high- and low-grade gliomas. (C) Box plots of AddModuleScore for the PERK pathway in different cell subpopulations. (D) t-Distributed Stochastic Neighbor Embedding (t-SNE) map of glioma cell distribution in high and low PERK pathways. (E) Gene Set Enrichment Analysis of differentially expressed genes between high- and low-PERK pathway glioma tumor cells. (F) A heatmap presenting key gene expression for selected EMT and hypoxia pathways in high- and low-PERK pathway groups. (G,H) A pseudotime heatmap of tumor cells in gliomas. (I) Changes in expression of selected EMT-related genes along the pseudotime developmental axis in gliomas. **** p 0.0001.

The role of PERK was further investigated to assess the invasion and migration ability of glioma cells under hypoxic conditions. We first evaluated the effect of PERK knockdown on LN229 cells under hypoxia. Western blot analysis indicated that the knockdown of PERK significantly inhibited the expression of phosphorylated PERK in glioma cells relative to control groups, whereas hypoxic conditions partially restored p-PERK levels (Fig. 4A). Our results further demonstrated a reduction in E-cadherin expression, but enhanced vimentin and N-cadherin expression in glioma cells exposed to hypoxic conditions (Fig. 4B). PERK knockdown also inhibited the invasion and migration ability of glioma cells under hypoxia (Fig. 4C,D). Thus, PERK silencing lowered the invasion and migration ability of in vitro glioma cells under hypoxia.

Fig. 4.

PKR-like endoplasmic reticulum kinase (PERK) silencing lowers the invasive and migratory ability of glioma cells under hypoxia. LN229 cells were transfected with negative control (NC) or PERK shRNA, and 48 h later, (A) the levels of PERK and phosphorylated (p)-PERK in glioma cells under hypoxia were evaluated by western blotting, (B) the protein levels of epithelial-mesenchymal transition (EMT) in glioma cells under hypoxia were evaluated by western blotting, (C) the migration of glioma cells under hypoxia was quantified using a wound closure rate. Scale bar: 200 µm, and (D) the number of invasive cells was examined by a transwell matrix penetration assay. Scale bar: 80 µm. *** p 0.001. All data are representative of three independent experiments.

The above findings indicate that PERK silencing reduced the invasion and migration ability of glioma cells, likely through modulation of EMT and ER stress-related pathways. However, our in vitro data need to be further validated in vivo. As glioma cells grow in a hypoxic microenvironment, we hypothesized that PERK silencing may also inhibit glioma invasion and migration in vivo, further influencing glioma growth. Consistently, the weights and volumes of the lenti-PERK shRNA tumours were significantly lower than those of control animals (Fig. 5A). Analysis of the expression of related tumor proteins revealed reduced p-PERK protein levels in the PERK-silenced gliomas than those in the control glioma cells (Fig. 5B). Furthermore, immunohistochemical analysis demonstrated decreased levels of vimentin and N-cadherin proteins, along with elevated levels of E-cadherin protein in PERK-silenced gliomas (Fig. 5C). These results suggest that PERK plays an important role in glioma invasion and migration and is associated with EMT and ER stress processes.

Fig. 5.

PKR-like endoplasmic reticulum kinase (PERK) silencing suppresses glioma growth and the epithelial-mesenchymal transition (EMT) process in vivo xenografts. (A) Lenti-PERK shRNA and lenti-NC LN229 cells were subcutaneously injected into BALB/c nude female mice. Tumour volume and tumour weights were measured using the in vivo proliferation assay. (B) Western blotting was performed for phosphorylated (p)-PERK and PERK in vivo. (C) Representative photomicrographs of immunohistochemical analyses for E-cadherin, N-cadherin, and vimentin in vivo. Scale bar: 50 µm. *** p 0.001. All data are representative of three independent experiments.

The potential role of PERK-high glioma cells in the tumor immune microenvironment was then explored through a CellChat analysis. Compared with the low PERK group, the high PERK group exhibited significantly enhanced signal outputs (Fig. 6A, Supplementary Fig. 2C,E), primarily concentrated in Macrophage migration inhibitory factor (MIF), SPP1, midkine (MK), and VEGF signaling pathways (Supplementary Fig. 2D). SPP1 signaling analysis revealed more active interactions between tumor cells with high PERK activity and macrophages than in the low PERK group (Fig. 6B, Supplementary Fig. 2F). Using publicly available spatial transcriptomics datasets, we analyzed two cases of low-grade and two of high-grade glioma samples. Hypoxic regions in the glioma were identified using the AddModuleScore, followed by differential gene expression and enrichment analyses (Supplementary Fig. 2G). These hypoxic regions were enriched for focal adhesion, endoplasmic reticulum lumen, and undold protein binding pathways, consistent with our bulk RNA-seq findings (Fig. 6C). Using SPOTlight, we deconvoluted single-cell data into spatial transcriptomics data and found that tumor cells with high PERK activation were more abundant in high-grade gliomas, whereas low PERK macrophages were more abundant in low-grade gliomas, even in the same GBM tissue, the tumor cells with high PERK expression were closer to the tumor core area, while the tumor cells with low PERK expression were located at the tumor border (Fig. 6D, Supplementary Fig. 2H). Spatial feature mapping of the EMT pathway using AddModuleScore and subsequent correlation analysis revealed that regions with denser distribution of high PERK activation glioma cells exhibited elevated EMT (Fig. 6E,F, cor = 0.6). Finally, we examined the spatial correlation between PERK-HIGH gliomas and TAMs and found a strong correlation between them (Fig. 6G, Supplementary Fig. 2I). Furthermore, regions enriched in PERK-high glioma cells also exhibited high expression of SPP1 and CD44, suggesting spatial consistency of this signaling axis (Fig. 6H).

Fig. 6.

Spatial transcriptomics and cell communication analyses suggest that tumor cells with high PKR-like endoplasmic reticulum kinase (PERK) activation may induce epithelial-mesenchymal transition (EMT) by enhancing secreted phosphoprotein 1 (SPP1)-CD44 signaling. (A) Heatmap of cell communication signaling via SPP1 in gliomas. (B) A Chordal diagram of SPP1 signal communication between glioma cells. (C) Enrichment analysis of hypoxic niche and other regional differential genes in the spatial transcriptome. (D) A feature map of high PERK activation glioma cells after SPOTLight deconvolution analysis. (E) Spatial feature map of EMT in glioma tissues. (F) Scatter plot showing the spatial correlation between high PERK activation glioma cells and the EMT pathway. (G) A heatmap of spatial cell communication correlations in glioma tissues. (H) Spatial distribution feature map of SPP1 and CD44 in glioma tissues. **** p 0.0001.