Skip to main content

Integrated transcriptomics uncovers an enhanced association between the prion protein gene expression and vesicle dynamics signatures in glioblastomas

Abstract

Background

Glioblastoma (GBM) is an aggressive brain tumor that exhibits resistance to current treatment, making the identification of novel therapeutic targets essential. In this context, cellular prion protein (PrPC) stands out as a potential candidate for new therapies. Encoded by the PRNP gene, PrPC can present increased expression levels in GBM, impacting cell proliferation, growth, migration, invasion and stemness. Nevertheless, the exact molecular mechanisms through which PRNP/PrPC modulates key aspects of GBM biology remain elusive.

Methods

To elucidate the implications of PRNP/PrPC in the biology of this cancer, we analyzed publicly available RNA sequencing (RNA-seq) data of patient-derived GBMs from four independent studies. First, we ranked samples profiled by bulk RNA-seq as PRNPhigh and PRNPlow and compared their transcriptomic landscape. Then, we analyzed PRNP+ and PRNP- GBM cells profiled by single-cell RNA-seq to further understand the molecular context within which PRNP/PrPC might function in this tumor. We explored an additional proteomics dataset, applying similar comparative approaches, to corroborate our findings.

Results

Functional profiling revealed that vesicular dynamics signatures are strongly correlated with PRNP/PrPC levels in GBM. We found a panel of 73 genes, enriched in vesicle-related pathways, whose expression levels are increased in PRNPhigh/PRNP+ cells across all RNA-seq datasets. Vesicle-associated genes, ANXA1, RAB31, DSTN and SYPL1, were found to be upregulated in vitro in an in-house collection of patient-derived GBM. Moreover, proteome analysis of patient-derived samples reinforces the findings of enhanced vesicle biogenesis, processing and trafficking in PRNPhigh/PRNP+ GBM cells.

Conclusions

Together, our findings shed light on a novel role for PrPC as a potential modulator of vesicle biology in GBM, which is pivotal for intercellular communication and cancer maintenance. We also introduce GBMdiscovery, a novel user-friendly tool that allows the investigation of specific genes in GBM biology.

Peer Review reports

Introduction

Glioblastoma (GBM) is an aggressive central nervous system tumor classified as a grade IV glioma by the World Health Organization (WHO) [1]. At the molecular level, the mutation status of isocitrate dehydrogenase (IDH) is an important diagnostic and prognostic biomarker. Survival expectancy for patients bearing IDH wild-type (IDHwt) GBM is dismal, and this type of cancer is more aggressive than the IDH-mutant counterpart. Given the differences between IDHwt and IDH-mutant tumors, a new classification defining GBM exclusively as IDHwt cancer has been proposed [2]. Transcriptionally, GBM is stratified into three molecular subtypes: proneural, classical, and mesenchymal [3]. Recently, cells resembling different molecular subtypes and spanning distinct genetic programs were shown to exist within a single tumor, highlighting the heterogeneous nature of this cancer [4]. In this context, GBM bears a subpopulation of glioma stem-like cells (GSCs) associated with its heterogeneity, aggressiveness, and therapy resistance [5].

In the past 20 years, standard treatment delivered to individuals with GBM has been surgical resection followed by radiotherapy and chemotherapy with temozolomide, providing patients with a median overall survival of 14-20 months [6,7,8]. Several aspects are proposed to be associated with GBM malignancy. Low therapy effectiveness is partially explained by resistance mechanisms to both radio and chemotherapy [9], and surgical resection is hindered by the invasive nature of GBM [7]. GSCs may play a role in therapy resistance by releasing extracellular vesicles (EVs), which carry bioactive molecules (protein, RNA, DNA, sugar, and lipids) able to manipulate the tumor microenvironment to support tumor growth and progression [10, 11].

Research efforts have focused on identifying novel potential molecular targets to be exploited in GBM therapeutics and to give patients better life expectancy and quality. Increasing evidence has demonstrated that the cellular prion protein (PrPC) plays a key role in GBM biology. PrPC is a glycosylphosphatidylinositol-anchored cell surface glycoprotein that faces the extracellular space and is encoded by the PRNP gene, mediating essential processes in mammalian nervous system [12,13,14]. This protein is a scaffolding molecule and orchestrates signaling complexes on the plasma membrane, interacting with a plethora of ligands [15]. Beyond its important physiological functions [16,17,18], PrPC is also crucial in different cancer types [19, 20], including brain tumors. PrPC is often found upregulated in GBM patient samples, and its silencing leads to diminished tumor growth and better survival expectancy in mice [21]. Moreover, PrPC modulates the expression of stemness and differentiation markers in GSCs, as well as cellular migration, proliferation, self-renewal, and tumorigenicity [22, 23]. Despite these recent advances, a deeper comprehension of the molecules and signaling pathways modulated by PrPC in GBM biology is still lacking and might elucidate novel therapeutic approaches against this glioma.

Herein, we analyzed four publicly available RNA sequencing (RNA-seq) datasets generated from patient-derived GBM at both bulk and single-cell (scRNA-seq) levels to unravel PrPC functions in the biology of this tumor. A fifth public dataset consisting of proteomics data also from patient-derived GBM confirmed and strengthened our findings at the protein level. We also developed GBMdiscovery, an R-based user-friendly software application that will help researchers understand the impact of their genes of interest in GBM biology through the analysis of all RNA-seq datasets investigated in this work.

Results

Pathways associated with membrane-enclosed organelles and secretion are enriched in GBM cells with high PRNP expression

To elucidate the repertoire of biological phenomena that might be modulated by PRNP/PrPC in GBM, we analyzed bulk RNA sequencing (RNA-seq) data of patient-derived primary GBM samples (n=157) from The Cancer Genome Atlas (TCGA) database (Fig. 1A). We initially verified that PRNP expression levels were increased in IDHwt GBM relative to IDH-mutant (Fig. 1B, left). Regarding molecular subtypes, PRNP presented similar expression levels between classical and mesenchymal GBM, and lower expression in the proneural subtype (Fig. 1B, right). In accordance with recent GBM classification guidelines, IDH-mutant (n=11) and samples of unknown subtype (n=22) or IDH status (n=4) were excluded prior to running any downstream analyses. Some unclassified samples in GBM subtype overlapped with unclassified IDH statuses, leading to 124 samples after filtering. We then normalized the count data of the remaining primary tumor samples according to log10(CPM+1) and ranked them by PRNP expression. Selecting GBM samples below the lower quartile to be the PRNPlow group (n=31), while samples above the upper quartile composed the PRNPhigh group (n=31) (Fig. 1C). Next, we assessed the composition of each group regarding GBM molecular subtypes (Fig. 1D). We found that both groups presented a similar composition of classical samples, while PRNPhigh consisted mainly of the mesenchymal subtype, with a decreased proportion of proneural samples relative to PRNPlow (Fig. 1D).

Fig. 1
figure 1

Impact of PRNP expression at bulk resolution in GBM samples from TCGA. A Schematic workflow of bulk RNA-seq data analyses of patient-derived primary GBM samples from TCGA (n=157). B Violin plots of PRNP normalized expression, according to log10 (counts per million [CPM]+1), in samples classified as IDH-mutant (n=11), IDHwt (n=142), or unclassified (n=4) (left, p=0.00016); and samples classified as classical (n=50), mesenchymal (n=67), proneural (n=18), or unclassified (n=22) (right, p=0.01). Kruskal-Wallis test. C GBM samples (n=124, after filtering) were ranked according to quartiles of PRNP expression. Those below the lower (PRNPlow, n=31) and above the upper quartile (PRNPhigh, n=31) were selected, as shown in the density plot of PRNP expression. D GBM molecular subtype composition of the PRNPlow and PRNPhigh groups. E Volcano plot of upregulated (red) and downregulated (blue) transcripts in PRNPhigh relative to PRNPlow. (The PRNP was removed from the plot). F Gene set enrichment analysis (GSEA) show enriched terms (Gene Ontology, GO) in the upregulated and downregulated transcripts of PRNPhigh GBM, sorted by normalized enrichment score (NES) and colored by adjusted p-value

