Abstract
Malignant tumors exhibit heterogeneous metabolic reprogramming, hindering the identification of translatable vulnerabilities for metabolism-targeted therapy. How molecular alterations in tumors promote metabolic diversity and distinct targetable dependencies remains poorly defined. Here we create a resource consisting of lipidomic, transcriptomic, and genomic data from 156 molecularly diverse glioblastoma (GBM) tumors and derivative models. Through integrated analysis of the GBM lipidome with molecular datasets, we identify CDKN2A deletion remodels the GBM lipidome, notably redistributing oxidizable polyunsaturated fatty acids into distinct lipid compartments. Consequently, CDKN2A deleted GBMs display higher lipid peroxidation, selectively priming tumors for ferroptosis. Together, this study presents a molecular and lipidomic resource of clinical and preclinical GBM specimens, which we leverage to detect a therapeutically exploitable link between a recurring molecular lesion and altered lipid metabolism in GBM.
Keywords: Glioblastoma, CDKN2A, Shotgun Lipidomics, Triacylglyceride, Ferroptosis, Lipid Peroxidation, Lipid Droplet, GPX4, RNA Sequencing
Graphical Abstrct
eTOC Blurb (~50 words)
Minami et al. integrate lipidomic, transcriptomic, and genomic profiling data to identify altered lipid metabolism in CDKN2A deleted glioblastoma (GBM). CDKN2A deletion remodels the distribution of polyunsaturated fatty acids into different lipid compartments, sensitizing GBMs with CDKN2A loss to lipid peroxidation and ferroptosis in vitro and in vivo.
INTRODUCTION
Molecular alterations in cancer rewire metabolism resulting in specific metabolic vulnerabilities that could be exploited for therapy.1–5 Large-scale multi-omic studies in cancer cell lines have further advanced our understanding of the intersection between molecular diversity and metabolic heterogeneity, indicating distinct metabolic vulnerabilities may exist within defined molecular subsets of cancer.6,7 However, the translatability of these findings into effective cancer treatments are potentially limited due to gaps in our understanding of the molecular and metabolic fidelity in conventional cell lines grown in vitro.8–12 Thus, accounting for both the molecular diversity and environment of human tumors is crucial to identify translatable metabolic targets for patients with cancer.
Glioblastomas (GBM) are universally lethal brain tumors that exhibit extensive genetic and transcriptional heterogeneity.13,14 Specific oncogenic alterations (e.g., EGFR mutation, IDH mutation) reprogram aspects of GBM cellular lipid metabolism, uncovering potential opportunities to exploit altered tumor metabolism for GBM therapy.15–19 Whether the diverse molecular landscape found in GBM leads to heterogeneous lipidomic phenotypes and distinct metabolic vulnerabilities is unknown. Here we developed a large library of molecularly annotated GBM tumors and derivative in vitro and in vivo orthotopic models to search for clinically relevant relationships between specific molecular lesions and lipid metabolic phenotypes.
RESULTS
A molecular and lipidomic resource of glioblastoma tumors and derivative pre-clinical models
84 tumors from patients with GBM were obtained from surgical resections and profiled using lipidomic, transcriptomic and whole exome sequencing methods (Figure 1A, S1A). We also analyzed 29 GBM orthotopic xenografts (PDX) and 43 gliomasphere (GS) cell cultures derived from our GBM tumor cohort (Figure 1A). For each sample, shotgun lipidomics quantified up to 1,020 distinct lipid species from 15 lipid subclasses, revealing considerable variation in the lipidome across GBM tumors as well as the derivative models (Figure 1B, S1B). Genomic analysis demonstrated that our tumor sample cohort contained the molecular lesions within the core pathways known to be altered in GBM (e.g., RTK/RAS/PI3K, P53, and RB signaling pathways),13 and encompass diverse patient characteristics that include tumor transcriptional subtypes (Classical, Mesenchymal, Proneural),14 recurrence, age, gender, and ethnicity (Figure 1B, S1A). Likewise, our PDX and GS models recapitulate both the recurring genetic lesions and transcriptional subtypes of GBM tumors (Figure 1B). Together, this large collection of clinical and preclinical GBM specimens capture and highlight the heterogeneity across both the molecular and lipidomic landscape of GBM.
Figure 1: A molecular and lipidomic resource of glioblastoma tumors and derivative pre-clinical models.
(A) Schematic representation of data collection workflow across glioblastoma (GBM) bulk tumor tissue (n=84) and patient derived models (orthotopic xenografts (n=29) and gliomaspheres (n=43)). (B) Plot describing sample type, subtype, lipid composition (class % of total), copy number alterations and mutations in GBM, tumor recurrence, patient age, and gender. Lipid composition is plotted as a z-score across samples within a tumor microenvironment type. (C) Correlation matrix of lipid species (species % of total, n=918 lipid species) and lipid related gene expression (n=999 genes) across GBM tumors. (D) Lipid class and saturation enrichment for each lipid cluster (significance calculated by hypergeometric test against total lipidome with FDR adjustment). Dot plot summarizing gene cluster-lipid cluster associations. Significant associations derived from Kolmogorov-Smirnov testing are outlined in black (see Methods). Gene clusters are labeled with summarized interpretations of enriched Gene Ontology terms. (E) Pie charts representing the log odds ratio scores of over-enriched lipid clusters for differentially abundant lipid species in altered samples relative to WT samples. See also Fig. S1, S2, and Table S1.
Lipids are structurally and functionally diverse molecules synthesized and regulated by multiple metabolic pathways and enzymes.20 To explore whether the heterogeneity observed in the GBM lipidome is linked to variation in lipid metabolic gene expression programs, we correlated lipid species compositions with the gene expression of lipid metabolism-related genes. Using our patient sample cohort, we found that lipid metabolic genes and lipid species clustered primarily into two groups separating triacylglycerides (TAGs) from other lipid species (Figure 1C, S2A). Lipid clusters within the TAG group further delineated monounsaturated (L1) and polyunsaturated (L2) TAGs, while clusters within the second group separated phospholipids (PC, PE), free fatty acids, and diacylglycerides (L3) from ether lipids and ceramides (L4) (Figure 1C, 1D, S2A). Gene clusters showed strong correspondence with associated lipid classes including TAGs (G1, G2), phospholipids (G3, G4), ether lipids (G3, G5), and ceramides (G4, G5) based on Gene Ontology (GO) term enrichment. Moreover, GO terms related to acyl tail desaturation were enriched in gene clusters (G1, G3) that correlated with lipid groups containing polyunsaturated fatty acids (PUFA) (L2 and L3) (Figure 1D, S2B, Table S1). We conclude that specific differences in lipid metabolic gene expression found across the somatic landscape of GBM are associated with distinct lipid clusters containing species representing specific subclasses and fatty acid saturation.
Next, to determine if specific oncogenic events align with distinct lipid metabolic clusters, we annotated the lipid species differentially abundant among previously defined oncogenic regulators of lipid metabolic reprogramming in GBM and/or other cancers (e.g., IDH mutation, TP53 mutation, EGFRvIII mutation, PTEN deletion/mutation, RB mutation, MDM2/4 amplification).18,19,21–24 Among the various oncogenic lesions represented in our GBM tumor cohort, we found significant heterogeneity both in the number of lipids that differed and in their associated enrichment for specific lipid clusters (Figure 1E, S2C). For example, IDH mutant GBMs were enriched in lipid clusters containing phospholipids and MUFAs (L3 and L4) (Figure 1E and S2C). This result confirms recent findings.18,19 Conversely, enriched lipids in EGFRvIII mutant GBMs consisted of TAGs and PUFAs (L2 cluster) in addition to phospholipids (L3 cluster). These observations support previous findings suggesting that mutant EGFR drives the synthesis of complex fatty acids containing PUFAs as well as phospholipids, both of which support oncogenic signaling.16,22 Notably, the distinct lipid composition signatures that distinguish IDH and EGFR mutant tumors were mirrored in their respective enrichment patterns observed among the lipid gene clusters (Figure S2D). Collectively, these data demonstrate that the heterogeneity in the GBM lipidome may arise from the diversity in oncogenic alterations among individual tumors.
CDKN2A deletion alters the GBM lipidome
To identify cell intrinsic regulators of GBM lipid metabolism, we next performed a correlation analysis from our patient-derived gliomaspheres between the lipidome (n=835 lipids) and protein-coding genes (n=10,778 genes). Genes were filtered to remove transcripts with low expression levels and both genes and lipid species were filtered to remove features with low variances (see Methods). Hierarchical clustering of genes with the most variable correlations (n=2,695 genes) identified two distinct gene-lipid groupings (Figure 2A, see Methods). To examine whether specific genomic alterations associated with these clusters, we scored individual samples for their cluster enrichment (Figure S3A) and then tested their association with gene-level somatic point mutations and locus-level copy number alterations (Figure 2B, S3B, see methods). While non-significant hits included mutant RB1 and TP53, both of which drive tumor lipid metabolic reprogramming,21,23 the only significant genetic alteration delineating these groups mapped to the chromosome 9p21 locus, the location of the tumor suppressor Cyclin Dependent Kinase Inhibitor 2A (CDKN2A). RNA expression confirmed CDKN2A as the most differentially expressed transcript at this locus among our gene-lipid groups (Figure 2B). CDKN2A, deleted in nearly 60% of GBM tumors, encodes the p14/16 proteins which canonically regulate cell cycle and cell fate pathways.25,26 A role for CDKN2A in regulating lipid composition in cancer has not been reported.
Figure 2: Comprehensive characterization of the GBM lipidome reveals the impact of CDKN2A deletion.
(A) Correlation matrix of lipid composition (species % of class, n=835) and gene expression of protein coding genes (n=2,695) in GS (n=43). (B) Scatterplot of gene- and locus-level somatic alterations (n=24 somatic events covering 109 genes) ranked by their association with the gene set composite score (see Methods). Scatterplot of the average Group 1 gene expression (log2 TPM) and shrunken log2 fold change (Group 2/Group 1) of genes at the chr9p21 locus identified. Shrunken fold changes and adjusted p-values from DESeq2 (See Methods). (C) Varimax rotated principal component analysis (PC1, PC2) of all lipids (species % of class, n=985 species) between CDKN2A WT (blue, n=12) and null (pink, n=31) GS. 50% data ellipses are shown. Quantification of difference in rotated component 1 scores between CDKN2A WT (n=12) and null (n=31) samples. Boxplots depict the median, interquartile range, and extrema (significance calculated by Student’s two-sample t-test). (D) Volcano plot of lipids (species % of class) altered by CDKN2A deletion in GS. Significantly altered lipids colored by species (Student’s t-test p-value < 0.05). Permutation p-value representing significance of differential abundance testing results displayed (see Methods) (E) Number of lipid species altered with CDKN2A deletion binned by previously defined lipid clusters (see Figure 1). (F) Acyl tail length and double bond characteristics of significantly altered lipid species identified in Figure 2D (significance calculated by Student’s two-sample t-test). (G) Acyl tail length and double bond characteristics of significantly altered lipid species in CDKN2A WT GS (GS104 and GS116) infected with: shRNA scramble control (shC) or two different shRNA targeting CDKN2A (KD-1 and KD-2). Each of these shRNAs target both p14 and p16 from the CDKN2A locus (significance calculated by Student’s two-sample t-test). For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****. See also Fig. S3 and S4.
Principal component analysis (PCA), followed by varimax rotation,27 of all measured lipid species differentiated CDKN2A wild-type (WT) from CDKN2A null gliomaspheres along the first rotated component (Figure 2C). Compared to CDKN2A WT gliomaspheres, CDKN2A deletion was associated with significant changes in the levels of 195 individual lipid species distributed across 15 detectable lipid subclasses (Figure 2D), as well as all four of the lipid clusters we defined above (Figure 2E). Closer inspection of these differentially enriched lipids revealed distinct features of acyl tail utilization. Lipid species more abundant in CDKN2A WT gliomaspheres comprised acyl tails that were highly desaturated (i.e., polyunsaturated fatty acids or PUFAs) and greater than 20-carbons in length (i.e., very long chain fatty acids) (Figure 2F). By contrast, the lipid species enriched in CDKN2A null gliomaspheres were comprised of long chain fatty acids that were generally shorter (i.e., 12–18 carbons) and saturated or monounsaturated. Interestingly, these distinct lipidomic characteristics (i.e., desaturation state or acyl tail length) were not generalizable to the broader lipidomes of CDKN2A WT and null cells (Fig. S3C and S3D), indicating that these characteristics were specific to a subset of lipid species and as such, likely a regulated process downstream of CDKN2A. In agreement with this idea, we found that shRNA-mediated silencing of p14/p16 expression in CDKN2A WT gliomaspheres shifted the composition and acyl tail characteristics of differentially abundant lipid species towards that of patient-derived CDKN2A deleted gliomaspheres (Figure 2G, S4A, S4B). Intriguingly, however, the link between CDKN2A deletion and altered lipid composition found in gliomaspheres could not be detected from a prior lipidomic analyses of established glioma cell lines (Figure S4C),6 emphasizing the value and importance of utilizing primary samples. Together, these data suggest that CDKN2A status has a profound effect on the tail length and saturation state of a subset of lipids in GBM.
Consistent with the diverse impact on lipid composition observed in CDKN2A null GBM, RNA sequencing (RNA-seq) analysis of our gliomaspheres revealed changes in the expression level of several lipid metabolic genes associated with CDKN2A deletion (Figure S4D). Transcripts representing all lipid gene clusters distinguished CDKN2A WT from CDKN2A null GBM, with over 160 genes associated with the synthesis of various lipid classes and modifications to acyl tails differentially expressed between the genotypes (Figure S4E). Moreover, the patterns of gene cluster enrichments observed between our CDKN2A WT and null GBM samples were preserved in RNA-seq data of GBM tumors obtained from The Cancer Genome Atlas (TCGA) dataset (n=139 samples with available RNA-seq and somatic alteration profiling) (Figure S4F). Thus, the considerable shifts in lipid composition with CDKN2A deletion are associated with significant changes in lipid gene expression programs.
CDKN2A deletion sensitizes GBM cells to lipid peroxidation and ferroptosis
CDKN2A deletion affects composition across all lipid subclasses, however triacylglycerides were the majority (58%) of the individual lipid species altered in CDKN2A null gliomaspheres (Figure 3A). While there were no differences in total TAG levels (Figure S5A), further examination revealed that the fraction of longer chain, highly desaturated TAGs were significantly reduced in CDKN2A deleted gliomaspheres. Conversely, shorter chain SFA-containing and MUFA-containing TAGs were enriched in CDKN2A deleted cells relative to CDKN2A WT counterparts (Figure 3A). Membrane phospholipids (e.g., phosphatidylethanolamines (PEs), phosphatidylcholines (PCs)) containing long acyl chain PUFAs are highly sensitive to lipid peroxidation, especially compared to phospholipid species with SFAs or MUFAs.28 By contrast, PUFAs localized in TAGs are sequestered within lipid droplets where they are protected against lipid oxidation.29 We hypothesized that the reduced abundance of long chain PUFA TAGs in CDKN2A null GBM tumors would increase lipid peroxidation. Consistent with this hypothesis, we observed greater basal lipid peroxidation, using the lipid oxidation probe C11 BODIPY 581/591, in CDKN2A null cells, when compared to their CDKN2A WT counterparts (Figure 3B, S5B). Likewise, silencing of p14/p16 expression in CDKN2A WT gliomaspheres decreased PUFA TAGs containing longer acyl chain tails (Figure 3C), which was associated with an increase in lipid peroxidation relative to non-silenced control cells (Figure 3D). Thus, we conclude that CDKN2A deletion shifts the acyl tail characteristics of the TAG pool resulting in high basal levels of lipid peroxidation in GBM.
Figure 3: CDKN2A deletion renders GBM susceptible to ferroptosis.
(A) Log2 fold changes of significantly different (Student’s t-test p-value < 0.05) TAG species in CDKN2A null/WT GS, ordered by directionality of change, total tail length, and total double bond number. Stacked bar chart indicating double bond and total carbon composition (Student’s two-sample t-test). (B) Confocal microscopy images showing basal C11-BODIPY lipid peroxidation in CDKN2A WT and CDKN2A null GS. Scale bar = 10 μm (60x magnification). Quantification of basal C11-BODIPY lipid peroxidation across CDKN2A WT and CDKN2A null GS plotted as a bar chart. Each point represents an individual cell. Mean +/− s.d. (Student’s two-sample t-test). (C) Log2 fold changes of significantly different TAG species with shCDKN2A (Student’s t-test p-value < 0.05) ordered by directionality of change, total tail length, and total double bond number. Stacked bar chart indicating double bond and total carbon composition (Student’s two-sample t-test). (D) Quantification of basal C11-BODIPY lipid peroxidation in GS116 (CDKN2A WT) infected with shC, KD-1, or KD-2. See (B). (E) Representative trace of cell death (%) under 2.5 μM RSL3 treatment over time (h) from the Incucyte live cell imaging system in a CDKN2A WT and null GS (two-way ANOVA followed by Tukey’s post-hoc test). (F) Heatmap of cell death (%) after 72 hours of treatment with RSL3 (2.5 μM and 5 μM) in CDKN2A WT and CDKN2A null GS (significance calculated by Student’s two-sample t-test). (G) Confocal microscopy images showing C11-BODIPY peroxidation induced after 24-hour treatment with 2.5 μM RSL3 in CDKN2A WT and CDKN2A null GS. Quantification of C11-BODIPY peroxidation plotted as fold change of RSL3 treatment/DMSO control. Each point represents an individual cell. Mean +/− s.d. (Student’s two-sample t-test). (H) Heatmap of cell death (%) after 72 hours of treatment with 2.5 μM RSL3 treatment with and without the addition of canonical ferroptosis inhibitors (100 μM DFO, 1 μM Ferrostatin-1, or 1 μM Liporoxstatin-1) in CDKN2A WT and null GS (Student’s two-sample t-test). (I) Confocal microscopy images showing C11-BODIPY peroxidation with DMSO control treatment or 2.5 μM RSL3 treatment after 24 hours in GS116 infected with shC, KD-1, or KD-2. See (G). (J) Heatmap of cell death (%) after 72 hours of treatment with 2.5 μM RSL3 treatment with and without the addition 1 μM Ferrostatin-1 in GS104 or GS116 infected with shC, KD-1, or KD-2 (Student’s two-sample t-test). For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****. See also Fig. S5.
Ferroptosis is a non-apoptotic cell death mechanism defined both by a dependence on iron and by the lethal accumulation of lipid peroxides.30 Having found that CDKN2A deletion is associated with increased basal lipid peroxidation, we posited that CDKN2A null gliomaspheres would be sensitized to induction of ferroptosis. Across a panel of CDKN2A null gliomaspheres (n=6), inhibition of glutathione peroxidase 4 (GPX4) with (1S,3R)-RSL3 (RSL3) or ML21031 resulted in pronounced cell death. Conversely, CDKN2A WT gliomaspheres (n=6) were largely insensitive to both compounds (Figure 3E, 3F, S5C). Consistent with the greater susceptibility of CDKN2A null GBMs to lethal lipid peroxidation, GPX4 inhibition increased lipid peroxidation only in CDKN2A null gliomaspheres and not in CDKN2A WT cells (Figure 3G). RSL3 or ML210-induced cell death was suppressed by the ferroptosis inhibitors ferrostatin-1 (Fer-1), liproxstatin-1 (Lip-1), the iron chelator deferoxamine (DFO), or by treatment with deuterated linolenic acid, which is far less sensitive to peroxidation than normal PUFAs (Figure 3H, S5D–F).32–34 Moreover, CDKN2A null gliomaspheres were more sensitive than their WT counterparts to Erastin-2, an inducer of ferroptosis via its ability to inhibit the cystine-glutamate antiporter xCT (Figure S5G). Finally, silencing of p14/p16 in CDKN2A WT gliomaspheres increased both lipid peroxidation and ferroptosis susceptibility (Figure 3I, 3J). Thus, disruption of CDKN2A, and specifically loss of p14/16 expression, renders GBM more susceptible to lipid peroxidation and ferroptosis.
CDKN2A deletion reduces oxidizable PUFA sequestration into lipid droplets
We sought to better understand how CDKN2A deletion selectively sensitizes GBM cells to ferroptosis. As noted, peroxidation of PUFA-containing phospholipids contributes to the execution of ferroptosis35 and sequestration of PUFAs into the TAG pool has been suggested to limit ferroptosis.36 Therefore, we posited that CDKN2A governs the partitioning of oxidizable PUFAs into the TAG pool, consequently modulating sensitivity to lipid peroxidation and ferroptosis. To begin testing this hypothesis, we provided the PUFA arachidonic acid (AA) to CDKN2A WT (n=3) and null (n=3) gliomaspheres and performed lipidomics; AA treatment significantly increased the sequestration of PUFAs into TAG species in CDKN2A WT gliomaspheres. By contrast, AA addition to CDKN2A null cells specifically expanded the AA-containing PC phospholipid pools (Figure 4A). Moreover, AA treatment markedly increased lipid droplet formation in CDKN2A WT gliomaspheres, but not in CDKN2A null counterparts (Figure 4B). Both the expanded PUFA TAG pool and the increase in lipid droplets observed with AA addition to CDKN2A WT GBM cells were blocked by pharmacologic inhibition of DGAT1 and DGAT2 – key enzymes in TAG biosynthesis (Figure 4A, 4B).37 Finally, we saw no differences in the relative abundance of cellular SFAs and MUFAs (e.g., 16:0, 16:1, 18:0 and 18:1) between CDKN2A WT and CDKN2A null gliomaspheres, consistent with the genotypes having similar de novo fatty acid synthetic capability (Figure S6A).
Figure 4: PUFA TAGs are protective against lipid peroxidation and ferroptosis in CDKN2A WT GBM.
(A) Stacked bar chart showing PUFA distribution (% of total) within TAGs, PCs, and PEs in CDKN2A WT (n=3) and null (n=3) GS under basal (AA −, DGATi −) conditions or with 24 hours 75 μM AA treatment without (AA+, DGATi−) or with co-treatment with 20 μM DGAT1 and 10 μM DGAT2 inhibitors (AA+, DGATi+). Mean +/− s.e.m. (Student’s two-sample t-test). (B) Confocal microscopy images showing accumulation of the LipidTOX neutral lipid stain in CDKN2A null and WT GS treated with DMSO control, 75 μM AA (AA), or 75 μM AA + 20 μM DGAT1 inhibitor + 10 μM DGAT2 inhibitor (AA+, DGATi+) for 24 hours. Scale bar = 10 μm (60x magnification). Quantification of lipid droplet area per cell across CDKN2A WT (n=4) and null (n=4) GS across these same conditions relative to its respective DMSO control. Each dot represents the average quantification per cell line for the specified condition. Mean +/− s.d. (Student’s two-sample t-test). (C) Confocal images showing C11-BODIPY peroxidation induced after 24-hour treatment with 75 μM AA (AA), 75 μM AA + 2.5 μM RSL3 (AA+, RSL3+) or 75 μM AA + 2.5 μM RSL3 + 20 μM DGAT1 inhibitor + 10 μM DGAT2 inhibitor (AA+, RSL3+, DGATi+) in CDKN2A null and WT GS. Quantification of C11-BODIPY peroxidation plotted as fold-change relative to 24 hours of 75 μM AA treatment. Each dot represents the average quantification per cell line. Mean +/− s.d. (Student’s two-sample t-test). (D) Cell death (%) after 72 hours of treatment with the same conditions shown in Fig. 4C in CDKN2A WT (n=3) and null (n=3) GS. Each dot represents the average cell death per cell line for the specified condition. Mean +/− s.d. (Student’s two-sample t-test). (E) Confocal images showing accumulation of the LipidTOX neutral lipid stain (L image set) with 24 hours of the same treatments indicated in (B) and C11-BODIPY peroxidation (R image set) with 24 hours treatment of the same conditions indicated in (C) in GS116 (CDKN2A WT) infected with shC, KD-1, or KD-2. Quantification of C11-BODIPY peroxidation plotted as fold-change relative to 24 hours of 75 μM AA treatment. Each point represents the average quantification per image for the specified condition and cell line. Mean +/− s.d. (Student’s two-sample t-test). (F) Cell death (%) after 72 hours of the same treatment indicated in (D) in GS116 (CDKN2A WT) infected with shC, KD-1, or KD-2. Each point represents a biological replicate. Mean +/− s.d. (significance calculated by Student’s two-sample t-test). For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****. See also Fig. S6.
Our results demonstrate that CDKN2A genetic status impacts intracellular PUFA distribution into distinct lipids classes (e.g., TAGs and phospholipids). Moreover, these CDKN2A mediated changes in lipid composition underlie susceptibility of GBM to lipid peroxidation and ferroptosis. In agreement with this idea, providing exogenous AA markedly increased RSL3-mediated lipid peroxidation and ferroptosis in CDKN2A null GBM cells. In distinction, DGAT1/2 inhibitors were required for AA-treated CDKN2A WT cells to increase RSL3-induced lipid peroxidation and ferroptosis to the levels observed in CDKN2A null gliomaspheres (Figure 4C, 4D, S6B). Taken together, these data support a model in which p14/16 signaling regulates the shunting of PUFAs into the TAG compartment, resulting in protection from lipid peroxidation and ferroptosis in GBM.
Loss of p16 is sufficient to render GBM sensitive to ferroptosis
CDKN2A encodes two proteins, p14 and p16. Canonically, p14 dictates cell fate by indirectly stabilizing p53, while p16 suppresses tumor formation by inhibiting CDK4/6.25,26 While a direct role for p14 or p16 in regulating cancer lipid metabolism has not been described, recent findings show that p16 contributes to hepatic lipid homeostasis.38 These observations led us to posit that p16 may regulate GBM lipid metabolism and consequent sensitivity to ferroptosis. To test this hypothesis, we first ablated CDKN2A using CRISPR-Cas9 in CDKN2A WT gliomaspheres using sgRNAs targeting exon 2 of CDKN2A (shared between p14 and p16) (Figure S7A, S7B).39 Knockout of exon 2 of CDKN2A (CDKN2A ex2 KO) abrogated p14 and p16 protein levels and, similar to our findings in CDKN2A null gliomaspheres, reduced the enrichment for desaturated long chain acyl tail lipid species (Figure S7C), with a notable decrease in the fraction of PUFA-containing TAGs (Figure 5A).
Figure 5: Loss of p16 is sufficient to render GBM sensitive to ferroptosis.
(A) Abundance of PUFA TAGs (n=260 lipid species, species % of class) between GS116 (CDKN2A WT) transfected with Control, CDKN2A ex2 KO, or CDKN2A ex2 KO + p16 addback. Z-scored across GS116 control, CDKN2A ex2 KO, and CDKN2A ex2 KO+p16. (B) Confocal microscopy images showing accumulation of the LipidTOX neutral lipid stain in GS116 Control, CDKN2A ex2 KO, or CDKN2A ex2 KO + p16 addback. Cells treated with DMSO control, 20 μM DGAT1 inhibitor + 10 μM DGAT2 inhibitor (DGATi), 75 μM AA (AA), or 75 μM AA + 20 μM DGAT1 inhibitor + 10 μM DGAT2 inhibitor (AA, DGATi) for 24 hours. Scale bar = 10 μm (60x magnification). Quantification of lipid droplet area/cell (FC relative to control). Each dot represents the average quantification per cell line for the specified condition. Mean +/− s.d. (Student’s two-sample t-test). (C) Confocal images showing basal MDA immunofluorescence staining. Quantification of MDA average fluorescent intensity/cell plotted. Each dot represents the average intensity per cell line. Mean +/− s.d. (significance calculated with parametric unpaired t-test). (D) Heatmap of cell death (%) induced by DMSO vehicle or 2.5μM RSL3 after 72 hours of treatment (Student’s two-sample t-test). For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****. See also Fig. S7.
Next, we reintroduced p16 into our CDKN2A ex2 KO cells to specifically restore p16 protein levels while maintaining reduced p14 expression (Figure S7A, S7B). p16 addback rescued many of the lipidomic features to that of CDKN2A WT control gliomaspheres, including enriched lipid species with very long chain and desaturated acyl tails (Figure S7C) and an increase in PUFA TAGs (Figure 5A). The rescue of PUFA TAG levels with p16 addback was coupled with a restored ability to shunt AA into lipid droplets (Figure 5B). These data further support a model where CDKN2A partitions oxidizable PUFA TAGs in lipid droplets. Correspondingly, we observed that addback of p16 expression could completely reverse basal lipid peroxidation levels, as determined by measuring both the lipid peroxidation byproduct malondialdehyde (MDA) (Figure 5C), and ferroptosis sensitivity relative to parental CDKN2A ex2 KO cells (Figure 5D). These results led us to conclude that p16 loss is sufficient to remodel components of the GBM lipidome, rendering GBM gliomaspheres susceptible to ferroptosis.
GPX4 inhibition decreases CDKN2A-deleted GBM tumor burden in vivo
We next explored whether CDKN2A genetic status influences the lipidome and ferroptotic potential of GBM tumors in vivo. In agreement with our cell culture findings, CDKN2A WT orthotopic PDX tumors (n=11) were enriched with TAGs containing longer acyl tails and more double bonds relative to CDKN2A null tumors (n=18) (Figure 6A). CDKN2A null PDX tumors also displayed higher levels of lipid peroxidation when compared to CDKN2A WT PDX tumors (Figure 6B, S8A). This result suggests that the CDKN2A null tumors would be more sensitive to the induction of ferroptosis. To test this hypothesis, we disrupted GPX4 expression using CRISPR-Cas9 in CDKN2A WT (n=2) and CDKN2A null (n=2) gliomaspheres (Figure 6C). As anticipated, genetic GPX4 ablation resulted in significantly higher ferroptosis in CDKN2A null gliomaspheres compared to WT gliomaspheres in vitro (Figure S8B–D). Next, we engineered a separate cohort of CDKN2A WT and CDKN2A null gliomaspheres with either non-targeting control (NT) sgRNAs or three different GPX4 sgRNAs. Transduced gliomaspheres were expanded in media containing Fer-1 (to prevent changes in viability with GPX4 ablation) for 72 hr and then implanted into the brains of NOD-scid gamma (NSG) mice. All three GPX4 sgRNAs significantly extended the survival of two different CDKN2A null PDX models (PDX025 and PDX187) relative to animals implanted with cells transduced with NT sgRNAs (Figure 6D). By contrast, disruption of GPX4 had no survival benefit for mice bearing CDKN2A WT PDX tumors (PDX005 and PDX208) (Figure 6D). Together, these results demonstrate that GPX4 inhibition significantly and selectively improves the outcomes of mice bearing CDKN2A null GBM tumors.
Figure 6. CDKN2A deletion induces lipid peroxidation and susceptibility to ferroptosis in GBM orthotopic xenografts.
(A) Heatmap showing enrichment of conserved TAGs between GS and PDX. TAGs enriched in CDKN2A WT tumors have a log2FC(null/WT)<0 (blue), while TAGs enriched in CDKN2A null tumors have a log2FC(null/WT)>0 (pink). Total double bonds and total carbon length annotated within heatmap. (B) Immunofluorescence microscopy of CDKN2A WT and CDKN2A null tumors stained for malondialdehyde (MDA) (orange), GFP (GFP tumor) (green) and DAPI (blue). Scale bar represents 50 μm (20x magnification). Quantification of MDA (mean fluorescent intensity) across CDKN2A WT (n=8) and null (n=6) PDX tissue. Box-plots show mean (middle-line) +/− s.d. Whiskers represent minimum and maximum data points (significance calculated by Mann-Whitney t-test). (C) Western blot of GPX4 and Actin in CDKN2A null gliomaspheres (GS025, GS187) and WT (GS005, GS208) infected with: nothing (parental), CRISPR non-targeting control (NT), or sgRNAs targeting GPX4 (sgGPX4-1, sgGPX4-2, sgGPX4-3). Molecular weight of closest ladder marker indicated to the right. (D) Kaplan-Meier plot showing probability of survival (%) in mice harboring CDKN2A null and WT orthotopic tumors infected with NT (n=6), sgGPX4-1 (n=6), sgGPX4-2 (n=6 for PDX025, PDX005, PDX208, n=10 for PDX187), or sgGPX4-3 (n=6). X-axis indicates days post tumor implantation. Significance calculated with Log-rank (Mantel-Cox) test between the survival curves for NT vs individual guides. (E) Heatmap showing enrichment of conserved TAGs between GS and tumors from patients with GBM. TAGs enriched in CDKN2A WT tumors have a log2FC(null/WT)<0 (blue), while TAGs enriched in CDKN2A null tumors have a log2FC(null/WT)>0 (pink). (F) Proposed mechanism of ferroptotic sensitivity in CDKN2A WT and CDKN2A null GBM. For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****. See also Fig. S8.
Finally, we asked whether the specific PUFA TAG lipid phenotype linked with CDKN2A genetic status was conserved in tumors from patients with GBM. Like our preclinical GBM models, CDKN2A deletion in GBM tumors (n=50) was associated with reduced fraction of TAGs containing highly desaturated long acyl tails relative to CDKN2A WT tumors (n=34) (Figure 6E). These data emphasize the clinical relevance of our cell culture and xenograft findings, and further supports the conclusion that lipid metabolic remodeling due to loss of CDKN2A increases membrane lipid peroxidation and consequently primes GBM for ferroptosis (Figure 6F).
DISCUSSION
The influences of molecular heterogeneity and environment on the metabolic features of cancer cells can confound our ability to identify translatable metabolic vulnerabilities that are exploitable in therapy.7,10,11 Here, we created a resource containing both molecular and lipidomic data from GBM tumor samples, together with derivative in vitro and in vivo models. This dataset will enable investigations into the genetic influences on GBM lipid metabolism and potential impacts of the tumor environment. We demonstrate the utility of this resource by integrating molecular and lipidomic features of our GBM sample cohort to discover that deletion of a recurring tumor suppressor, CDKN2A, remodels the GBM lipidome. Crucially, GBM lipid composition is similarly remodeled with CDKN2A deletion in the human brain, as a xenograft in the mouse brain, and when grown as gliomaspheres in cell culture. This concordance across environments (in vitro and in vivo) increases the translational potential of this metabolic vulnerability in patients with GBM having CDKN2A deletion. We anticipate and encourage others to use this dataset to uncover other potential therapeutically actionable approaches for this deadly disease.
How CDKN2A regulates the GBM lipidome remains less clear. Our gene expression data indicate large changes in lipid metabolic programing with deletion of CDKN2A in GBM, and we posit it unlikely that the molecular underpinning of how CDKN2A influences the lipidome can be attributed to a single lipid metabolic gene or pathway. Moreover, there is tight substrate and feedback regulation of lipid enzymes, so the lipidome of CDKN2A null GBM likely represents a composite of transcriptional and post-translation regulation. Our data show that p16 encoded by CDKN2A acts as a fundamental regulator of lipid composition and protection from ferroptotic priming in GBM through its ability to alter acyl tail composition and shunting of fatty acids into TAGs within lipid droplets. Interestingly, CDKN2A has been previously linked to metabolic syndrome, to hepatic fatty acid oxidation and to ketogenesis, as well as to adipose tissue differentiation.38,40,41 Although the molecular mechanisms linking CDKN2A with lipid metabolism in normal tissues are also not well understood, these observations, combined with our studies, suggest a fundamental role for CDKN2A in regulating cellular lipid composition in normal and neoplastic tissues. Undoubtedly, deep mechanistic work deconvoluting how a cell cycle regulatory protein influences lipid composition of normal and malignant cells will be an important goal of future studies. Finally, as the expression of several other genes adjacent to CDKN2A are similarly downregulated with deletion of this locus, including those encoding metabolic enzymes (e.g., MTAP, HACD4), we cannot exclude the possibility that their loss also contributes to altered GBM lipid metabolism in CDKN2A deleted GBMs.
The ferroptosis pathway is emerging as a candidate therapeutic target in several cancers.42,43 We find that loss of CDKN2A primes a large subset of patient-derived GBM models to ferroptosis both in vitro and in orthotopic xenografts. One important concept in ferroptosis is that sequestration of oxidizable PUFAs into specific lipid species or subcellular compartments can buffer susceptibility to oxidative stress and induction of cell death.29,36 Supporting this mechanism, our data demonstrate that the sensitivity of CDNK2A null GBM to ferroptosis is associated with differences in the partitioning of oxidizable PUFAs into TAGs and lipid droplets. Our results provide proof-of-concept evidence that GPX4 may be a relevant therapeutic target for a large, genetically stratified subset of patients with GBM. There are currently no pharmacological GPX4 inhibitors suitable for use in the brain; our findings provide a strong rationale for developing brain penetrant small molecules that induce lipid peroxidation-mediated GBM cell death via inhibition of GPX4 or other ferroptotic targets. In summary, our study provides a framework for the combined use of preclinical and clinical tumor samples to investigate the impact of molecular diversity on metabolic heterogeneity, with the goal of identifying clinically relevant metabolic vulnerabilities.
STAR Methods
RESOURCE AVAILABILITY
Lead Contact
Further information and requests for resources and reagents should be directed to and will be fulfilled by the lead contact, David Nathanson ([email protected]).
Materials Availability
This study did not generate new unique reagents.
Data and Code: Whole exome sequencing and mutation data has been deposited at the database of Genotypes and Phenotypes (dbGaP) and is publicly available as of the date of publication. Accession numbers are listed in the key resources table. RNA-seq data has been deposited at dbGaP and is publicly available as of the date of publication. Accession numbers are listed in the key resources table. Shotgun lipidomics data has been deposited at Mendeley and is publicly available as of the date of publication. Accession numbers are listed in the key resources table. Lastly, the data generated in Fig. 1C (Lipid-Gene Correlation Matrix) has been deposited to Mendeley and is publicly available as of the date of publication. The DOI is listed in the key resources table.
This paper does not report original code.
Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.
Key resources table
| REAGENT or RESOURCE | SOURCE | IDENTIFIER |
|---|---|---|
| Antibodies | ||
| P16 INK4A (D3W8G) Rabbit mAb | Cell Signaling | 92803 |
| P14 ARF (4C6/4) Mouse mAb | Cell Signaling | 2407 |
| GPX4 Mouse mAb | R&D Systems | MAB5457-SP |
| β-Actin mAb | Cell Signaling | 3700S |
| Anti-Mouse IgG HRP-linked Antibody | Cell Signaling | 7076 |
| Anti-Rabbit IgG HRP-linked Antibody | Cell Signaling | 7074 |
| MDA Goat pAb | Abcam | ab27644 |
| Anti-Goat Donkey pAb Alexa Fluor 647 | Abcam | ab150135 |
| Bacterial and virus strains | ||
| One Shot Stbl3 Chemically Competent E. coli | Invitrogen | C737303 |
| Biological samples | ||
| Patient Derived Brain Tumor Tissue | UCLA Health | N/A |
| Patient Derived GBM Orthotopic Xenografts | This paper | N/A |
| Patient Derived GBM Cell Lines (Gliomasheres) | This paper | N/A |
| Chemicals, peptides, and recombinant proteins | ||
| RSL3 | Selleckchem | S155 |
| Ferrostatin-1 | Cayman Chemicals | 17729 |
| Liproxstatin-1 | Cayman Chemicals | 17730 |
| Desferrioxamine (DFO) | Cayman Chemicals | 14595 |
| DGAT1 and 2 inhibitors | DGAT1: Santa Cruz Biotechnology DGAT2: Millipore Sigma |
DGAT1: 959122-11-3 DGAT2: PZ0233 |
| ML210 | Selleckchem | S0788 |
| Arachidonic Acid (AA) | Sigma Aldrich | 10931–250mg |
| Buthionine Sulfoximine (BSO) | Selleckchem | S9728 |
| Erastin2 | Cayman Chemicals | HY-139087 |
| Deposited data | ||
| Lipidomics Data (Bulk patient, orthotopic xenograft, gliomasphere, isogenics) | This paper | Mendeley doi:10.17632/kjtdgk3f25.1 |
| Standardized WES Copy Number Data (Patient, orthotopic xenograft, gliomasphere) | This paper | dbGaP (NCBI) Accession: phs003286 |
| Standardized WES Mutation Data (Patient, orthotopic xenograft, gliomasphere) | This paper | dbGaP (NCBI) Accession: phs003286 |
| Standardized RNA-Sequencing Data (Patient, orthotopic xenograft, gliomasphere) | This paper | dbGaP (NCBI) Accession: phs003286 |
| Patient Lipid Gene Expression-Lipid Species Abundance (Species % Total) Correlation | This paper | Mendeley doi:10.17632/kjtdgk3f25.1 |
| Experimental models: Cell lines | ||
| Human: 293FT cells | Invitrogen | R70007 |
| Experimental models: Organisms/strains | ||
| Mouse: NOD.Cg-Prkdcscid Il2rgtm1Wjl/SzJ | The Jackson Laboratory | RRID:IMSR_JAX:005557 |
| Oligonucleotides | ||
| sgRNA targeting sequences for GPX4: See methods | N/A | |
| shRNA targeting p14/p16 | TRCN0000255853 | |
| shRNA targeting p14/p16 | TRCN0000255849 | |
| Recombinant DNA | ||
| LentiCRISPR v2 plasmid | Addgene | Plasmid 52961 (depositing lab: Feng Zhang) |
| pLKO plasmid | Addgene | Plasmid 10878 (depositing lab: David Root) |
| MISSION® pLKO.1-puro Non-Mammalian shRNA Control Plasmid DNA | Sigma-Aldrich | SHC002 |
| pHIV-ARF-mOrange2 | Addgene | Plasmid 110731 (depositing lab: Robert Judson-Torres |
| pEN35-CDKN2A-Ex2-R | Addgene | Plasmid 110737 (depositing lab: Robert Judson-Torres |
| pEN35-CDKN2A-Ex2-L | Addgene | Plasmid 110736 (depositing lab: Robert Judson-Torres |
| pEN35-CKDN2A-Ex2-CMV-EGFP | Addgene | Plasmid 110734 (depositing lab: Robert Judson-Torres |
| Software and algorithms | ||
| ImageJ Fiji | Schindelin et al. (2012) | https://imagej.net/software/fiji/ |
| Zen Blue 3.4 | Zeiss Group | https://www.micro-shop.zeiss.com/en/us/softwarefinder/software-categories/zen-blue/ |
| Prism | GraphPad | https://www.graphpad.com/ |
| BioRender | BioRender | https://biorender.com/ |
| Incucyte | Sartorius | https://www.sartorius.com/en/products/live-cell-imaging-analysis |
| Samtools (v1.10) | Li et al (2011) | http://www.htslib.org/ |
| CutAdapt (v2.8) | Martin (2011) | https://cutadapt.readthedocs.io/en/stable/installation.html |
| BBsplit (June 11, 2018) | BBTools (v38.58) | https://www.sourceforge.net/projects/bbmap/ |
| BWA-MEM (0.7.17-r1188) | Li and Durbin (2012) | https://github.com/lh3/bwa |
| GATK (v4.2.0.0) | Broad Institute | https://gatk.broadinstitute.org/hc/en-us |
| Seal | BBTools (v38.58) | https://www.sourceforge.net/projects/bbmap/ |
| Toil RNA-Seq (v3.3.5) | Vivian et al. 2017 | https://github.com/BD2KGenomics/toil-rnaseq |
| Mutect2 v4.2.0.0 | Broad Institute | https://gatk.broadinstitute.org/hc/en-us/articles/360056969692-Mutect2 |
| VarScan2 (v2.4.3) | Koboldt et al. (2012) | https://varscan.sourceforge.net/ |
| MuSE v1.0rc | MD Anderson Cancer Center | https://bioinformatics.mdanderson.org/public-software/muse/ |
| CNVkit (v0.99) | Talevich et al. (2016) | https://cnvkit.readthedocs.io/en/stable/ |
| R (v4.1.0) | CRAN | https://www.r-project.org |
| GSVA (v1.40.1) | Hänzelmann et al. (2013) | https://bioconductor.org/packages/release/bioc/html/GSVA.html |
| EnrichR | Ma’ayan Lab | https://maayanlab.cloud/Enrichr/ |
| Gene Ontology (GO) Resource | Ashburner et al. (2000) | http://geneontology.org/ |
| Kyoto Encylclopedia of Genes and Genomes (KEGG) | Kanehisa Laboratories | https://www.genome.jp/kegg/ |
| Reactome | Jassal et al, 2020 | https://reactome.org/ |
| Other | ||
| ProLong Gold antifade reagent with DAPI | Invitrogen | P36935 |
| HCS LipidTOX™ Deep Red neutral lipid dye | Invitrogen | H34477 |
| CellTOX cell death green fluorescent dye | Promega | G8741 |
| Incucyte nuclear red fluorescent dye | Sartorius | 4717 |
| C11 BODIPY 581/591 | Thermo Scientific | D3861 |
| Hoechst nuclear stain | Thermo Scientific | 62249 |
| Tumor dissociation kit, human | Miltenyi Biotec | 130-094-929 |
| Brain Tumor Dissociation Kit (P) | Miltenyi Biotec | 130-095-942 |
| Mouse Cell Depletion Kit | Miltenyi Biotec | 130-104-694 |
| Myelin Removal Beads II, h&m&r, 2×4ml | Miltenyi Biotec | 130-096-433 |
| CD45 MicroBeads, human | Miltenyi Biotec | 130-045-801 |
| Debris Removal Solution | Miltenyi Biotec | 130-109-398 |
| Laminin, Mouse | Corning | 354232 |
| Heparin sodium salt from porcine intestinal mucosa | Sigma | H3149 |
| Gibco FGF-Basic (AA 1–155) Recombinant Human Protein | Thermo Scientific | PHG0263 |
| Gibco EGF Recombinant Human Protein | Thermo Scientific | PHG0313 |
| Gibco TrypLE Express Enzyme (1X), phenol red | Thermo Scientific | 12605028 |
| Gibco DMEM/F-12, HEPES | Thermo Scientific | 11330032 |
| Gibco DMEM, high glucose, pyruvate, no glutamine | Thermo Scientific | 10313039 |
| Mycoalert Detection Kit | Lonza | LT07-318 |
| BenchMark™ Fetal Bovine Serum | Gemini Bio | 100-106-500 |
EXPERIMENTAL MODEL AND SUBJECT DETAILS
Patient-derived GBM tumors and gliomaspheres
After explicit informed consent was obtained from patients, all patient-derived tumor tissue was obtained through the UCLA Institutional Review Board (IRB) protocol 10-000655. Sample size for shotgun lipidomics, RNA-seq, and WES were determined based on tissue availability. Patient sex for all derivative samples and models are annotated in Fig. 1B. Tumors were mechanically and enzymatically dissociated using the Miltenyi Biotec Human tumor dissociation kit. Red blood cell lysis buffer removed red blood cells. Antibody-conjugated magnetic beads removed CD45+ cells and myelinated cells in two column-based filtration steps. Primary GBM cells were established and maintained in gliomasphere conditions consisting of DMEM/F12 (Gibco), B27 (Invitrogen), penicillin–streptomycin (Invitrogen), and GlutaMAX (Invitrogen) supplemented with heparin (5 μg/mL, Sigma), EGF (20 ng/mL, Sigma), and FGF (20 ng/mL, Sigma). All cells were grown under 37 °C, 20% O2, and 5% CO2 and were routinely monitored and tested negative for the presence of mycoplasma with a commercially available kit (MycoAlert, Lonza). Gliomasphere cell lines were used at fewer than 15 passages. All cells were authenticated by short tandem repeat (STR) analysis.
Mice
Female, immunocompromised NOD scid gamma (NSG) mice, 6–8 weeks of age, were purchased from the University of California Los Angeles (UCLA) Medical Center animal-breeding facility and Jackson Laboratories. All mice were kept under defined pathogen-free conditions at the AAALAC-approved animal facility of the Division of Laboratory Animals (DLAM) at UCLA. Mice did not receive any prior treatments and were drug and test naïve prior to injection of tumor cells. All animal experiments were performed with the approval of the UCLA Office of Animal Resource Oversight (OARO).
Patient-derived orthotopic xenografts
The majority of the patient-derived orthotopic xenografts (n=27) were developed as previously described.44 Briefly, gliomaspheres were transduced with secreted Gaussia luciferase (sGluc)-GFP reporter to enable non-invasive quantification of tumor burden45 as well as endpoint GFP-guided microdissection of the tumor tissue from the surrounding normal brain. Gliomaspheres were dissociated and injected (5 × 105 cells per injection) into the right striatum of the brain in female NSG mice (8–9 weeks old). Injection coordinates were 2 mm lateral and 1mm posterior to bregma, at a depth of 2mm. For two models, PDX074 and PDX152, freshly isolated and purified patient tumor cells were transduced with sGluc and 5 × 105 cells and inoculated into the brains of mice per the same protocol as gliomaspheres implantation.
Tumor burden was monitored based on 1–2x weekly measurements by sGluc. Upon endpoint, mice were euthanized and tumors were dissected from the brains. The isolated tumor tissue was then submitted for shotgun lipidomics analysis. For WES and RNA-seq analysis, tumor cells were purified via mechanical and enzymatic dissociation using the Miltenyi Biotec Mouse tumor dissociation kit. Antibody-conjugated magnetic beads removed myelinated and mouse cells in a two column-based filtration step. Sample size for shotgun lipidomics, RNA-seq, and WES were determined based on model availability. For in vivo studies involving genetically manipulated tumors (see Fig. 6), sample size was estimated based on previous studies.
METHOD DETAILS
Secreted Gaussia luciferase measurements
Cells were infected with a lentiviral vector containing a secreted Gaussia luciferase (sGluc)-encoding reporter gene (Targeting Systems no. GL-GFP) and intracranially implanted into the right striatum in mice (4 × 105 cells/mouse). To measure the levels of sGluc, 6 μL of blood was collected from the tail vein and immediately mixed with 50 mM EDTA to prevent coagulation. sGluc activity was obtained by measuring chemiluminescence after injection of 100 μL of 100 μM coelentarazine (Nanolight) in a 96-well plate, as described before.45
Genetic manipulations
Lentivirus particles used for genetic manipulation were produced by transfection of 293-FT cells (Thermo) with Lipofectamine 2000 (Invitrogen). Virus particles were collected 48 h after transfection. The CDKN2A shRNA guide sequences used to knockdown p14 and p16 TRCN0000255853 and TRCN0000255849 in the Vector backbone pLKO (Addgene, plasmid 10878). Mission pLKO.1-puro Non-Mammalian shRNA Control (Sigma-Aldrich, SHC002) was used as the shRNA scramble control (shC). GPX4 CRISPR-Cas9 gene disruption was performed with the LentiCRISPR V2 vector (Addgene, plasmid 52961) using the following oligonucleotide sequences:
| Forward (5’ > 3’) | Reverse (5’ > 3’) | |
| sgGPX4-1 | caccgttaacctggacaagtaccgg | aaacccggtacttgtccaggttaac |
| sgGPX4-2 | caccgtgcacgagttttccgccaag | aaaccttggcggaaaactcgtgcac |
| sgGPX4-3 | caccgaagtaaactacactcagctc | aaacgagctgagtgtagtttacttc |
Cells were spinfected with lentivirus generated from LentiCRISPR V2 non-targeting (NT) control, sgGPX4-1, sgGPX4-2, or sgGPX4-3 (500 uL lentivirus, 1 million cells plated in 1 mL media in a 6-well plate, 1 ug/mL Polybrene) at 800g for 1 hour and 30 minutes (32C). Cells were immediately transferred into normal gliomasphere media and subjected to puromycin selection (1μM puromycin). Post-spinfection, cells were continuously cultured in 1μM ferrostatin-1 (refreshed every 2–3 days) until experimentation or implantation into mice.
CRISPR-Cas9 mediated CDKN2A disruption, as well as addback of p16, were performed as previously described.39 Cells were infected with pEN35-CDKN2A-Ex2-R, pEN35-CDKN2A-Ex2-L, and pHDR-CDKN2A-Ex2-CMV-EGFP cDNA (1 ug of cDNA per construct per 1,000,000 cells). Cells were allowed to expand and were then sorted via FACS isolation. After validation of CDKN2A knock-out with Western blotting, cells were then infected with pHIV-INK4A-mOrange2 cDNA (1 ug of cDNA per 1,000,000 cells), expanded, and then sorted via FACS isolation.
Immunoblotting
Cells were collected and lysed in RIPA buffer (Boston BioProducts) containing Halt Protease and Phosphatase Inhibitor (Thermo Fisher Scientific). Lysates were centrifuged at 14,000 × g for 15 min at 4°C. Protein samples were then boiled in NuPAGE LDS Sample Buffer (Invitrogen) and NuPAGE Sample Reducing Agent (Invitrogen), separated with SDS–PAGE on 12% Bis-Tris gels (Invitrogen) and transferred to nitrocellulose membranes (GE Healthcare). Immunoblotting was performed per the antibody manufacturers’ specifications, as mentioned previously. Membranes were developed with the SuperSignal system (Thermo Fisher Scientific).
Exogenous fatty acid addition – shotgun lipidomics
For Figure 4A, cells were plated at a density of 5 million cells in 50 mL media and treated with the following conditions (in triplicate) for 24 hours prior to collection for shotgun lipidomics analysis: Vehicle (DMSO), 75 μM arachidonic acid (AA), or 20 μM DGAT1 inhibitor+10 μM DGAT2 inhibitor+75 μM AA for 24 hours.
Shotgun lipidomics
After removal of cell culture media, cells were resuspended in 200 uL of ice-cold PBS and collected in glass tubes. Samples were prepared for shotgun lipidomic analysis at the UCLA Lipidomics Core (detailed protocol published by Hsieh W.Y. et al).46 In short, a modified Bligh and Dyer extraction was performed on all samples.47 Prior to extraction, a 13-class lipid class Lipidizer Internal Standard (AB Sciex, 5040156) was added to each sample. After two extractions, the pooled organic layers were dried down in a Genevac EZ-2 Elite, and lipid samples were resuspended in 1:1 methanol/dichloromethane with 10 mM Ammonium Acetate. Samples were analyzed on the Sciex Lipidyzer Platform for targeted quantitative measurement of 1100 lipid species across 13 classes. Differential Mobility Device on Lipidyzer was tuned with SelexION tuning kit (Sciex 5040141). Instrument settings, tuning settings, and MRM list available upon request. Data analysis was performed according to the workflow published by Su, B. et al (2021).48
Annotating lipid species characteristics
For results examining acyl tail characteristics of lipid species, monoacyl species were annotated based off their one fatty acyl tail (e.g., LPC 14:0 – carbon number: 14; number of double bonds: 0). Diacyl species were annotated based off their two fatty acyl tails (e.g., DG18:0_18:3 counted as two different acyl tails, one with 18 carbons and 0 double bonds, and another with 18 carbons and 3 double bonds). Triacyl species are fragmented into multiple species, with each species representing a different acyl tail, and were annotated based off their one detectable fatty acyl tail per species (e.g., TG 49:2-16:1 counted as having one acyl tail with 16 carbons and 1 double bond). For results examining lipids at the molecular species level, TAG length and double bond numbers were based on the single detected tail and saturation was determined based on information on both the number of double bonds on the detected tail and total number of double bonds detected across all three tails (e.g., TG52:5-16:0 counted as PUFA because there must be at least two double bonds on one of the other acyl tails). For results examining TAG species characteristics, double bonds and total carbons were determined by looking at the overall number of double bonds and overall number of carbons present in the TAG species.
Fatty acid analysis
Fatty Acid analysis was carried out as described previously.49 Briefly, cells were collected and counted for normalization. Nonadecanoic methyl ester internal standard (N-19-M, Nu-Chek Prep) was added to derivatization mixture. Total fatty acids methyl esters were prepared by acid methanolysis reaction carried out for 12–16 hours at 45C. Resulting FAMEs were extracted in hexane. Samples were analyzed with an Agilent 7890B/5977A GC-MS using a 30m UI-DB-WAX column (Agilent 122-7032UI). GC oven and MSD settings available upon request. Fatty acids were identified by retention time and parental ion m/z. Relative amounts were normalized to internal standard.
Quantifying lipid peroxidation with BODIPY 581/591 C11
Lipid peroxidation was assessed using the lipid peroxidation probe C11 BODIPY 581/591 from Thermo Fisher (D3861). To assess lipid peroxidation in vitro, gliomaspheres (non-GFP transduced) were split and cultured for 24 h in previously described gliomasphere medium except the B-27 supplement is without antioxidants (Gibco, 10889038). After 24 h, the cells were collected, dissociated with trypsin, and re-suspended as single cells in Hank’s Buffered Salt Solution (HBSS) containing the C11 BODIPY 581/591 (1:1000) and Hoechst nuclear stain (1:1000). Cells were incubated for 15 minutes in the dark in solution on a Cell-tak treated coverslip to enable cell adhesion. Coverslips were then washed once in HBSS and inverted onto a slide with 100 μl HBSS and parafilm spacers, and then sealed with Valap to enable live cell confocal microscopy. Images were acquired on a confocal LSM880 microscope at 63x magnification to capture DAPI (405/462), non-oxidized C11 (594/667), and oxidized C11 (488/554). Quantification was performed as follows: .czi files were opened in ImageJ FIJI and the red and green channels were merged. ROIs were drawn around 30 representative cells, and values were generated on a cell-by-cell basis for both red and green channels. Oxidation was quantified as (green average intensity)/ (red average intensity+green average intensity). Representative images were generated by merging red, green, and blue channels in ImageJ and adjusting brightness/contrast for visualization (same parameters applied to all images being compared to one another).
Cell death kinetics – Incucyte
Opaque-sided clear bottomed 96-well plates or 384-well plates were laminin-coated for 1 hour and then washed two times with neurobasal media. Gliomaspheres were dissociated into single cell suspension, filtered, plated at 20,000 cells per well (96-well) or 5,000 cells per well (384-well), and allowed to attach for 4 to 12 h. Once cells were attached, media was briefly removed and replaced with medium containing antioxidant-free B27 and CellTOX cell death green fluorescent dye (1:1000) (Promega, G8741) and nuclear red fluorescent dye (1:1000) (Sartorius, 4632) as well as the required compounds. The plate was then set in an Incucyte live-cell imaging system (Sartorius) which acquires images every 2 h for up to five days. Incucyte software image analysis was performed using top-hat image normalization (Satorius Incucyte User Manual), and percent cell death was calculated using percentage of green objects/mm2 over red objects/mm2.
Cell death kinetics – Trypan blue exclusion
Gliomaspheres were dissociated into single cell suspension, filtered, and plated at 200,000 cells per well in a 12-well plate. Cells plated in triplicate per condition in DMEM/F12+antioxidant free B27+glutamax+Pen/strep with either DMSO vehicle or 2.5μM RSL3. After 72 hours, cells were collected, split with TrypLE for 3 minutes, and resuspended in 1 mL media. Cells were then counted with the Countess 3 Automated Cell Counter at a 1:1 ratio of trypan blue: cell suspension. Brightfield was used to capture total cell count, Dead Cells (%), and Live Cells (%).
Quantifying MDA levels with immunofluorescent microscopy - Tissue
Tumor tissues were fixed in 4% PFA and paraffin-embedded, and slices were mounted onto slides. Slides were incubated at 60C for 30 min, then deparaffinized in xylene (100%, 3×5 min). Slides were then rehydrated in graded alcohols (100%, 100%, 95%, 80%, 70%, 50%, 2 min each). Antigen retrieval was performed in antigen retrieval solution (100 mM Tris, 5% (w/v) urea, pH 9.5). Blocking was done with 10% BSA in TBST for 5 minutes. Slides were incubated overnight at 4°C with 1:500 goat anti-MDA primary antibody (Abcam ab27644) followed by secondary incubation with 1:1000 donkey anti-goat Alexa 647 (Abcam ab150131). Coverslips were mounted to the slides using the ProLong Gold antifade reagent with DAPI (Invitrogen, P36935). Confocal microscopy on the LSM880 microscope at 20x magnification was used to capture tissue sections (MDA (633/697), GFP (488/546), and DAPI (405/498)). Average intensity of MDA was quantified in ImageJ Fiji, and images were converted to JPEG files using the Zen Blue software. To quantify MDA intensity, the threshold was set based off of a representative, high MDA image. This threshold was then used for all images being compared to one another. For image visualization, the Alexa-647 far red channel was changed to orange.
Quantifying MDA levels with immunofluorescent microscopy - Cells
Cells were plated at a density of 200K per well in a 12-well plate and incubated for 24 hours in DMEM/F12+glutamax+P/S+B27 without antioxidants. Cells were collected, pelleted, washed with PBS, pelleted, then resuspended in 100 uL PBS and added to Cell-Tak’d coverslips. After 1 hour, cells were fixed with 4% PFA and incubated for an additional 15 minutes. PFA was washed three times with PBS and were permeabilized with 0.5% Triton X-100 in PBS. After 15 minutes, cells were washed three times with PBS and were blocked with 10% BSA in PBS for an hour. After blocking, cells were washed three times with PBS and incubated in primary antibody overnight at 4C (goat anti-MDA primary antibody diluted 1:500 in 1% BSA in PBS). The next day, cells were washed three times with PBS to remove primary antibody. Cells were then incubated for 1 hour in the dark with secondary antibody (donkey anti-goat Alexa-647 secondary antibody diluted 1:1000 in 1% BSA in PBS). Cells were washed three times with PBS to remove excess secondary, washed one time with distilled water to remove PBS, then mounted to slides with ProLong Gold antifade reagent with DAPI. Confocal microscopy on the LSM880 microscope at 60x magnification was used to capture at least 30 cells (MDA (633/697), GFP (488/546), and DAPI (405/498)). Average intensity of MDA was quantified in ImageJ Fiji, and images were converted to JPEG files using the Zen Blue software. To quantify MDA intensity, the threshold was set based off a representative, high MDA image. This threshold was then used for all images being compared to one another. For image visualization, the Alexa-647 far red channel was changed to orange.
Lipid droplet imaging and quantification
Lipid droplets were assessed using HCS LipidTOX™ Deep Red neutral lipid dye. Cells were collected and trypsinized following standard procedure for gliomasphere cultures. The cell pellet was then re-suspended in media containing the LipidTOX™ dye (1:1000). Cells were distributed onto a coverslip pre-treated with Cell-tak to enable rapid cell adhesion. Cells were allowed to attach to the coverslip for 30 min before the addition of 4% PFA to fix the cells. Coverslips were washed in HBSS and then water and mounted onto a slide using Invitrogen Prolong Gold Antifade mounting reagent with DAPI. Confocal microscopy on the LSM880 microscope at 63x magnification was used to capture the lipid droplet (488/563) and DAPI (405/453). Lipid droplet area and number per cell was quantified using ImageJ FIJI. To quantify lipid droplet area, the threshold was set based off of a representative, high LD image (only masking punctate signal). This threshold was then used for all images being compared to one another. Particles were then analyzed, summed per image, then divided by the number of cells (at least 30 cells) captured in that image. For image visualization, the color of the far-red channel was switched to green on ImageJ, and the far red and blue channels were merged into a composite image.
Whole exome and RNA sequencing
Prior to January 21st 2021, whole exome libraries were prepared using the SeqCap EZ Human Exome Library V3 kit and sequenced on the Illumina HiSeq 3000 or Novaseq 6000 systems (paired-end, 150bp). Afterwards, libraries were constructed using the KAPA HyperExome enrichment kit and sequenced on the NovaSeq 6000 system. On average, 120-fold (100~150 fold) exome-wide target coverage was achieved for all sequenced tumor and normal blood and peripheral blood mononuclear cell (PBMC) samples. Transcriptome libraries were constructed from poly(A) selected RNA using the Nugen Universal Plus mRNA-seq library prep kit and sequenced on the HiSeq 3000 or Novaseq 6000 systems (paired-end, 150bp). Gliomas were RNA sequenced and achieved an average read depth of approximately 35 million mapped deduplicated reads per sample.
Whole exome QC and mapping
For WES, adapter trimming was performed with Cutadapt,50 with minimum read length after trimming 20 and mean quality value before trimming of 20. To remove potential mouse contamination from sequenced library, we mapped trimmed reads to mouse (mm10) and human (GRCh38) genomes simultaneously with BBsplit from the BBtools package.51 After removing contaminant and ambiguous reads, alignment to the hg38 human reference genome was carried out with BWA 0.7.17-r1188.52 We followed GATK v4.2.0.0 best practices to mark duplicate reads and perform base quality score recalibration with ApplyBQSR.53
RNA-seq QC and mapping
For RNA sequencing, reads unambiguously aligning to the mouse genome were removed with the same mapping strategy as for WES, using Seal from the BBtools package. Filtered reads were then processed through the UCSC Toil RNA sequencing pipeline54 for quality control, adapter trimming, sequence alignment, and expression quantification. Briefly, sequence adapters are trimmed and sequences are then aligned to human reference genome (GRCh38) using STAR v2.4.2a55 and gene expression quantification is performed using RSEM v1.2.25.56 Expression values are output as raw counts and transcripts per million (TPM) for downstream analysis.
Somatic mutation and copy number alteration calling
For tumor samples with a sequenced matched normal, Mutect2 v4.2.0.057 MuSE v1.0rc,58 and Varscan259 were used to call mutations. To avoid potential false positive variants, variants with minimum coverage of 20 reads and identified by at least two mutation callers were selected for further analysis. For samples without matched normal sequencing, variants were called using Mutect2 in tumor-only mode. Variants were compared against the matched normal sample when available, as well as a constructed panel of normal samples following GATK best practices. Variants were then filtered based on occurrence frequency in COSMIC database60 (>7 occurrences in CNS tumors) or annotated as “likely oncogenic” or “oncogenic” by OncoKB.61 CNVkit was used to detect copy number changes from whole-exome sequencing data.62
EGFRvIII calls were generated based on alternative transcript splicing of EGFR detected in RNA sequencing data. EGFRvIII calls were determined by calculating the fraction of reads mapping to the junction between exons 7 and 8 and the aberrant junction between exons 1 and 8 (indicative of an in-frame deletion of exons 2–7) to generate EGFRvIII transcript allele frequencies (TAFs). Samples with TAFs > 10% were considered positive for EGFRvIII variant.
Integrated CDKN2A status based on whole exome and RNA sequencing
To generate a final call of CDKN2A status for all patient tumors, we combined information from our whole exome and RNA sequencing pipelines. First, segment level copy number values were manually reviewed for consistency with copy number ratios calculated within each 5 kilobase window overlapping CDKN2A. Cases with discordant bin-level and segment-level copy number estimates were manually adjusted based on the average of bins. Patient tumor samples with CDKN2A deep deletions based on WES were typed as CDKN2A null. Samples with CDKN2A shallow deletions or diploid copy number based on WES combined with greatly decreased CDKN2A RNA expression (log2 TPM < 2.5) and copy number deletions detected in their matched patient-derived models were typed as CDKN2A null. All other samples with shallow deletions, diploid copy number, or copy number gains of CDKN2A based on WES without further evidence of deletion were CDKN2A wild type.
Lipidome-transcriptome correlation analysis
For patient tumors, matched lipidomic profiling (% of total, n=918 species overlapping in patient and gliomasphere profiling) and gene expression data (log2 TPM+1) from 84 patient tumor samples were integrated to evaluate relationships between variation in gene expression and lipid composition. Gene expression data were filtered to genes annotated in a Gene Ontology, KEGG, or Reactome term related to lipid metabolism with average expression above 0.5 log2 TPM+1 (n = 999).63–65 Pearson’s correlation coefficients were then calculated for each lipid-gene pair. Correlation values were hierarchically clustered based on Euclidean distance using Ward’s method. Clusters of lipids and genes were defined using static tree cutting of hierarchical clustering results.
For gliomaspheres, lipidomic profiling (% of lipid class, n=985 species) and gene expression data (log2 TPM+1) from 43 samples were used. Gene expression data were first filtered to protein coding genes with average expression above 0.5 log2 TPM+1. Genes and lipids were further filtered to remove the bottom 20% and 10% of low variance features respectively (10,778 remaining genes, 835 remaining lipid species). Pearson’s correlation coefficients were then calculated and the top 25% of genes (n=2,695 genes) with most variable correlation values were selected for downstream analysis and plotting based on their ability to preserve patterns observed globally.
Correlation cluster annotation and lipid-gene cluster co-enrichment
Lipid clusters defined in the patient lipid-gene correlation analysis were annotated based on their enrichment for specific lipid classes, saturation levels, and chain length categories based on hypergeometric p-values calculated through comparison of each cluster to the background of all analyzed lipid species. Hypergeometric p-values for each category were adjusted using the Benjamini-Hochberg procedure. Gene clusters were annotated with their enriched lipid metabolic features based on geneset over-representation analysis run on Gene Ontology Biological Process (GO BP) terms using the PANTHER Classification System.66 GO BP results were filtered to terms involved in the metabolism or biosynthesis of detectable lipid classes (i.e. excluding terms including “negative regulation” or “catabolism”) with adjusted p-value < 0.05 and then sorted based on fold enrichment to identify enriched lipid metabolic terms.
Co-enrichment scores quantifying the relationships between pairs of lipid and gene clusters were calculated based on the average of two two-sample Kolmogorov-Smirnov statistics: one statistic calculated on genes ranked based on their correlation with lipids of a given lipid cluster, and another statistic calculated on lipids ranked based on their correlation with genes of a given gene cluster. The first statistic quantifies the specificity of gene cluster correlations for each lipid cluster while the second statistic quantifies the specificity of lipid cluster correlations for each gene cluster. P-values were adjusted using the Benjamini-Hochberg procedure and cluster pairs with adjusted p-values < 0.05 for both gene-wise and lipid-wise enrichment as well as congruent directionality (positive or negative enrichment) were deemed significant co-enrichment results.
Genetic alteration - lipid cluster enrichment scoring
To quantify whether lipid species in each lipid cluster are enriched positively or negatively with a given genetic alteration, a log odds ratio was calculated between the fraction of significantly over-abundant lipid species within a lipid cluster and the fraction of significantly under-abundant lipid species within the same cluster (Student’s t-test p-value < 0.05). Fisher’s exact test, with Haldane-Anscombe correction (adding 0.5 to each cell to avoid zeroes), was carried out on counts of the significantly different lipid species inside or outside of a given lipid cluster for a given genetic alteration. This test was performed for each pair of lipid clusters and genetic alterations and p-values were adjusted using the Benjamini-Hochberg procedure.
Gene cluster scoring and genomic association testing
Single sample GSEA (ssGSEA)67 was used to score individual samples for their enrichment for gene sets defined in correlation analyses in patient tumors and gliomaspheres. For gliomaspheres, a composite score was defined as the difference between the enrichment scores for the two gene clusters and was used to test for associations with genomic alterations.
For genome-wide testing, gene-level copy number amplifications (8+ copies), deep deletions (< 0.5 copies), and somatic mutations detected via WES were individually considered for testing. For each gene, samples were grouped by their binary alteration status and Welch’s t-test was used to quantify their association with ssGSEA composite scores. Gene-level results for somatic amplifications and deletions were then clustered into locus-level results by grouping genes with genomic coordinates within 1 Megabase of each other and taking the most extreme test statistic from the set of clustered genes to represent the locus. Locus-level results for amplifications and deletions were then merged with gene-level results for point mutations and p-values were adjusted using the Benjamini-Hochberg procedure.
Differential expression analysis
Differential gene expression analysis was carried out using DESeq2 on expected count data from RSEM rounded to whole numbers.68 Genes were pre-filtered to remove genes with average counts below 100. Fold change estimates were shrunken using the adaptive shrinkage option lfcShrink(…, type = “ashr”) provided in DESeq2.68
Principal component analysis and Varimax rotation
Principal component analysis (PCA) was carried out using the prcomp() function in R v4.1.0 without scaling features.69 PCA loadings for PC1 and PC2 were then rotated using varimax rotation without Kaiser normalization using the R package varimax.27 Rotated component scores were calculated by multiplying the original PC scores by the varimax rotation matrix.
Differential lipid abundance permutation analysis
To determine the statistical significance of differential lipid abundance results obtained using Student’s t-test, permutation analysis was performed by randomly shuffling the group labels without replacement and re-running the differential lipid abundance testing 100,000 times. The number of differentially abundant lipids were then counted for various significance levels (α = [0.25, 0.15, 0.05, 0.01, 0.001]) and quantiles of the permutation results (median, top 10%, top 5%, top 1%) for comparison to the observed number of differentially abundant lipids at the same significance levels. A final permutation p-value was calculated as the fraction of random permutations with equal or more significantly different lipid species than the observed result at α = 0.05.
TCGA GBM multi-omics data and analysis
Somatic point mutations and copy number alterations for TCGA GBM samples was accessed through CBioPortal (Glioblastoma, Cell 2013).13,70,71 Gene expression results from RSEM were accessed from the Toil RNA-seq Recompute72 project publicly available through UCSC XenaBrowser.73 Datasets were then filtered to only include samples present in all datasets (n = 142).
Samples were then scored for their enrichment for each gene cluster defined in the patient correlation analysis using ssGSEA. ssGSEA scores were then stratified by CDKN2A alteration status (defined as either deletion or somatic mutation with loss of heterozygosity) and statistical testing was performed using Welch’s t-test for both our cohort and TCGA GBM.
Gene expression subtyping
Tumors from patients with GBM, orthotopic xenografts, and gliomasphere cultures were scored for their enrichment of previously defined TCGA GBM subtypes using ssGSEA.74 An empirical distribution of expected ssGSEA scores for each subtype signature were generated from permutations of the patient tumor expression values generating 1 million simulated ssGSEA scores per subtype. These simulated distributions were used to convert observed ssGSEA scores into permutation z-scores and empirical p-values indicating significance of subtype enrichment as done previously.74 In cases where a sample showed enrichment for multiple subtypes, the sample was called based on its maximum z-score.
QUANTIFICATION AND STATISTICAL ANALYSIS
Comparisons were made with two-tailed unpaired Student’s t-tests, and unless otherwise stated, unadjusted p-values <0.05 were considered statistically significant. Welch’s t-test was used when comparison groups were observed to have disparate variances. No further statistical methods or tests assessing the distribution of the data were used. All data from multiple independent experiments were assumed to be t-distributed. For each experiment, replicates are as noted in the figure legends. Data represent mean ± s.d. values unless otherwise indicated in the figure legends. All statistical analyses were calculated in Prism 9.0 (GraphPad). Sample size was selected based on previous studies. No samples were excluded. For all figures: p >0.05=ns; p <0.05 = *; p<0.01 = **; p<0.001 = ***; p<0.0001 = ****.
Supplementary Material
Supplemental Table Titles:
Table S1: GO biological processes results for gene clusters generated from correlation matrix of lipid species across tumors from patients with GBM. Related to Figure 1C.
Highlights (<85 characters including spaces, 3–4 bullet points).
Molecular and lipidomic resource of over 150 GBM tumors and derivative models
Unbiased multi-omic analysis identifies CDKN2A as a regulator of lipid metabolism
CDKN2A deletion reduces oxidizable PUFA sequestration into lipid droplets
GBM with CDKN2A deletion are susceptible to lipid peroxidation and ferroptosis
ACKNOWELDGEMENTS
The authors thank Dr. Harvey Herschman, Dimitri Cadet, Dr. Jonathan Tsang, and other past and present members of the Nathanson lab for their valuable critiques and discussions. We thank Dr. Mikhail Shchepiniov for deuterated PUFAs. We further acknowledge Dr. Michaela Veliova for feedback regarding experimental planning. We thank the Brain Tumor Translational Resource (BTTR) and Dr. Xinmin Li with the Center for Genomics and Bioinformatics (TCGB) for recruitment of tumor specimens and sequencing, respectively. All graphics created with Biorender (Biorender.com). Finally, we are grateful to the patients with brain cancer and their families who consented to donating tumor tissue for this research. This work was funded by NIH grants R01NS121319 (DAN, SJB, SJD), P50CA211015 (DAN, TFC, SJB, TGG, LML), P50CA211015 (DAN, TFC, SJB, TGG, LML), R01CA227089 (DAN, TGG, TFC, LML), P01HL146358 (SJB) and the Department of Defense grant W81XWH-20-1-0453 (JKM, DM, DAN).
Inclusion and Diversity
We support inclusive, diverse and equitable conduct of research.
Footnotes
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
DECLARATION OF INTERESTS
D.A.N. is a co-founder of Trethera Corporation and has equity in the company. D.A.N. and T.F.C. are co-founders of Katmai Pharmaceuticals and have equity in the company. S.J.D. is co-founder of Prothegen Inc., serves on the scientific advisory board of and Hillstream BioPharma, and holds patents related to ferroptosis. T.G.G. has consulting and equity agreements with Auron Therapeutics, Boundless Bio, Coherus BioSciences and Trethera Corporation.
References:
- 1.Ying H, Kimmelman AC, Lyssiotis CA, Hua S, Chu GC, Fletcher-Sananikone E, Locasale JW, Son J, Zhang H, Coloff JL, et al. (2012). Oncogenic Kras maintains pancreatic tumors through regulation of anabolic glucose metabolism. Cell 149, 656–670. 10.1016/J.CELL.2012.01.058. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 2.Garcia-Cao I, Song MS, Hobbs RM, Laurent G, Giorgi C, De Boer VCJ, Anastasiou D, Ito K, Sasaki AT, Rameh L, et al. (2012). Systemic Elevation of PTEN Induces a Tumor-Suppressive Metabolic State. Cell 149, 49–62. 10.1016/J.CELL.2012.02.030. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 3.Dang L, White DW, Gross S, Bennett BD, Bittinger MA, Driggers EM, Fantin VR, Jang HG, Jin S, Keenan MC, et al. (2009). Cancer-associated IDH1 mutations produce 2-hydroxyglutarate. Nature 2009 462:7274 462, 739–744. 10.1038/nature08617. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 4.Maddocks ODK, Berkers CR, Mason SM, Zheng L, Blyth K, Gottlieb E, and Vousden KH (2012). Serine starvation induces stress and p53-dependent metabolic remodelling in cancer cells. Nature 2012 493:7433 493, 542–546. 10.1038/nature11743. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 5.Kim J, Lee HM, Cai F, Ko B, Yang C, Lieu EL, Muhammad N, Rhyne S, Li K, Haloul M, et al. (2020). The hexosamine biosynthesis pathway is a targetable liability in KRAS/LKB1 mutant lung cancer. Nature Metabolism 2020 2:12 2, 1401–1412. 10.1038/s42255-020-00316-0. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 6.Li H, Ning S, Ghandi M, Kryukov GV, Gopal S, Deik A, Souza A, Pierce K, Keskula P, Hernandez D, et al. (2019). The landscape of cancer cell line metabolism. Nature Medicine 2019 25:5 25, 850–860. 10.1038/s41591-019-0404-8. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 7.Chen PH, Cai L, Huffman K, Yang C, Kim J, Faubert B, Boroughs L, Ko B, Sudderth J, McMillan EA, et al. (2019). Metabolic Diversity in Human Non-Small Cell Lung Cancer Cells. Mol Cell 76, 838–851.e5. 10.1016/J.MOLCEL.2019.08.028. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 8.Domcke S, Sinha R, Levine DA, Sander C, and Schultz N (2013). Evaluating cell lines as tumour models by comparison of genomic profiles. Nature Communications 2013 4:1 4, 1–10. 10.1038/ncomms3126. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 9.Lee J, Kotliarova S, Kotliarov Y, Li A, Su Q, Donin NM, Pastorino S, Purow BW, Christopher N, Zhang W, et al. (2006). Tumor stem cells derived from glioblastomas cultured in bFGF and EGF more closely mirror the phenotype and genotype of primary tumors than do serum-cultured cell lines. Cancer Cell 9, 391–403. 10.1016/J.CCR.2006.03.030. [DOI] [PubMed] [Google Scholar]
- 10.Marin-Valencia I, Yang C, Mashimo T, Cho S, Baek H, Yang XL, Rajagopalan KN, Maddie M, Vemireddy V, Zhao Z, et al. (2012). Analysis of tumor metabolism reveals mitochondrial glucose oxidation in genetically diverse human glioblastomas in the mouse brain in vivo. Cell Metab 15, 827–837. 10.1016/J.CMET.2012.05.001. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 11.Davidson SM, Papagiannakopoulos T, Olenchock BA, Heyman JE, Keibler MA, Luengo A, Bauer MR, Jha AK, O’Brien JP, Pierce KA, et al. (2016). Environment impacts the metabolic dependencies of ras-driven non-small cell lung cancer. Cell Metab 23, 517–528. 10.1016/J.CMET.2016.01.007. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 12.Barekatain Y, Ackroyd JJ, Yan VC, Khadka S, Wang L, Chen KC, Poral AH, Tran T, Georgiou DK, Arthur K, et al. (2021). Homozygous MTAP deletion in primary human glioblastoma is not associated with elevation of methylthioadenosine. Nature Communications 2021 12:1 12, 1–13. 10.1038/s41467-021-24240-3. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 13.Brennan CW, Verhaak RGW, McKenna A, Campos B, Noushmehr H, Salama SR, Zheng S, Chakravarty D, Sanborn JZ, Berman SH, et al. (2013). The somatic genomic landscape of glioblastoma. Cell 155, 462. 10.1016/j.cell.2013.09.034. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 14.Verhaak RGW, Hoadley KA, Purdom E, Wang V, Qi Y, Wilkerson MD, Miller CR, Ding L, Golub T, Mesirov JP, et al. (2010). Integrated Genomic Analysis Identifies Clinically Relevant Subtypes of Glioblastoma Characterized by Abnormalities in PDGFRA, IDH1, EGFR, and NF1. Cancer Cell 17, 98–110. 10.1016/J.CCR.2009.12.020. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 15.Guo D, Prins RM, Dang J, Kuga D, Iwanami A, Soto H, Lin KY, Huang TT, Akhavan D, Hock MB, et al. (2009). EGFR signaling through an Akt-SREBP-1-dependent, rapamycin-resistant pathway sensitizes glioblastomas to antilipogenic therapy. Sci Signal 2. 10.1126/SCISIGNAL.2000446/SUPPL_FILE/2_RA82_SM.PDF. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 16.Gimple RC, Kidwell RL, Kim LJY, Sun T, Gromovsky AD, Wu Q, Wolf M, Lv D, Bhargava S, Jiang L, et al. (2019). Glioma Stem Cell-Specific Superenhancer Promotes Polyunsaturated Fatty-Acid Synthesis to Support EGFR Signaling. Cancer Discov 9, 1248–1267. 10.1158/2159-8290.CD-19-0061. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 17.Williams KJ, Argus JP, Zhu Y, Wilks MQ, Marbois BN, York AG, Kidani Y, Pourzia AL, Akhavan D, Lisiero DN, et al. (2013). An essential requirement for the SCAP/SREBP signaling axis to protect cancer cells from lipotoxicity. Cancer Res 73, 2850–2862. 10.1158/0008-5472.CAN-13-0382-T/651001/AM/AN-ESSENTIAL-REQUIREMENT-FOR-THE-SCAP-SREBP. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 18.Lita A, Pliss A, Kuzmin A, Yamasaki T, Zhang L, Dowdy T, Burks C, de Val N, Celiku O, Ruiz-Rodado V, et al. (2021). IDH1 mutations induce organelle defects via dysregulated phospholipids. Nature Communications 2021 12:1 12, 1–16. 10.1038/s41467-020-20752-6. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 19.Fack F, Tardito S, Hochart G, Oudin A, Zheng L, Fritah S, Golebiewska A, Nazarov PV, Bernard A, Hau A-C, et al. Altered metabolic landscape in IDH-mutant gliomas affects phospholipid, energy, and oxidative stress pathways. 10.15252/emmm.201707729. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 20.Jain JL et al. (1979). Fundamentals of Biochemistry. In (S. Chand & Company; ), pp. 594–639. [Google Scholar]
- 21.Freed-Pastor WA, Mizuno H, Zhao X, Langerød A, Moon SH, Rodriguez-Barrueco R, Barsotti A, Chicas A, Li W, Polotskaia A, et al. (2012). Mutant p53 disrupts mammary tissue architecture via the mevalonate pathway. Cell 148, 244–258. 10.1016/j.cell.2011.12.017. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 22.Bi J, Ichu T-A, Zanca C, Furnari FB, Cravatt BF, and Correspondence PSM (2019). Oncogene Amplification in Growth Factor Signaling Pathways Renders Cancers Dependent on Membrane Lipid Remodeling. Cell Metab 30, 525–538.e8. 10.1016/j.cmet.2019.06.014. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 23.Muranaka H, Hayashi A, Minami K, Kitajima S, Kohno S, Nishimoto Y, Nagatani N, Suzuki M, Kulathunga LAN, Sasaki N, et al. (2017). A distinct function of the retinoblastoma protein in the control of lipid composition identified by lipidomic profiling. Oncogenesis 2017 6:6 6, e350–e350. 10.1038/oncsis.2017.51. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 24.Liu Y, He Y, Jin A, Tikunov AP, Zhou L, Tollini LA, Leslie P, Kim TH, Li LO, Coleman RA, et al. (2014). Ribosomal protein-Mdm2–p53 pathway coordinates nutrient stress with lipid metabolism by regulating MCD and promoting fatty acid oxidation. Proc Natl Acad Sci U S A 111, E2414–E2422. 10.1073/PNAS.1315605111/SUPPL_FILE/PNAS.201315605SI.PDF. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 25.Bates S, Phillips AC, Clark PA, Stott F, Peters G, Ludwig RL, and Vousden KH (1998). p14ARF links the tumour suppressors RB and p53. Nature 1998 395:6698 395, 124–125. 10.1038/25867. [DOI] [PubMed] [Google Scholar]
- 26.Zhao R, Choi BY, Lee MH, Bode AM, and Dong Z (2016). Implications of Genetic and Epigenetic Alterations of CDKN2A (p16INK4a) in Cancer. EBioMedicine 8, 30–39. 10.1016/j.ebiom.2016.04.017. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 27.Kaiser HF (1958). The varimax criterion for analytic rotation in factor analysis. Psychometrika 23, 187–200. 10.1007/BF02289233/METRICS. [DOI] [Google Scholar]
- 28.Ayala A, Muñoz MF, and Argüelles S (2014). Lipid Peroxidation: Production, Metabolism, and Signaling Mechanisms of Malondialdehyde and 4-Hydroxy-2-Nonenal. Oxid Med Cell Longev 2014. 10.1155/2014/360438. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 29.Bailey AP, Koster G, Lechene CP, Postle AD, and Gould Correspondence AP (2015). Antioxidant Role for Lipid Droplets in a Stem Cell Niche of Drosophila. Cell 163, 340–353. 10.1016/j.cell.2015.09.020. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 30.Dixon SJ, Lemberg KM, Lamprecht MR, Skouta R, Zaitsev EM, Gleason CE, Patel DN, Bauer AJ, Cantley AM, Yang WS, et al. (2012). Ferroptosis: An Iron-Dependent Form of Nonapoptotic Cell Death. Cell 149, 1060–1072. 10.1016/J.CELL.2012.03.042. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 31.Eaton JK, Furst L, Ruberto RA, Moosmayer D, Hilpmann A, Ryan MJ, Zimmermann K, Cai LL, Niehues M, Badock V, et al. (2020). Selective covalent targeting of GPX4 using masked nitrile-oxide electrophiles. Nature Chemical Biology 2020 16:5 16, 497–506. 10.1038/s41589-020-0501-5. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 32.Dixon SJ, Winter GE, Musavi LS, Lee ED, Snijder B, Rebsamen M, Superti-Furga G, and Stockwell BR (2015). Human Haploid Cell Genetics Reveals Roles for Lipid Metabolism Genes in Nonapoptotic Cell Death. ACS Chem Biol 10, 1604–1609. 10.1021/ACSCHEMBIO.5B00245/SUPPL_FILE/CB5B00245_SI_002.PDF. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 33.Miotto G, Rossetto M, Roveri A, Venerando R, Vučković A-M, Di Paolo ML, Bosello-Travain V, Zaccarin M, Maiorino M, Toppo S, et al. (2018). Insight the mechanism of ferroptosis inhibition by ferrostatin-1. Free Radic Biol Med 120, S120–S121. 10.1016/J.FREERADBIOMED.2018.04.397. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 34.Yang WS, Kim KJ, Gaschler MM, Patel M, Shchepinov MS, and Stockwell BR (2016). Peroxidation of polyunsaturated fatty acids by lipoxygenases drives ferroptosis. Proc Natl Acad Sci U S A 113, E4966–E4975. 10.1073/PNAS.1603244113/SUPPL_FILE/PNAS.1603244113.SD01.XLSX. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 35.Kagan valerian, Mao G, Qu F, pedro Friedmann angeli J, doll sebastian, st Croix C, Hussain dar H, Liu B, tyurin vladimir, ritov, vladimir B, et al. (2017). Oxidized arachidonic and adrenic pes navigate cells to ferroptosis. nature CHeMICaL BIOLOGY | 13. 10.1038/nCHeMBIO.2238. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 36.Dierge E, Debock E, Guilbaud C, Corbet C, Mignolet E, Mignard L, Bastien E, Dessy C, Larondelle Y, and Feron O (2021). Peroxidation of n-3 and n-6 polyunsaturated fatty acids in the acidic tumor environment leads to ferroptosis-mediated anticancer effects. Cell Metab 33, 1701–1715.e5. 10.1016/J.CMET.2021.05.016. [DOI] [PubMed] [Google Scholar]
- 37.Yen CLE, Stone SJ, Koliwad S, Harris C, and Farese RV (2008). DGAT enzymes and triacylglycerol biosynthesis. J Lipid Res 49, 2283–2301. 10.1194/JLR.R800018-JLR200. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 38.Deleye Y, Cotte AK, Hannou SA, Hennuyer N, Bernard L, Derudas B, Caron S, Legry V, Vallez E, Dorchies E, et al. (2020). CDKN2A/p16INK4a suppresses hepatic fatty acid oxidation through the AMPKα2-SIRT1-PPARα signaling pathway. J Biol Chem 295, 17310–17322. 10.1074/JBC.RA120.012543. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 39.Zeng H, Jorapur A, Shain AH, Lang UE, Torres R, Zhang Y, McNeal AS, Botton T, Lin J, Donne M, et al. (2018). Bi-allelic Loss of CDKN2A Initiates Melanoma Invasion via BRN2 Activation. Cancer Cell 34, 56–68.e9. 10.1016/J.CCELL.2018.05.014. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 40.Rabhi N, Hannou SA, Gromada X, Salas E, Yao X, Oger F, Carney C, Lopez-Mejia IC, Durand E, Rabearivelo I, et al. (2018). Cdkn2a deficiency promotes adipose tissue browning. Mol Metab 8, 65–76. 10.1016/J.MOLMET.2017.11.012. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 41.Hannou SA, Wouters K, Paumelle R, and Staels B (2015). Functional genomics of the CDKN2A/B locus in cardiovascular and metabolic disease: what have we learned from GWASs? Trends Endocrinol Metab 26, 176–184. 10.1016/J.TEM.2015.01.008. [DOI] [PubMed] [Google Scholar]
- 42.Tsoi J, Robert L, Paraiso K, Galvan C, Sheu KM, Lay J, Wong DJL, Atefi M, Shirazi R, Wang X, et al. (2018). Multi-stage Differentiation Defines Melanoma Subtypes with Differential Vulnerability to Drug-Induced Iron-Dependent Oxidative Stress. Cancer Cell 33, 890–904.e5. 10.1016/J.CCELL.2018.03.017. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 43.Badgley MA, Kremer DM, Carlo Maurer H, DelGiorno KE, Lee HJ, Purohit V, Sagalovskiy IR, Ma A, Kapilian J, Firl CEM, et al. (2020). Cysteine depletion induces pancreatic tumor ferroptosis in mice. Science (1979) 368. 10.1126/SCIENCE.AAW9872/SUPPL_FILE/AAW9872S3.MOV. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 44.Mai WX, Gosa L, Daniels VW, Ta L, Tsang JE, Higgins B, Gilmore WB, Bayley NA, Harati MD, Lee JT, et al. (2017). Cytoplasmic p53 couples oncogene-driven glucose metabolism to apoptosis and is a therapeutic target in glioblastoma. Nature Medicine 2017 23:11 23, 1342–1351. 10.1038/nm.4418. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 45.Tannous BA (2009). Gaussia luciferase reporter assay for monitoring biological processes in culture and in vivo. Nature Protocols 2009 4:4 4, 582–591. 10.1038/nprot.2009.28. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 46.Hsieh WY, Williams KJ, Su B, and Bensinger SJ (2021). Profiling of mouse macrophage lipidome using direct infusion shotgun mass spectrometry. STAR Protoc 2, 100235. 10.1016/J.XPRO.2020.100235. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 47.BLIGH EG, and DYER WJ (1959). A rapid method of total lipid extraction and purification. Can J Biochem Physiol 37, 911–917. 10.1139/O59-099. [DOI] [PubMed] [Google Scholar]
- 48.Su B, Bettcher LF, Hsieh W-Y, Hornburg D, Pearson MJ, Blomberg N, Giera M, Snyder MP, Raftery D, Bensinger SJ, et al. (2021). A DMS Shotgun Lipidomics Workflow Application to Facilitate High-Throughput, Comprehensive Lipidomics. 10.1021/jasms.1c00203. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 49.Williams KJ, and Bensinger SJ (2020). Cellular Fatty Acid Analysis in Macrophage Using Stable Isotope Labeling. Methods Mol Biol 2184, 47–60. 10.1007/978-1-0716-0802-9_4. [DOI] [PubMed] [Google Scholar]
- 50.Martin M (2011). Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J 17, 10–12. 10.14806/EJ.17.1.200. [DOI] [Google Scholar]
- 51.Bushnell B (2014). BBMap: A Fast, Accurate, Splice-Aware Aligner. In (Conference: 9th Annual Genomics of Energy & Environment Meeting). [Google Scholar]
- 52.Li H, and Durbin R (2009). Fast and accurate short read alignment with Burrows–Wheeler transform. Bioinformatics 25, 1754–1760. 10.1093/BIOINFORMATICS/BTP324. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 53.Van der Auwera GA, and O’Connor BD (2020). Genomics in the Cloud: Using Docker, GATK, and WDL in Terra 1st ed. (O’Reilly Media; ). [Google Scholar]
- 54.Vivian J, Rao AA, Nothaft FA, Ketchum C, Armstrong J, Novak A, Pfeil J, Narkizian J, Deran AD, Musselman-Brown A, et al. (2017). Toil enables reproducible, open source, big biomedical data analyses. Nature Biotechnology 2017 35:4 35, 314–316. 10.1038/nbt.3772. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 55.Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, and Gingeras TR (2013). STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21. 10.1093/BIOINFORMATICS/BTS635. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 56.Li B, and Dewey CN (2011). RSEM: Accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics 12, 1–16. 10.1186/1471-2105-12-323/TABLES/6. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 57.Benjamin D, Sato T, Cibulskis K, Getz G, Stewart C, and Lichtenstein L Calling Somatic SNVs and Indels with Mutect2. 10.1101/861054. [DOI] [Google Scholar]
- 58.Fan Y, Xi L, Hughes DST, Zhang J, Zhang J, Futreal PA, Wheeler DA, and Wang W (2016). MuSE: accounting for tumor heterogeneity using a sample-specific error model improves sensitivity and specificity in mutation calling from sequencing data. Genome Biol 17, 178. 10.1186/S13059-016-1029-6/FIGURES/3. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 59.Koboldt DC, Zhang Q, Larson DE, Shen D, McLellan MD, Lin L, Miller CA, Mardis ER, Ding L, and Wilson RK (2012). VarScan 2: Somatic mutation and copy number alteration discovery in cancer by exome sequencing. Genome Res 22, 568–576. 10.1101/GR.129684.111. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 60.Tate JG, Bamford S, Jubb HC, Sondka Z, Beare DM, Bindal N, Boutselakis H, Cole CG, Creatore C, Dawson E, et al. (2018). COSMIC: the Catalogue Of Somatic Mutations In Cancer. Nucleic Acids Res 47, 941–947. 10.1093/nar/gky1015. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 61.Chakravarty D, Gao J, Phillips S, Kundra R, Zhang H, Wang J, Rudolph JE, Yaeger R, Soumerai T, Nissan MH, et al. (2017). OncoKB: A Precision Oncology Knowledge Base. 10.1200/PO.17.00011, 1–16. 10.1200/PO.17.00011. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 62.Talevich E, Shain AH, Botton T, and Bastian BC (2016). CNVkit: Genome-Wide Copy Number Detection and Visualization from Targeted DNA Sequencing. PLoS Comput Biol 12, e1004873. 10.1371/JOURNAL.PCBI.1004873. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 63.Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT, et al. (2000). Gene Ontology: tool for the unification of biology. Nature Genetics 2000 25:1 25, 25–29. 10.1038/75556. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 64.Kanehisa M, Furumichi M, Tanabe M, Sato Y, and Morishima K (2017). KEGG: new perspectives on genomes, pathways, diseases and drugs. Nucleic Acids Res 45, D353–D361. 10.1093/NAR/GKW1092. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 65.Jassal B, Matthews L, Viteri G, Gong C, Lorente P, Fabregat A, Sidiropoulos K, Cook J, Gillespie M, Haw R, et al. (2020). The reactome pathway knowledgebase. Nucleic Acids Res 48, D498–D503. 10.1093/NAR/GKZ1031. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 66.Mi H, Muruganujan A, Ebert D, Huang X, and Thomas PD (2019). PANTHER version 14: more genomes, a new PANTHER GO-slim and improvements in enrichment analysis tools. Nucleic Acids Res 47, D419–D426. 10.1093/NAR/GKY1038. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 67.Barbie DA, Tamayo P, Boehm JS, Kim SY, Moody SE, Dunn IF, Schinzel AC, Sandy P, Meylan E, Scholl C, et al. (2009). LETTERS Systematic RNA interference reveals that oncogenic KRAS-driven cancers require TBK1. 10.1038/nature08460. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 68.Love MI, Huber W, and Anders S (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 15, 1–21. 10.1186/S13059-014-0550-8/FIGURES/9. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 69.R Core Team (2021). R: A language and environment for statistical computing. https://www.r-project.org/. [Google Scholar]
- 70.Cerami E, Gao J, Dogrusoz U, Gross BE, Sumer SO, Aksoy BA, Jacobsen A, Byrne CJ, Heuer ML, Larsson E, et al. (2012). The cBio Cancer Genomics Portal: An Open Platform for Exploring Multidimensional Cancer Genomics Data. Cancer Discov 2, 401–404. 10.1158/2159-8290.CD-12-0095. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 71.Gao J, Aksoy BA, Dogrusoz U, Dresdner G, Gross B, Sumer SO, Sun Y, Jacobsen A, Sinha R, Larsson E, et al. (2013). Integrative Analysis of Complex Cancer Genomics and Clinical Profiles Using the cBioPortal. Sci Signal 6, 1–1. 10.1126/SCISIGNAL.2004088. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 72.Vivian J, Rao AA, Nothaft FA, Ketchum C, Armstrong J, Novak A, Pfeil J, Narkizian J, Deran AD, Musselman-Brown A, et al. (2017). Toil enables reproducible, open source, big biomedical data analyses. Nature Biotechnology 2017 35:4 35, 314–316. 10.1038/nbt.3772. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 73.Goldman MJ, Craft B, Hastie M, Repečka K, McDade F, Kamath A, Banerjee A, Luo Y, Rogers D, Brooks AN, et al. (2020). Visualizing and interpreting cancer genomics data via the Xena platform. Nature Biotechnology 2020 38:6 38, 675–678. 10.1038/s41587-020-0546-8. [DOI] [PMC free article] [PubMed] [Google Scholar]
- 74.Wang Q, Hu B, Hu X, Kim H, Squatrito M, Scarpace L, deCarvalho AC, Lyu S, Li P, Li Y, et al. (2017). Tumor Evolution of Glioma-Intrinsic Gene Expression Subtypes Associates with Immunological Changes in the Microenvironment. Cancer Cell 32, 42–56.e6. 10.1016/J.CCELL.2017.06.003. [DOI] [PMC free article] [PubMed] [Google Scholar]
Associated Data
This section collects any data citations, data availability statements, or supplementary materials included in this article.
Supplementary Materials
Supplemental Table Titles:
Table S1: GO biological processes results for gene clusters generated from correlation matrix of lipid species across tumors from patients with GBM. Related to Figure 1C.