We used raw counts as input in DESeq2 [24] to identify differentially expressed transcripts (DETs) between the groups, using PRNPlow as control. In total, 2634 significant DETs were identified, of which 1749 were upregulated and 885 downregulated (Fig. 1E, Table S1). We subjected DETs to gene set enrichment analysis (GSEA) and observed that several terms (Gene Ontology, GO) associated with cell cycle and DNA damage response (including negative regulation of cell cycle process, cell cycle checkpoint signaling, cell division, double-strand break repair, and DNA repair) were downregulated in PRNPhigh (Fig. 1F, Table S2). On the other hand, terms related to different aspects of membrane-enclosed organelles (including inner dynein arm assembly, microtube bundle formation, external encapsulating structure), plasma membrane events (cell projection, cell periphery, cell projection assembly, plasma membrane-bound cell projection assembly, plasma membrane) and secretion (extracellular region, extracellular matrix, collagen-containing extracellular matrix, extracellular space, extracellular matrix structural constituent) were upregulated in PRNPhigh (Fig. 1F, Table S2). Additionally, over-representation analysis (ORA) of upregulated transcripts showed enrichment for extracellular transport and the trans-Golgi network (Table S3). This caught our attention because, in addition to being a critical scaffolding protein on the plasma membrane [25], PrPC has been shown to be involved in intracellular trafficking and exosome biogenesis [26]. Moreover, many processes associated with the immune response (immunoglobulin complex, immune response, antigen receptor-mediated signaling pathway, among other terms) were substantially enriched in the PRNPhigh group, in line with results that show an increase in the mesenchymal molecular subtype in these samples [27]. Overall, the initial bulk RNA-seq data assessment suggests that high PRNP expression might be associated with intra- and extracellular transport pathways.

PRNP-positive cells are associated with vesicle dynamics in patient-derived GBM at single-cell level

While bulk RNA-seq can provide remarkable insights to biological questions, it has the limitation that information on individual cells is lost, and only the average transcriptional profile of the whole sample is obtained [28]. Our observations from bulk transcriptomics prompted us to ascertain the extent of our findings at the single-cell level, exploring three publicly available patient-derived GBM single-cell RNA-seq (scRNA-seq) datasets from independent studies [4, 29, 30] (Fig. 2A). For each dataset, we first isolated the malignant from the non-neoplastic cells based on each study’s metadata, filtering out IDH-mutant, recurrent, pediatric, and unclassified samples. Then, we visualized how PRNP expression was distributed across the tumors (Figs. 2B and S1A and B), finding that PRNP is widely expressed in cancer cells, and its levels are heterogeneous among different samples.

Fig. 2
figure 2

Integrated analysis of single-cell GBM datasets shows the impact of PRNP expression. A Schematic workflow of single-cell RNA-seq (scRNA-seq) data analyses carried out on three independent and publicly available datasets. B UMAP representation of PRNP expression in GBM cells from Darmanis et al., Neftel et al., and Richards et al. (C) Venn diagram of common marker genes of PRNP+ cells in all scRNA-seq datasets. D Functional profiling (ORA, GO) of the 840 common genes found in (C), ranked by -log10 (Adjusted p-value)

Using a similar approach as in our bulk RNA-seq analysis, we classified the single-cells as PRNP=0 (PRNP-) or PRNP>0 (PRNP+). Marker genes of PRNP+ cells were identified and subjected to ORA for each dataset, resulting in several over-represented vesicle-related terms (Fig. S2). Next, we compared common PRNP+ marker genes among the studies. Our results show that 840 genes lay at the intersection of the three datasets (Fig. 2C, Table S4), and, once again, numerous vesicle-associated terms were enhanced among the intersecting genes (Fig. 2D, Table S5). In addition, endoplasmic reticulum (ER)- and Golgi-associated terms, structures intrinsically involved in intracellular traffic, were also enriched (Fig. 2D, Table S5). Altogether, our findings from single-cell transcriptomics provide an in-depth insight into the association of PRNP and vesicle biology in GBM.

PRNP positively correlates with vesicle-associated genes in GBM

Subsequently, we aimed to understand how PRNP expression correlates with the levels of other genes expressed by cancer cells (Fig. 3A). We calculated Pearson correlation between PRNP, and all the genes identified in the sequencing of each study (Fig. 3B), and compared the positively correlated genes between the datasets, finding 2541 common genes (Fig. 3C, Table S6). Besides vesicle biogenesis, processing, intra- and intercellular trafficking, secretion, budding, fission and fusion, a plethora of components of the endocytic pathway are also enriched among the common positively correlated genes (Fig. 3D, Table S7). These data reinforce the putative function of PRNP/PrPC in the interplay between membrane-bound compartments and cellular trafficking.

Fig. 3
figure 3

PRNP positively correlates with vesicle-associated genes at single-cell resolution. A Schematic workflow of the correlation analysis between PRNP expression and all genes identified in the scRNA-seq of the three independent patient-derived GBM datasets. B Pearson correlation shows significant positively (red) and negatively (blue) correlated genes with PRNP for the three scRNA-seq datasets. C Venn diagram of common positively correlated genes identified in (B). D ORA analysis (gProfiler2 and GO) of genes positively correlated with PRNP

To strengthen our conclusions, we combined the findings from the analyzed bulk- and single-cell transcriptome data. We identified which genes were markers of PRNP+ GBM cells at the intersection among all scRNA-seq datasets and were also upregulated in the PRNPhigh samples from TCGA bulk RNA-seq data. Our approach yielded a panel of 73 genes, one of which was PRNP (Fig. 4A, Table S8). In general, the genes from our panel displayed a higher expression in tumor samples (primary and recurrent) compared to non-neoplastic tissues (Fig. S3), and increased levels in IDHwt GBM compared to IDH-mutant (Figure S4). Among molecular subtypes, the genes from our 73-gene panel had a more heterogeneous expression level but presented a prominent enhancement in the mesenchymal subtype (Figure S5). Confirming the observations from previous analyses, functional profiling of this gene panel demonstrated that vesicle biology terms (vesicle, extracellular vesicle, extracellular exosomes, intracellular vesicle, among others), were enriched (Fig. 4B, Table S9). These findings further corroborate the association between PRNP expression and vesicle dynamics in GBM biology, highlighting a set of genes whose expression could be affected by PRNP levels in the tumor.

Fig. 4
figure 4

Identification of a 73-gene panel enriched in traffic-related structures and vesicle dynamics in PRNPhigh samples and PRNP+ GBM cells. A Venn diagram of the 73 common genes with increased expression levels in PRNPhigh GBM cells (TCGA bulk RNA-seq) and upregulated in PRNP+ cells (Darmanis et al., Neftel et al. and Richards et al.). B Over-representation analysis (gProfiler2 - GO, KEGG, Reactome, and WikiPathways) of common genes found in (A). C Schematic workflow of the analysis of an in-house cohort of GBM patient-derived samples. D Relative expression of PRNP, SYPL1, DSTN, ANXA1 and RAB31 in patient-derived GBM samples expressing Low (n=26) or High (n=26) PRNP levels by RT-qPCR and 2-ΔΔCt, normalized by HPRT, GUSB and TBP. Student t-test, ****p<0.0001

As a proof-of-concept, we used an in-house collection of patient-derived GBM samples to validate the expression levels of a set of genes from our 73-gene panel. For this selection, we searched the literature for previous experimental evidence indicating that the genes from our panel were related to intracellular traffic, endocytic or exocytic vesicles, and cytoskeleton maintenance, which is also relevant for organelle transport. Besides PRNP, we selected synaptophysin Like 1 (SYPL1) [31] and destrin (DSTN) [32] for gene expression quantification. Additionally, we also selected annexin A1 (ANXA1) [33] and Ras-related protein Rab-31 (RAB31), genes that are not shown in our panel but are noteworthy vesicle components and regulators of transport vesicles, including in GBM [34, 35]. We performed RT-qPCR on 104 patient-derived GBM, dividing them according to their PRNP expression. Briefly, we used the Δ-Cq values for normalized PRNP expression and stratified the samples into quartiles according to these values, using the samples above the upper quartile and below the lower quartile to define, respectively, high (n=26) and low (n=26) PRNP levels (Fig. 4C). Expression of SYPL1, DSTN, ANXA1 and RAB31 was increased in samples with high PRNP expression compared to samples with low PRNP expression (Fig. 4D). Therefore, our experimental results using an additional cohort of patient-derived samples corroborate our findings from RNA-seq analyses.

Patient-derived whole-proteome analysis shows that PrPC expression is associated with vesicle dynamics

To investigate whether our results would be corroborated beyond the transcriptional level, we analyzed a publicly available patient-derived GBM proteomics dataset [36]. We divided tumor samples into groups with either high (PrPC-high; n=13) or low PrPC (PrPC-low; n=13) levels (Figure S6A) and identified differentially expressed proteins (DEPs) between these conditions, using PrPC-low as control (Figure S6B, Table S10). Protein classification shows membrane traffic proteins and cytoskeleton proteins among the upregulated DEPs (Figure S6C). Functional profiling revealed that several vesicle dynamics and traffic terms were enriched in the upregulated proteins in PrPC-high tumors (Figure S6D), which was in accordance with our previous transcriptomics data. Together, our findings advocate for an association of PrPC with vesicle phenomena not only at transcriptional but also at protein levels.

Vesicle dynamics signatures are positively correlated with PRNP expression in other patient-derived solid tumors

To evaluate if our findings extended to other cancers, we analyzed the impact of PRNP expression in other tumor types. We calculated the correlation between PRNP expression levels and vesicle dynamics signatures in different stages of patient-derived solid tumors from TCGA from: adrenal gland, bladder, breast, gallbladder, colon, head and neck, kidney, lung, ovary, pancreas, prostate, rectum, skin, stomach, esophagus, thyroid, and uterus (Figure S7, Table S11).

Our results showed that, for most cancer types, PRNP displayed significantly positive correlations with several vesicle dynamics signatures in at least one tumor stage. In particular, bladder, head and neck, kidney, and lung tumors were the cancer types with the highest number of positive correlations between PRNP expression and vesicle dynamics signatures among distinct tumor stages. These findings demonstrate, therefore, that PRNP expression is correlated with vesicle-related events in solid tumors from different sites, suggesting a conserved role of PRNP modulating these processes in cancer.

Increased vesicle dynamics signatures associated with shortened overall survival of glioblastoma patients

To understand the clinical relevance of our findings, our next step was to explore the association between PRNP expression and vesicle dynamics processes on patients’ overall survival, using TCGA-GBM survival data. We observed that a high PRNP expression is associated with a worse prognosis of GBM patients (Fig. 5A). Furthermore, we selected representative trafficking and vesicle dynamics gene sets from GO, based on our previous results. Subsequentially, we calculated gene signatures in IDHwt primary GBMs from TCGA, determining the optimal cutoff to separate samples with either low or high signatures, and performing Kaplan-Meier survival estimates (Fig. 5B).

Fig. 5
figure 5

Impact of PRNP and vesicle dynamics signatures on patients’ overall survival. A Kaplan-Meier curves demonstrate the overall survival of GBM patients either below or above the optimal cutoff, calculated for PRNP expression levels. A risk table is shown at the left. Statistical significance was assessed using a log-rank test, and p-value is stated in the graph. B Kaplan-Meier curves demonstrate the overall survival of GBM patients either below or above the optimal cutoff, calculated for gene signatures of traffic- or vesicle-related processes (GO). Statistical significance was assessed using log-rank tests, and p-values are stated in the graphs. Inserts show PRNP normalized expression in samples above and under the cutoff of each gene signature

The results demonstrate that enrichment for transport vesicle, exocytosis, exocytic process, vesicle docking, intercellular transport, vesicle lumen, membrane invagination, regulation of vesicle fusion and secretion is linked to shorter patients’ overall survival (Fig. 5B). Concomitantly with the survival curves, PRNP expression was evaluated in samples above and below cutoff for each signature and was significantly increased in exocytic process (Fig. 5B, inserts).

These results highlight the importance of our findings since PrPC mRNA and protein levels are associated with biological phenomena that are likely to be responsible for poor survival expectancy for patients with GBM.

GBMdiscovery as a novel platform to unravel the impact of specific genes in glioblastoma biology

We found a potential implication of PrPC in GBM biology, demonstrating that its mRNA and protein levels are associated with intracellular traffic and, more specifically, vesicle dynamics in this tumor. Since transcriptomics provided such interesting biological insights, we developed GBMdiscovery, an R-based shiny app that allows users to discover the implications of their own genes of interest in GBM biology (Fig. 6). The app is user-friendly, allowing researchers from any area with no bioinformatics background to use it, and it can be downloaded from https://github.com/marilenehohmuth/GlioblastomaDiscovery.

Fig. 6
figure 6

Leading page and example analysis of GBMdiscovery. GBMdiscovery is an R-based shiny app that allows users to discover the implications of their genes of interest in GBM biology, combining the four publicly available patient-derived bulk- and scRNA-seq datasets analyzed in the present study

GBMdiscovery analyzes all publicly available RNA-seq datasets presented in this study to find which genes and biological phenomena might be differentially regulated upon distinct levels of a specific gene. First, the app analyzes each dataset individually, displaying both genes and processes that are altered in each study, based on differential levels of the gene-of-interest of the user. GBMdiscovery then compares the results of all datasets to provide a list that contains the genes that are commonly down or upregulated upon distinct levels of the gene-of-interest in a ubiquitous manner across the four RNA-seq datasets. In addition, our app also calculates the correlation between this gene and all genes expressed by malignant cells identified in the studies. Therefore, GBMdiscovery aims to facilitate uncovering potentially interesting genes that may contribute to GBM biology, ultimately leading to the discovery of novel therapeutic targets against the tumor.

Discussion

Intra- and intercellular trafficking are broad terms comprising diverse biological phenomena such as endocytosis, protein sorting, vesicle transport, and secretion, which play a crucial role in cell signaling and communication. While critical for physiological functions and development [37], these mechanisms are hijacked by neoplastic cells to promote tumor survival, progression, and spreading [38, 39]. These processes are known to contribute to GBM, the most aggressive form of brain tumor [40, 41]. Therefore, identification of molecules that orchestrate such pathways is imperative for developing new treatment strategies against this glioma. Herein, we described for the first time that intracellular traffic processes, particularly vesicle dynamics, are upregulated in GBM patient-derived samples/cells with high PRNP/PrPC expression (Fig. 7). Moreover, we described a positive correlation between PRNP expression and the expression of genes associated with cellular trafficking and vesicle-related processes, proposing a putative role for this molecule in the control of GBM vesicle biogenesis, processing, and intra/inter-cellular transport.

Fig. 7
figure 7

Schematic diagram of the mechanisms involved in extra- and intracellular trafficking by which PrPC may affect GBM biology. PrPC is transiently found in the ER and Golgi during its biogenesis and is internalized by a vesicle-mediated process to be recycled or degraded. Vesicular biogenesis, endocytic path and exocytosis are processes enriched in GBM cells with high PRNP expression

Growing evidence has demonstrated that intracellular transport is essential for cancer metabolic reprogramming, modulating the energetic balance of cancer cells and their proliferative and invasive behaviors [39]. Regarding cell-cell interactions, exosomes have been implicated in exchanging a range of different cargoes between cancer cells and neoplastic and stromal cells [42]. The signal molecules carried from one cell to another ultimately lead to the regulation of pivotal mechanisms for tumor maintenance [38]. Exosomes have broad implications for cancer progression, as they contribute to forming pre-metastatic and metastatic niches [38, 42]. In GBM, EVs secreted by the tumor are enriched in pro-angiogenic factors. EVs derived from hypoxic GBM cells can promote tumor growth, vascularization, and proliferation, in addition to acidifying the tumor microenvironment, which leads to augmented EV trafficking both inwards and outwards the cells [43]. Small EVs released by GBM cells may be used for cancer subtyping since their cargo reflects the tumor phenotype and has been associated with enhanced glioma growth, proliferation, survival, migration, invasion, metabolic reprogramming, drug resistance, angiogenesis, and immunomodulation of the microenvironment [44]. In gliomas, microvesicles cargo can enhance tumor growth [45] transport molecules that modify gene expression in recipient cells [46], and trigger the proliferation of endothelial cells [47]. In the current study, we propose PrPC, a protein that regulates tumorigenesis and tumor growth [15, 21, 23], as an orchestrator of vesicle-related pathways in GBM.

In physiological conditions, previous studies have demonstrated that PrPC regulates exosome secretion in astrocytes and fibroblasts [26] and participates in synaptic vesicle release [48]. PrPC was implicated in the secretion of soluble factors by astrocytes and in their communication with neurons and extracellular matrix organization [49]. Evidence suggests this glycoprotein is sorted into exosomes and enriched in EVs [50]. Moreover, PrPC-expressing exosomes are tumorigenic and highly secreted by drug-resistant colorectal cancer cells upon hypoxic conditions in the tumor microenvironment [51]. The cellular location of PrPC is dynamic. While it is mainly located on the cell surface, PrPC is also found in other cellular compartments, such as ER and the Golgi apparatus, which are directly involved in intracellular traffic pathways [52, 53]. Additionally, this protein is observed in intracellular vesicles in different stages of endocytosis after being internalized [52]. Therefore, throughout its dynamic trafficking within the cell and leveraging its scaffolding features, PrPC may interact with orchestrators of intracellular traffic and modulate vesicle dynamics mechanisms (Fig. 7). PRNP/PrPC may also affect the expression of downstream genes and proteins involved in intracellular traffic phenomena, ultimately regulating vesicle-related processes, but this remains to be further investigated in vitro.

In this context, our work identified several genes whose expression is positively correlated with PRNP in GBM. From our 73-gene panel (PRNP + 72 genes), we assessed the expression of PRNP, SYPL1 and DSTN, on a collection of patient-derived GBMs. SYPL1 encodes synaptophysin-like 1 (SYPL1), a protein that has emerged as a prognostic marker and potential target in pancreatic ductal adenocarcinoma [31], hepatocellular carcinoma [54], and colorectal cancer [55]. Additionally, SYPL1 was described in the transport of vesicles [31]. DSTN, in turn, encodes destrin, a regulator of actin depolymerization [32]. The function of destrin in actin cytoskeleton maintenance is relevant from the vesicular biology perspective since the cytoskeleton plays key roles in vesicle transport [56, 57]. We also assessed the expression of ANXA1 and RAB31, since these are relevant molecules already involved in vesicle phenomena in GBM. Specifically, ANXA1 encodes the membrane protein annexin A1 (ANXA1), which is involved in tethering and EV aggregation [33]. ANXA1 is a marker gene of mesenchymal states in GBM [4], which is consistent with, as we showed, PRNPhigh/PRNP+ cells being mainly classified as mesenchymal. ANXA1 is overexpressed in high grade GBM samples and is a pivotal molecule to regulate GBM cell proliferation, migration, and invasion [34]. Moreover, increased expression level of ANXA1 in gliomas was associated with worse prognosis [58]. The protein encoded by RAB31 acts both as a driver of intraluminal vesicle formation and as a suppressor of multivesicular body degradation, ultimately regulating the exosome pathway independent of the endosomal sorting complex required for transport [59]. Interestingly, RAB31 knockdown significantly affected the capacity of GBM cell line growth in animal models [60]. Furthermore, a recent study reports that RAB31 silencing attenuates glioma invasion, process mediated by EVs from glioma-derived endothelial cells [35].

Interestingly, we stablished that our findings present clinical relevance, since we showed worse prognosis in patient-derived samples with both enhanced vesicle biology signatures and high PRNP/PrPC expression. Moreover, as reported here, the association between PRNP and vesicular dynamics has a potential impact in the biology of several types of solid tumors, strengthening this pathway as putative novel target for cancer therapeutic strategies.

Altogether, our findings provide valuable resources to further clarify molecular mechanisms governed by PrPC in tumor biology and open new paths for future GBM research. Of great interest, we make available a novel user-friendly tool that will expand our knowledge on this glioma, allowing users to evaluate the impact of their genes of interest in GBM biology in a robust manner across various RNA-seq datasets at both bulk and single-cell resolutions.

Methods

Bulk RNA sequencing data analysis

Data retrieval

GBM bulk RNA-seq data and corresponding metadata were retrieved from TCGA through TCGAbiolinks (v2.30.0) in R (v4.3.2). We downloaded raw count data of 175 samples of the TCGA-GBM project, of which 157 corresponded to primary tumor samples, 13 to recurrent samples, and 5 to solid normal tissue. To assess PRNP expression across different IDH statuses and transcriptional subtypes, the raw counts of all primary tumors were normalized according to counts per million (CPM) using edgeR (v3.42.4) [61] and converted to the log10 (CPM+1) scale. After that, IDH-mutant samples, and samples whose IDH status or transcriptional subtype was unknown were removed from the cohort; 124 samples remained in the dataset.

Establishment of groups with distinct PRNP expression levels

Following the same procedure described above, we obtained log10 (CPM+1)-normalized counts for all 124 primary IDHwt GBMs. Thereafter, we found the quartiles of PRNP expression across the samples, classifying those below the lower quartile – which displayed decreased PRNP levels – as PRNPlow (n=31) and the ones above the upper quartile – with increased PRNP expression – as PRNPhigh (n=31).

Identification of differentially expressed transcripts between PRNP high and PRNP low samples

After obtaining groups with different PRNP expression levels, we used the raw count data of the PRNPhigh and PRNPlow groups as input in DESeq2 (v1.40.2) [24] for differential expression analysis. PRNPlow was used as the control group. Log fold change shrinkage was performed with the shrinkage estimators of a normal distribution. Differentially expressed transcripts (DETs) satisfying padjusted<=0.05 were considered significant and selected for downstream analysis.

Functional profiling

We subjected the significant DETs to functional profiling with clusterProfiler (v4.10.0). Gene set enrichment analysis (GSEA) was performed with clusterProfiler’s gseGO function using the Gene Ontology (GO) database. All GO’s categories ("BP”, biological processes; “CC”, cellular components; and “MF”, molecular functions) were considered. The correction of p-values was done with the Benjamini-Hochberg (BH) method and GO terms whose padjusted<=0.05 were considered significant. Over-representation analysis (ORA) was carried out with clusterProfiler’s enrichGO function with the same settings and significance threshold specified for GSEA.

Single-cell RNA sequencing data analysis

Data retrieval

GBM single-cell RNA-seq (scRNA-seq) data were retrieved from three distinct sources. First, we downloaded raw count data and corresponding metadata of 3589 single cells from the work of Darmanis et al. (http://gbmseq.org); second, log-transformed, transcripts per million (TPM)-normalized data and associated information of 7930 single cells were obtained from the study of Neftel et al. via the Broad Institute Single-Cell Data Portal (https://singlecell.broadinstitute.org/single_cell/study/SCP393/single-cell-rna-seq-of-adult-and-pediatric-glioblastoma); and third, raw count data and related metadata of 44712 single cells were downloaded from the work of Richards et al. through Broad Institute’s Portal as well (https://singlecell.broadinstitute.org/single_cell/study/SCP503/gradient-of-developmental-and-injury-reponse-transcriptional-states-define-functional-vulnerabilities-underpinning-glioblastoma-heterogeneity).

Cell filtering

Cells from all three works had already been filtered by the authors of the respective studies regarding scRNA-seq data quality control metrics. Cell type classification of each cell was provided by the authors in the studies’ associated metadata. Here, we filtered cells according to their cell type assignments. As we aimed to investigate PRNP expression implications only in neoplastic cells, we isolated those malignant cells in the three datasets, resulting in 1091 cells for Darmanis et al. dataset; 6896 cells for Neftel et al. dataset; and 14207 cells for Richards et al. dataset. Recurrent tumors were excluded from the datasets, as well as pediatric samples and IDH-mutant tumors.

Dimensionality reduction

We processed the data of malignant cells from each of the three datasets using FUSCA (v1.3.1) [62]. For the Darmanis et al. and Richards et al. datasets, we normalized and scaled count data according to all genes, subsequently performing principal component analysis (PCA). In the case of the Neftel et al. dataset, as count data were already normalized, we skipped the normalization step and only scaled the data according to all genes, performing PCA afterwards. Uniform manifold approximation and projection (UMAP) was then employed for all three datasets, considering the 15 principal components with the highest standard deviations.

Correlation between PRNP and all genes expressed by glioblastoma cells

We calculated the Pearson correlation between PRNP expression and each gene’s expression levels for each dataset individually. Statistical significance was defined as p<=0.05. We compared the correlation results of each dataset and selected the genes that were positively correlated with PRNP in all three datasets (n=2541) for functional profiling. We performed ORA with clusterProfiler’s enrichGO. All GO’s categories were considered, and p-values were corrected with the BH method. Terms whose padjusted<=0.05 were considered significant.

Classification of glioblastoma cells according to PRNP expression

Using normalized count data, we classified the malignant cells of each dataset according to their PRNP levels: those which displayed expression of this gene were assigned as PRNP+, while the ones that did not express PRNP were considered PRNP-.

Identification of marker gene signatures in PRNP + and PRNP - glioblastoma cells

Once we established groups of GBM cells either expressing PRNP or not in each dataset, we proceeded to the identification of marker gene signatures that characterized PRNP+ and PRNP- cells. For each case, we utilized FUSCA’s findSignatures function to find marker genes using Wilcox statistical tests applying a cut-off of p<=0.05. Marker genes of PRNP+ and PRNP- GBM cells of each dataset were subjected to GSEA and ORA with gseGO and enrichGO, respectively, using all categories of the GO database, with p-values being corrected with the BH method. Terms whose padjusted<=0.05 were considered significant.

Comparison of the results from each RNA sequencing dataset (bulk and single-cell)

We compared the upregulated transcripts in the PRNPhigh group from TCGA bulk RNA-seq analysis with marker gene signatures of PRNP+ GBM cells from each scRNA-seq dataset. We subjected the 73 (PRNP + 72) genes that had higher expression in PRNPhigh/PRNP+ in all RNA-seq datasets to ORA with enrichGO with the same settings specified above. Since there were only a few genes in the gene panel, we also used gProfiler2 (v0.2.2) to perform ORA including GO and other databases, such as KEGG, Reactome, and WikiPathways.

Survival analysis

Impact of PRNP expression on overall survival

We stratified the 124 primary IDHwt GBM samples into two groups based on PRNP expression levels using the optimal cut-off determined by maxstat.test (v0.7-25) in R. The optimal cut-off was calculated with log-rank statistics and the HL method for p-value approximation. Samples above the optimal cut-off were considered to have high levels of the gene of interest, whilst samples below that threshold were considered to have low levels. Finally, for Kaplan-Meier survival estimates, we used the survival (v3.5-7) and survminer (v0.4.9) packages.

Impact of gene set signatures on overall survival

We used log10(CPM+1)-normalized counts of all 124 primary IDHwt GBM samples to calculate gene set signature scores with gene set variation analysis (GSVA) [63]. GO gene sets were obtained from the Molecular Signatures Database (MSigDB) (v.7.5.1). We used the maxstat.test (v0.7-25) to calculate the optimal cut-off based on log-rank statistics and the HL method for p-value approximation. Samples above the optimal cut-off were considered to have high levels of the gene set of interest, whilst samples below that threshold were considered to have low levels. Finally, for Kaplan-Meier survival estimates, we used the survival (v3.5-7) and survminer packages (v0.4.9).

Proteome analysis

Data retrieval, processing, and establishment of groups with different PrPC levels

GBM proteomics data were obtained from the Proteomics Identifications Database (PRIDE). The selected dataset consisted of patient-derived GBMs (PRIDE Project PXD015545) [36], containing 12 sets with 6 samples each, tagged with 6 different tandem mass tags (TMT). The GBM Global sets 1 to 6 and 7 to 12 were selected for further analysis. Samples were filtered for contaminants, for reverse database, and the ones only identified by site. Next, samples were normalized by dividing each sample from each set by the Global Internal Standard (GIS) of the respective set. Following this, each ratio obtained by the normalization was log-transformed [log2(ratio)]. Then, the median of the samples was subtracted to centralize the normals in each sample, and with that, the groups of sets (1 to 6 and 7 to 12) were joined. Samples of normal tissue, which were 4, and the 12 GIS were removed, remaining 50 GBM samples (10 proneural, 10 neural, 14 classical, and 16 mesenchymal). Following this, samples were filtered according to PrPC log2 (ratio) and separated into the PrPC-High (n=13) and PrPC-Low (n=13) expression groups, and subjected to differential expression analysis, using PrPC-Low as controls.

Differential expression analysis, overrepresentation analysis, and enrichment map visualization

For the comparison between the PrPC-High and PrPC-Low groups, the Benjamini-Hochberg method was used to determine differentially expressed proteins (DEPs) using PrPC-Low as control. In this analysis, we employed FDR≤0.05 and s0=0.1 as cut-offs. Upregulated DEPs in PrPC-High samples were subjected to ORA in g:Profiler using the g:SCS test to determine statistical significance and GO, Reactome, KEGG, and WikiPathways as gene set databases. A cut-off of p<0.05 was employed. ORA results were visualized through an enrichment map [64], which was built using the EnrichmentMap plugin in Cytoscape. Overrepresented terms were clustered with the AutoAnnotate plugin. Protein classes were analyzed using Panther database [65].

Analysis of other solid tumors

Using TCGAbiolinks (v2.30.0), we retrieved the data of adrenocortical carcinoma (ACC), bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC), cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), esophageal carcinoma (ESCA), head and neck squamous cell carcinoma (HNSC), kidney renal clear cell carcinoma (KIRC), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), ovarian serous cystadenocarcinoma (OV), pancreatic adenocarcinoma (PAAD), pheochromocytoma and paraganglioma (PCPG), prostate adenocarcinoma (PRAD), rectum adenocarcinoma (READ), skin cutaneous melanoma (SKCM), stomach adenocarcinoma (STAD), and uterine corpus endometrial carcinoma (UCEC) samples from TCGA. In total, we retrieved the data of 19 solid tumor types other than GBM, spanning different organs and distinct locations of the human body. For each cancer type, we subset only primary tumors and normalized raw count data according to log10(CPM+1), calculating vesicle dynamics signature scores with GSVA based on selected GO gene sets from MSigDB (v7.5.1) afterwards. Then, for each tumor stage, we computed the Spearman correlation between PRNP expression and those vesicle dynamics signatures, considering p<=0.05 statistically significant.

Experimental procedures

Samples

GBM samples were collected during the surgical procedure and immediately frozen in liquid nitrogen after resection by the group of Brain Tumors and Metastases of the Division of Neurosurgical Clinic of the Central Institute of the Department of Neurology of the Hospital das Clinicas of the School of Medicine of the University of Sao Paulo (HC-FMUSP). The casuistry of this project consisted of 103 GBMs confirmed by histopathological analysis of the HC-FMUSP Pathological Anatomy Division and are part of the biorepository samples collected during the Clinical Genome Project, approved by the National Ethics Commission (CONEP) and the Ethics Committee for Analysis of Research Projects (CAPPesq) of HC-FMUSP, under protocol number 830/01. Post-informed consents were obtained from all patients with tumor and epilepsy or from their legal guardians.

Total RNA extraction and cDNA synthesis

Frozen tumor samples were analyzed in 4μm cryosections and stained with hematoxylin and eosin for quality verification. Total RNA from tumors were also extracted with RNeasy Mini Kit. The quality and concentration of RNA were determined by measuring the absorbance in NanoDrop Spectrophotometers (Thermo Fisher Scientific) at 260 and 280 nm, and A260/A280 ratios greater than 1.8 were considered satisfactory for purity. cDNA synthesis was performed by reverse transcription with SuperScript III reverse transcriptase, RNase inhibitor (RNaseOUT), random oligonucleotides, and oligodT, according to the manufacturer's recommendations (Thermo Fisher Scientific). cDNA was treated with 1U RNase H (Thermo Fisher Scientific) at 37°C for 30 minutes and at 72°C for 10 minutes, diluted in Tris-EDTA buffer, and stored at -20°C for subsequent analysis of gene expression.

Quantitative real-time PCR

Gene expression levels in tissue samples were analyzed by reverse transcription quantitative real-time PCR (RT-qPCR). Reactions were performed by the incorporation method of SYBR Green on QuantStudio 3 (Thermo Fisher Scientific). The reactions were performed in duplicates in a final reaction volume of 10 μl, containing 5 μl of Power SYBR Green PCR Master Mix (Thermo Fisher Scientific), 2.5 μl of cDNA, and 2.5 μl of primers in a pre-standardized concentration. The primers used for the qPCR reaction were the following: ANXA1 Forward: GCGGTGAGCCCCTATCCTA; Reverse: TGATGGTTGCTTCATCCACAC; RAB31 Forward: GGGGTTGGGAAATCAAGCATC; Reverse: GCCAATGAATGAAACCGTTCCT; DSTN Forward: ATTTTGTGGGAATGCTTCCTGA; Reverse: GCATCCTTGGAGCTTGCATAG; SYPL1 Forward: AAGATTACGTCCTCATAGGCGA; Reverse: TCGTGTAGCCAACATAAAGCAG; PRNP Forward: AGTCAGTGGAACAAGCCGAG; Reverse: CTGCCGAAATGTATGATGGGC. Amplification conditions were: initial incubation at 50°C for 5 minutes and at 95°C for 10 minutes, followed by 40 cycles of 95°C for 15 seconds and 60°C for 60 seconds. The expression value of the gene was normalized with reference genes as internal control: HPRT Forward: TGAGGATTT GGAAAGGGTGT; Reverse: GAGCACACAGAGGGCTACAA; GUSB Forward: GAAAATACGTGGTTGGAGAGCTCATT; Reverse: CCGAGTGAAGATCCCCTTTTTA; and TBP Forward: AGGATAAGAGAGCCACGAACCA; Reverse: CTTGCTGCCAGTCTGGACTGT. Single-product amplification was confirmed by analyzing its dissociation curve.

RT-qPCR analysis

One hundred three GBM samples were separated from non-tumoral tissue and classified based on Δ-Cq values for normalized PRNP expression, using the geometric mean of TBP, HPRT, and GUSB as endogenous controls. Specifically, samples were stratified into quartiles according to the normalized values, and samples lying above the upper quartile and below the lower quartile, respectively, high (n=26) and low (n=26) PRNP levels, were used for relative gene expression analysis, with PRNP low as the control group. Relative gene expression was assessed by the 2-ΔΔCt method [66], followed by Student's t-test to assess the significance of the results.

Availability of data and materials

All RNA-seq datasets analyzed in this study are publicly available for download. GBM single-cell RNA-seq (scRNA-seq) data were retrieved from: http://gbmseq.org; https://singlecell.broadinstitute.org/single_cell/study/SCP393/single-cell-rna-seq-of-adult-and-pediatric-glioblastoma; https://singlecell.broadinstitute.org/single_cell/study/SCP503/gradient-of-developmental-and-injury-reponse-transcriptional-states-define-functional-vulnerabilities-underpinning-glioblastoma-heterogeneity.

Supplementary tables will be made available before publication. The scripts that were used for the analyses and the GBMdiscovery app will be made publicly available at https://github.com/marilenehohmuth/GlioblastomaDiscovery.

References

  1. Wesseling P, Capper D. WHO 2016 Classification of gliomas. Neuropathol Appl Neurobiol. 2018;44:139–50.

    Article  PubMed  CAS  Google Scholar 

  2. Louis DN, Perry A, Wesseling P, Brat DJ, Cree IA, Figarella-Branger D, et al. The 2021 WHO classification of tumors of the central nervous system: a summary. Neuro Oncol. 2021;23:1231–51.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  3. Wang Q, Hu B, Hu X, Kim H, Squatrito M, Scarpace L, et al. Tumor evolution of glioma-intrinsic gene expression subtypes associates with immunological changes in the microenvironment. Cancer Cell. 2017;32:42–56.e6.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Neftel C, Laffy J, Filbin MG, Hara T, Shore ME, Rahme GJ, et al. An Integrative Model of Cellular States, Plasticity, and Genetics for Glioblastoma. Cell. 2019;178:835–849.e21.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  5. Sidaway P. CNS cancer: Glioblastoma subtypes revisited. Nat Rev Clin Oncol. 2017;14:587.

    Article  PubMed  Google Scholar 

  6. Stupp R, Mason WP, van den Bent MJ, Weller M, Fisher B, Taphoorn MJ, Belanger K, Brandes AA, Marosi C, Bogdahn U, Curschmann J, Janzer RC, Ludwin SK, Gorlia T, Allgeier A, Lacombe D, Cairncross JG, Eisenhauer E, Mirimanoff RO; European Organisation for Research and Treatment of Cancer Brain Tumor and Radiotherapy Groups; National Cancer Institute of Canada Clinical Trials Group. Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma. N Engl J Med. 2005;352(10):987–96.

  7. Hatoum A, Mohammed R, Zakieh O. The unique invasiveness of glioblastoma and possible drug targets on extracellular matrix. Cancer Manag Res. 2019;11:1843–55.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  8. Lakomy R, Kazda T, Selingerova I, Poprach A, Pospisil P, Belanova R, Fadrus P, Vybihal V, Smrcka M, Jancalek R, Hynkova L, Muckova K, Hendrych M, Sana J, Slaby O, Slampa P. Real-World Evidence in Glioblastoma: Stupp’s Regimen After a Decade. Front Oncol. 2020;10:840.

  9. Ou A, Alfred Yung WK, Majd N. Molecular mechanisms of treatment resistance in glioblastoma. Int J Mol Sci. 2021;22:1–24.

    Google Scholar 

  10. Wang X, Yu X, Xu H, Wei K, Wang S, Wang Y, Han J. Serum-derived extracellular vesicles facilitate temozolomide resistance in glioblastoma through a HOTAIR-dependent mechanism. Cell Death Dis. 2022;13(4):344.

  11. Ma C, Nguyen HPT, Jones JJ, Stylli SS, Whitehead CA, Paradiso L, Luwor RB, Areeb Z, Hanssen E, Cho E, Putz U, Kaye AH, Morokoff AP. Extracellular Vesicles Secreted by Glioma Stem Cells Are Involved in Radiation Resistance and Glioma Progression. Int J Mol Sci. 2022;23(5):2770.

  12. Zanata SM, Lopes MH, Mercadante AF, Hajj GNM, Chiarini LB, Nomizo R, et al. Stress-inducible protein 1 is a cell surface ligand for cellular prion that triggers neuroprotection. EMBO J. 2002;21:3307–16.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  13. Lopes MH, Hajj GNM, Muras AG, Mancini GL, Castro RMPS, Ribeiro KCB, et al. Interaction of cellular prion and stress-inducible protein 1 promotes neuritogenesis and neuroprotection by distinct signaling pathways. J Neurosci. 2005;25:11330–9.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  14. Coitinho AS, Freitas ARO, Lopes MH, Hajj GNM, Roesler R, Walz R, et al. The interaction between prion protein and laminin modulates memory consolidation. Eur J Neurosci. 2006;24:3255–64.

    Article  PubMed  Google Scholar 

  15. Linden R. The biological function of the prion protein: A cell surface scaffold of signaling modules. Front Mol Neurosci. 2017;10:77.

  16. Lebreton S, Zurzolo C, Paladino S. Organization of GPI-anchored proteins at the cell surface and its physiopathological relevance. Crit Rev Biochem Mol Biol. 2018;53:403–19.

    Article  PubMed  CAS  Google Scholar 

  17. Castle AR, Gill AC. Physiological functions of the cellular prion protein. Front Mol Biosci. 2017;4:19.

  18. Mehrpour M, Codogno P. Prion protein: From physiology to cancer biology. Cancer Lett. 2010;290:1–23.

    Article  PubMed  CAS  Google Scholar 

  19. Ding M, Chen Y, Lang Y, Cui L. The Role of Cellular Prion Protein in Cancer Biology: A Potential Therapeutic Target. Front Oncol. 2021;11:742949.

  20. Yang X, Zhang Y, Zhang L, He T, Zhang J, Li C. Prion protein and cancers. Acta Biochimica et Biophysica Sinica. 2014;46:431–40.

    Article  PubMed  CAS  Google Scholar 

  21. Lopes MH, Santos TG, Rodrigues BR, Queiroz-Hazarbassanov N, Cunha IW, Wasilewska-Sampaio AP, et al. Disruption of prion protein-HOP engagement impairs glioblastoma growth and cognitive decline and improves overall survival. Oncogene. 2015;34:3305–14.

    Article  PubMed  CAS  Google Scholar 

  22. Corsaro A, Bajetto A, Thellung S, Begani G, Villa V, Nizzari M, Pattarozzi A, Solari A, Gatti M, Pagano A, Würth R, Daga A, Barbieri F, Florio T. Cellular prion protein controls stem cell-like properties of human glioblastoma tumor-initiating cells. Oncotarget. 2016;7(25):38638–57.

  23. Iglesia RP, Prado MB, Cruz L, Martins VR, Santos TG, Lopes MH. Engagement of cellular prion protein with the co-chaperone Hsp70/90 organizing protein regulates the proliferation of glioblastoma stem-like cells. Stem Cell Res Ther. 2017;8(1):76.

  24. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

  25. Prado MB, Escobar MIM, Alves RN, Coelho BP, Fernandes CF de L, Boccacino JM, et al. Prion protein at the leading edge: Its role in cell motility. Int J Mol Sci. 2020;21:1–18.

  26. Dias MVS, Teixeira BL, Rodrigues BR, Sinigaglia-Coimbra R, Porto-Carreiro I, Roffé M, et al. PRNP/prion protein regulates the secretion of exosomes modulating CAV1/caveolin-1-suppressed autophagy. Autophagy. 2016;12:2113–28.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  27. Kaffes I, Szulzewsky F, Chen Z, Herting CJ, Gabanic B, Velázquez Vega JE, Shelton J, Switchenko JM, Ross JL, McSwain LF, Huse JT, Westermark B, Nelander S, Forsberg-Nilsson K, Uhrbom L, Maturi NP, Cimino PJ, Holland EC, Kettenmann H, Brennan CW, Brat DJ, Hambardzumyan D. Human Mesenchymal glioblastomas are characterized by an increased immune cell presence compared to Proneural and Classical tumors. Oncoimmunology. 2019;8(11):e1655360.

  28. Longo SK, Guo MG, Ji AL, Khavari PA. Integrating single-cell and spatial transcriptomics to elucidate intercellular tissue dynamics. Nat Rev Genet. 2021;22:627–44.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  29. Richards LM, Whitley OKN, MacLeod G, Cavalli FMG, Coutinho FJ, Jaramillo JE, et al. Gradient of Developmental and Injury Response transcriptional states defines functional vulnerabilities underpinning glioblastoma heterogeneity. Nat Cancer. 2021;2:157–73.

    Article  PubMed  CAS  Google Scholar 

  30. 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:1399–410.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  31. Song Y, Sun X, Duan F, He C, Wu J, Huang X, Xing K, Sun S, Wang R, Xie F, Mao Y, Wang J, Li S. SYPL1 Inhibits Apoptosis in Pancreatic Ductal Adenocarcinoma via Suppression of ROS-Induced ERK Activation. Front Oncol. 2020;10:1482.

  32. Verdoni AM, Smith RS, Ikeda A, Ikeda S. Defects in actin dynamics lead to an autoinflammatory condition through the upregulation of CXCL5. PLoS One. 2008;3(7):e2701.

  33. Rogers MA, Buffolo F, Schlotter F, Atkins SK, Lee LH, Halu A, Blaser MC, Tsolaki E, Higashi H, Luther K, Daaboul G, Bouten CVC, Body SC, Singh SA, Bertazzo S, Libby P, Aikawa M, Aikawa E. Annexin A1-dependent tethering promotes extracellular vesicle aggregation revealed with single-extracellular vesicle analysis. Sci Adv. 2020;6(38):eabb1244.

  34. Chen R, Chen C, Han N, Guo W, Deng H, Wang Y, et al. Annexin-1 is an oncogene in glioblastoma and causes tumour immune escape through the indirect upregulation of interleukin-8. J Cell Mol Med. 2022;26:4343–56.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  35. Suo J, Wang Y, Wang L, Qiu B, Wang Z, Yan A, Qiang B, Han W, Peng X. RAB31 in glioma-derived endothelial cells promotes glioma cell invasion via extracellular vesicle-mediated enrichment of MYO1C. FEBS Open Bio. 2024;14(1):138–47. 

  36. Oh S, Yeom J, Cho HJ, Kim JH, Yoon SJ, Kim H, et al. Integrated pharmaco-proteogenomics defines two subgroups in isocitrate dehydrogenase wild-type glioblastoma with prognostic and therapeutic opportunities. Nat Commun. 2020;11(1):3288. https://doi.org/10.1038/s41467-020-17139-y.

  37. York HM, Coyle J, Arumugam S. To be more precise: The role of intracellular trafficking in development and pattern formation. Biochem Soc Trans. 2020;48:2051–66.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  38. Becker A, Thakur BK, Weiss JM, Kim HS, Peinado H, Lyden D. Extracellular vesicles in cancer: cell-to-cell mediators of metastasis. Cancer Cell. 2016;30:836–48.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  39. Sneeggen M, Guadagno NA, Progida C. Intracellular Transport in Cancer Metabolic Reprogramming. Front Cell Dev Biol. 2020;8:597608.

  40. Russo MN, Whaley LA, Norton ES, Zarco N, Guerrero-Cázares H. Extracellular vesicles in the glioblastoma microenvironment: A diagnostic and therapeutic perspective. Mol Aspects Med. 2023;91:101167.

  41. Fernandes CFL, Coelho BP, Souza MCDS, Boccacino JM, Soares SR, Araújo JPA, Melo-Escobar MI, Lopes MH. Extracellular vesicles throughout development: A potential roadmap for emerging glioblastoma therapies. Semin Cell Dev Biol. 2023;133:32–41.

  42. Xu R, Rai A, Chen M, Suwakulsiri W, Greening DW, Simpson RJ. Extracellular vesicles in cancer — implications for future improvements in cancer care. Nat Rev Clin Oncol. 2018;15:617–38.

    Article  PubMed  CAS  Google Scholar 

  43. Vader P, Breakefield XO, Wood MJA. Extracellular vesicles: Emerging targets for cancer therapy. Trends Mol Med. 2014;20:385–93.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  44. Yekula A, Yekula A, Muralidharan K, Kang K, Carter BS, Balaj L. Extracellular Vesicles in Glioblastoma Tumor Microenvironment. Front Immunol. 2020;10:3137.

  45. Skog J, Würdinger T, van Rijn S, Meijer DH, Gainche L, Curry WT, et al. Glioblastoma microvesicles transport RNA and proteins that promote tumour growth and provide diagnostic biomarkers. Nat Cell Biol. 2008;10:1470–6.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  46. Li CCY, Eaton SA, Young PE, Lee M, Shuttleworth R, Humphreys DT, et al. Glioma microvesicles carry selectively packaged coding and noncoding RNAs which alter gene expression in recipient cells. RNA Biol. 2013;10:1333–44.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  47. Liu S, Sun J, Lan Q. Glioblastoma microvesicles promote endothelial cell proliferation through Akt/beta-catenin pathway. Int J Clin Exp Pathol. 2014;7(8):4857–66.

  48. Robinson SW, Nugent ML, Dinsdale D, Steinert JR. Prion protein facilitates synaptic vesicle release by enhancing release probability. Hum Mol Genet. 2014;23:4581–96.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  49. Lima FRS, Arantes CP, Muras AG, Nomizo R, Brentani RR, Martins VR. Cellular prion protein expression in astrocytes modulates neuronal survival and differentiation. J Neurochem. 2007;103:2164–76.

    Article  PubMed  CAS  Google Scholar 

  50. Hartmann A, Muth C, Dabrowski O, Krasemann S, Glatzel M. Exosomes and the prion protein: More than one truth. Front Neurosci. 2017;11:194.

  51. Yun CW, Lee JH, Go G, Jeon J, Yoon S, Lee SH. Prion protein of extracellular vesicle regulates the progression of colorectal cancer. Cancers (Basel). 2021;13(9):2144. Erratum in: Cancers (Basel). 2021;13(14).

  52. Alves RN, Iglesia RP, Prado MB, Escobar MIM, Boccacino JM, Fernandes CF de L, et al. A new take on prion protein dynamics in cellular trafficking. Int J Mol Sci. 2020;21:1–19.

  53. Iglesia RP, Prado MB, Alves RN, Escobar MIM, Fernandes CFL, Fortes ACDS, Souza MCDS, Boccacino JM, Cangiano G, Soares SR, de Araújo JPA, Tiek DM, Goenka A, Song X, Keady JR, Hu B, Cheng SY, Lopes MH. Unconventional Protein Secretion in Brain Tumors Biology: Enlightening the Mechanisms for Tumor Survival and Progression. Front Cell Dev Biol. 2022;10:907423.

  54. Chen DH, Wu QW, Li XD, Wang SJ, Zhang ZM. SYPL1 overexpression predicts poor prognosis of hepatocellular carcinoma and associates with epithelial-mesenchymal transition. Oncol Rep. 2017;38:1533–42.

    Article  PubMed  CAS  Google Scholar 

  55. Liu L, He Q, Li Y, Zhang B, Sun X, Shan J, et al. Serum SYPL1 is a promising diagnostic biomarker for colorectal cancer. Clinica Chimica Acta. 2020;509:36–42.

    Article  CAS  Google Scholar 

  56. Schuh M. An actin-dependent mechanism for long-range vesicle transport. Nat Cell Biol. 2011;13:1431–6.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  57. Tabb JS, Molyneaux BJ, Cohen DL, Kuznetsov SA, Langford GM. Transport of ER vesicles on actin filaments in neurons by myosin V. J Cell Sci. 1998;111(Pt 21):3221–34.

  58. Zhang D, Wang W, Zhou H, Su L, Han X, Zhang X, Han W, Wang Y, Xue X. ANXA1: An Important Independent Prognostic Factor and Molecular Target in Glioma. Front Genet. 2022;13:851505.

  59. Wei D, Zhan W, Gao Y, Huang L, Gong R, Wang W, et al. RAB31 marks and controls an ESCRT-independent exosome pathway. Cell Res. 2021;31:157–77.

    Article  PubMed  CAS  Google Scholar 

  60. Pan Y, Zhang Y, Chen L, Liu Y, Feng Y, Yan J. The critical role of Rab31 in cell proliferation and apoptosis in cancer progression. Mol Neurobiol. 2016;53:4431–7.

    Article  PubMed  CAS  Google Scholar 

  61. Robinson MD, McCarthy DJ, Smyth GK. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2009;26:139–40.

    Article  PubMed  PubMed Central  Google Scholar 

  62. Lummertz da Rocha E, Rowe RG, Lundin V, Malleshaiah M, Jha DK, Rambo CR, Li H, North TE, Collins JJ, Daley GQ. Reconstruction of complex single-cell trajectories using CellRouter. Nat Commun. 2018;9(1):892.

  63. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-Seq data. BMC Bioinformatics. 2013;14:7.

  64. Merico D, Isserlin R, Stueker O, Emili A, Bader GD. Enrichment map: A network-based method for gene-set enrichment visualization and interpretation. PLoS One. 2010; 5(11):e13984.

  65. Thomas PD, Ebert D, Muruganujan A, Mushayahama T, Albou LP, Mi H. PANTHER: Making genome-scale phylogenetics accessible to all. Protein Sci. 2022;31:8–22.

    Article  PubMed  CAS  Google Scholar 

  66. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2-ΔΔCT method. Methods. 2001;25:402–8.

    Article  PubMed  CAS  Google Scholar 

Download references

Acknowledgements

We are grateful to the authors who generated each dataset analyzed in the present work for making the data publicly available to the whole scientific community. Part of the results shown here (bulk analysis) was based upon data generated by the: TCGA Research Network: https://www.cancer.gov/tcga.

Funding

This research was funded by Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP; 2020/07450-5 to J.M.B., 2020/10671-3 to G.C., 2019/14952-0 to B.P.C., 2019/14741-9 to C.F.L.F., 2017/26158-0 to M.B.P., 2019/11097-1 to M.I.M.E., 2023/12717-9 to B.P.S., 2017/20271-0 and 2018/15557-4 to M.H.L., and 2015/03614-5, 2016/0577-1 and 2021/01207-4 to S.M.O.S and S.K.N.M.), Faculdade de Medicina da USP (FMUSP), Conselho Nacional de Desenvolvimento Científico e Tecnológico (101796/2020-0 to J.M.B., 317214/2021-7 to S.M.O.S., 304541/2020-6 to S.K.N.M., and 409941/2021-2, 307677/2022-2 and 406258/2022-8 to M.H.L), Instituto Serrapilheira (R-2111-39726 to E.L.R) and Fundação de Amparo à Pesquisa e Inovação de Santa Catarina (FAPESC, 2022TR000123 to E.L.R).

Author information

Authors and Affiliations

Authors

Contributions

J.M.B. and M.H.L. conceived the study; J.M.B. and R.S.P. retrieved, processed, and analyzed bulk and single-cell RNA sequencing data; R.S.P. developed the web application; G.C., C.F.L.F. retrieved, processed, and analyzed proteome data; J.M.B., E.L.R., C.F.L.F., B.P.C., and M.H.L. discussed approaches to analyze RNA sequencing data; E.L.R. contributed with the analyses; J.M.B., G.C., B.P.C., C.F.L.F., E.L.R., and M.H.L. analyzed and discussed the results; J.M.B., M.I.M.E., and M.H.L. conceptualized the illustrations; M.I.M.E. designed the illustrations; J.M.B. and M.B.P. selected the primers for qRT-PCR; S.M.O.S. and S.K.N.M. provided patient-derived normal, and tumor samples for qRT-PCR; P.R.S. performed qRT-PCR; P.R.S., S.M.O.S., C.F.L.F., M.B.P. and J.M.B., discussed and analyzed qRT-PCR results; B.P.S., assisted with shiny app hosted at Institution server., J.M.B., C.F.L.F., B.P.C., and M.H.L. wrote the manuscript; J.M.B., R.S.P., C.F.L.F., G.C., P.R.S., M.B.P., B.P.C., M.I.M.E., S.A., G.D.B., S.M.O.S., S.K.N.M., E.L.R., and M.H.L. revised the manuscript. All authors have read and agreed to the final version of the manuscript.

Corresponding authors

Correspondence to Edroaldo Lummertz da Rocha or Marilene Hohmuth Lopes.

Ethics declarations

Ethics approval and consent to participate

For RT-qPCR validation analysis, the casuistry of this project consisted of 103 GBMs and are part of the biorepository samples collected during the Clinical Genome Project, approved by the National Ethics Commission (CONEP) and the Ethics Committee for Analysis of Research Projects (CAPPesq) of the Hospital das Clínicas da Faculdade de Medicina da Universidade de São Paulo (HC-FMUSP), under protocol number 830/01.

Consent for publication

All authors approved the manuscript

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.

Supplementary Information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, 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 changes were made. 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/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Boccacino, J.M., dos Santos Peixoto, R., Fernandes, C.F. et al. Integrated transcriptomics uncovers an enhanced association between the prion protein gene expression and vesicle dynamics signatures in glioblastomas. BMC Cancer 24, 199 (2024). https://doi.org/10.1186/s12885-024-11914-6

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12885-024-11914-6

Keywords