Sexual dimorphism in glioma glycolysis underlies sex differences in survival

Throughout each stage of life, males and females differ on multiple levels as a consequence of differences in sex chromosome complement, as well as the organizational and activational effects of sex hormones (1). Net sex differences in phenotype arise through the combined effects of genetics, epigenetics, and metabolism from the time of fertilization. Multiple studies have demonstrated that male embryos grow significantly faster than female embryos (2, 3). Among the mechanisms that underlie these sex differences in growth rates are enhanced central carbon and amino acid metabolism in males and X-linked drivers of metabolism including the glycolytic enzyme phosphoglycerate kinase (PGK) (4). These metabolic effects may persist into adulthood, as recent evidence suggests that healthy males at rest and in exercise rely more on carbohydrate utilization than females who are more reliant on lipid metabolism (5). However, it is currently unclear whether these differences in metabolic mechanisms contribute to sex differences in human disease, specifically cancer.

Numerous epidemiological studies have identified the sex of the patient as a significant factor that impacts cancer incidence and survival. In general, males have a higher incidence and mortality relative to females in cancers throughout the body, with the exception of a few cancers (e.g., gallbladder and anorectal cancers) (6–8). Similar patterns are present in intracranial tumors. Although extra-axial meningiomas are more commonly diagnosed in females, males have an increased incidence of parenchymal brain tumors relative to females that is independent of age (9). For example, boys constitute the majority of group 3 and 4 medulloblastomas, with peak incidences less than 5 years of age and in the peripubertal period (10, 11). Moreover, boys with high-grade gliomas have significantly reduced overall and progression-free survival relative to girls (12). These findings suggest that existing sex disparities, at least in the case of brain cancers, are not completely attributable to the actions of circulating sex hormones. In fact, a recent multidimensional analysis of molecular differences in cancer between males and females determined strong sex effects in the presence of mutations, DNA methylation, transcript, and protein expression among renal, bladder, liver, lung, head and neck squamous, and thyroid cancers, but a weak sex effect in low-grade and high-grade gliomas (13). These findings suggest that there are alternative molecular mechanisms, perhaps related to metabolism, that may explain the sex disparity seen in brain cancers that are otherwise invisible to global molecular profiling studies.

Metabolism is a critical factor that is required for tumor survival and tumorigenicity. One of the hallmarks of cancer metabolism is enhanced glucose uptake and its conversion to lactate despite the presence of available oxygen (i.e., the Warburg effect or aerobic glycolysis) (14). This process not only results in the rapid generation of ATP for cellular energetics through glycolysis, but can also contribute to biosynthetic pathways required for proliferation (15, 16). Fluorodeoxyglucose (FDG) PET is a well-established tool for quantifying glucose uptake in tumors. Not only does FDG uptake positively correlate with glioma grade, but it inversely correlates with survival (17). Together, these findings suggest that sexual dimorphism in nutrient utilization may exist in brain cancers and that this might contribute to sex differences in survival, as well as require sex-specific interpretations of diagnostic tests like FDG-PET. Therefore, we sought to determine if the normal sex differences in glucose metabolism would have correlates in glioma metabolism and whether there would be opportunities for refined risk stratification by incorporating sex-specific analysis of glycolysis.

Glycolytic gene expression stratifies risk in gliomas. To determine if there was a sexual dimorphism in glioma glycolysis that could explain differences in survival, we investigated the lower-grade glioma (LGG) dataset in The Cancer Genome Atlas (TCGA) (18, 19). This dataset included transcriptomic and genomic data in a nearly equal number of male (n = 285) and female (n = 228) patients with grade 2 and 3 gliomas. No significant difference in overall survival (OS) between males and females existed in the LGG patients (Supplemental Figure 1; supplemental material available online with this article; https://doi.org/10.1172/jci.insight.92142DS1). Next, we assessed RNA sequencing (RNA-Seq) expression data of 36 transcripts encoding hexose transporters, glycolytic enzymes, and monocarboxylate (i.e., lactate and pyruvate) transporters (MCTs) in male versus female LGG samples. Overall, there were minimal but significant differences in only 2 of the 36 genes, with lactate dehydrogenase B (LDHB) exhibiting a significant but minimal expression increase in males compared with females (1.1-fold, P = 0.02) (Supplemental Figure 2). Hexokinase 1 (HK1), conversely, was slightly increased in females relative to males (1.1-fold, P = 0.02). Together, this was consistent with previous findings of a weak sex effect on transcript expression in the TCGA LGG dataset (13). However, we were interested in whether there might be glycolytic subtypes within each sex that correlated with survival, specifically. We hypothesized that subgroups within a sex with enhanced glycolytic gene expression would manifest decreased OS.

To determine whether there were glycolytic subgroups in males and females, we performed an unsupervised analysis using 36 glycolytic genes. We stratified and Z score–normalized the gene expression data by sex, and applied a K-means clustering analysis to separate the male and female LGG samples each into 2 clusters (i.e., high versus low glycolytic expression). Thirty-two male and 27 female samples were distinguished by their increased glycolytic gene expression. These were denoted as cluster 2, and the majority of male and female samples, which did not overexpress these transcripts, were denoted as cluster 1. Male and female cluster 2 was characterized by a total of 14 transcripts and 10 transcripts, respectively, with a mean expression Z-score value greater than 1 relative to male and female cluster 1 that had 0 transcripts (Supplemental Table 1). The dissimilarity of cluster 2 relative to cluster 1 in both males and females was further confirmed with multidimensional scaling analysis (MDS) (Figure 1). Despite the glycolytic overexpression seen in both male and female cluster 2, survival analyses of these clusters identified a sex difference in survival where cluster 2 males performed poorly compared with cluster 1 males and all females. Cluster 2 males had a median OS of 41.46 months compared with 98.16 months for cluster 1 males (P = 0.0005). No statistically significant glycolytic cluster–specific differences in OS were seen for females; cluster 2 had a median OS of 146.02 months compared with a cluster 1 median OS of 78.15 months (P = 0.3113) (Figure 1). Unbiased K-means clustering analyses using glycolytic gene expression led to 2 potentially significant discoveries: (a) a glycolytic gene expression threshold could exist above which males but not females are defined by decreased OS and (b) decreased male OS could be driven by a subset of these 36 glycolytic transcripts.

Figure 1

K-means clustering identifies sex differences in glycolysis. (A) Heatmap generated from the K-means (K = 2) clustering analysis identifies a cluster of males characterized by high glycolytic gene expression. (B) Multidimensional scaling (MDS) analysis demonstrates dissimilarity of the 2 clusters. (C) Survival analysis demonstrates that the cluster of males with glycolytic gene overexpression have significantly shorter survival than the remainder of males. (D–F) Same analyses performed for females, but no significant differences in overall survival were present. P values were calculated using the log-rank test. Numbers in parentheses refer to number of deaths/total patients in that group.

To optimally define glycolytic subgroups and determine which glycolytic transcripts contribute to survival differences, we developed a TCGA data mining algorithm that extracted survival information as a function of transcript level on a sex-specific basis using RNA-Seq data (Figure 2). First, we defined the optimal glycolytic gene expression threshold for stratifying survival differences in males. We applied an unbiased sliding Z-score threshold (range 0–3 in 0.25-unit increments; note that all genes have similar range after Z-score normalization regardless of sex) to glycolytic gene expression in both male and female LGG samples. Using the log-rank test to assess statistical significance in OS differences between the male subgroups, we determined that a Z score of 1.75 maximized male differences in survival (median OS difference = 75.99 months, hazard ratio [HR] 2.46, P = 0.0018). As expected, no Z-score threshold was able to identify female glycolytic subgroups showing a statistically significant OS difference (P = 0.9541) (Supplemental Table 2).

Figure 2

Survival-based identification of sex-specific differences in gene expression. RNA-Seq data from lower-grade gliomas were interrogated for males or females whose overexpression of a specific glycolytic transcript beyond a specific Z-score threshold conferred significantly different (P < 0.05) overall survival compared with other individuals within that sex. Significance calculated using the log-rank test.

Next, we used this optimized Z-score threshold to identify which of the 36 glycolytic transcripts were driving the survival differences in the male LGG samples. The Z-score threshold of 1.75 included 11 genes (GAPDH, LDHA, PGK1, HK3, PFKL, GCK, GPI, PGAM2, SLC2A5, SLC16A3, and SLC16A8) whose overexpression was associated with significantly decreased OS in males (Figure 3 and Supplemental Table 3). The male high-glycolytic group was defined as any male who overexpressed at least 1 of the 11 genes that was associated with significantly decreased survival, resulting in a total of 63 males. All other males were defined as male low-glycolytic. A total of 77 females overexpressed any 1 of the 11 genes and were assigned to the female high-glycolytic group. However, the only gene that was associated with decreased OS in females was hexokinase 3 (HK3). It was overexpressed in only 9 of the 77 high-glycolytic females (Supplemental Table 3). Male high-glycolytic patients had significantly decreased median OS relative to male low-glycolytic patients (36.33 months vs. 105.12 months, P < 0.0001). However, there were no differences seen in the female patient groups (75.1 months vs. 87.39 months, P = 0.5183). Similarly, male high-glycolytic patients had significantly decreased median disease-free survival (DFS) (22.17 months vs. 72.01 months, P < 0.0001) and females exhibited no differences (38.9 months vs. 45.14 months, P = 0.7396) (Figure 3).

Figure 3

Male gliomas are uniquely stratified by glycolytic gene overexpression. (A) Glycolytic pathway demonstrating 11 transcripts (red) whose overexpression confer decreased overall survival within a given sex. (B) Overall survival and (C) disease-free survival analyses reveal male-specific stratification based upon expression of 11 glycolytic transcripts. Any sample with overexpression of at least 1 glycolytic transcript was placed in the high-glycolytic group. All other samples in that sex were placed in the low-glycolytic group. P values were calculated using the log-rank test. Numbers in parentheses refer to number of deaths/total patients in that group.

We questioned whether outliers within the male high-glycolytic group might be skewing the survival differences. We investigated the heterogeneity of the high-glycolytic samples by looking at coexpression of the significant glycolytic transcripts within the male and female high-glycolytic groups. Interestingly, no single patient sample was associated with overexpression of all 11 genes. In fact, the majority of male samples (38 of 63 samples) were associated with overexpression of only 1 transcript, with a maximum of 8 coexpressed transcripts. Similar observations were seen with females, with the majority of samples overexpressing a single transcript (38 of 77 samples) (Supplemental Table 4).

To examine if there was a dose-dependent response of glycolytic transcripts on patient survival, we binned the males and females into 3 subtypes. Metabolic subtype 1 was characterized by overexpression of zero transcripts, subtype 2 by overexpression of 1–3 transcripts, and subtype 3 by overexpression of 4–8 transcripts. Thus, subtype 1 is analogous to the low-glycolytic group and subtypes 2 and 3 make up the high-glycolytic group. Interestingly, males were stratified into 3 groups with distinct survival profiles. Male metabolic subtype 3 had a median OS of 16.16 months, male subtype 2 had a median OS of 41.1 months, and male subtype 1 had a median OS of 105.1 months. These metabolic subtypes, however, failed to stratify female patients (Figure 4). Interestingly, comparison of the male glycolytic subtypes with individual glycolytic transcript expression data and the results of the initial K-means analysis demonstrated that not only was metabolic subtype 3 characterized by a majority of K-means cluster 2, but that subtype 3 was driven by LDHA, GAPDH, PGK1, SLC16A3, PFKL, and GPI expression. These findings were less evident for females (Figure 4). It is noteworthy that 3 individual approaches for detecting male-specific sex differences in glycolytic gene expression (unbiased K-means, optimal Z-score threshold, and dose dependence) did not incorporate any clinical, pathologic, or molecular data other than the sex of the patient and the corresponding gene expression data from that patient’s tumor. Thus, these findings could have immediately translatable applications for clinical practice that could be used to stratify patients in addition to conventional methods. This prompted us to investigate the effects of sex-specific glycolytic stratification on conventional glioma classification that uses histopathology and genomics.

Figure 4

Metabolic subtyping of glycolytic gene expression further stratifies male lower-grade glioma. (A) Male high-glycolytic and low-glycolytic groups were further stratified into 3 metabolic subtypes based upon the number of coexpressed glycolytic transcripts. Metabolic subtype 1 was characterized as zero coexpressed transcripts (i.e., the low-glycolytic group), subtype 2 was defined as coexpression of 1–3 transcripts, and subtype 3 was defined as coexpression of 4–8 transcripts. No sample coexpressed more than 8 of the 11 transcripts. (B) Metabolic subtyping does not stratify females. (C) Visualization of these metabolic subtypes demonstrate that metabolic subtype 3 has the shortest overall survival (OS) in males and is driven by a cohort of glycolytic transcripts including LDHA, GAPDH, PGK1, SLC16A3, PFKL, and GPI. Red boxes indicate patient samples whose overexpression of that transcript (Z > 1.75) resulted in significantly decreased survival for that sex. Blue boxes indicate samples whose overexpression of that transcript (Z > 1.75) did not result in significantly different survival for that sex. P values were calculated using the log-rank test. Numbers in parentheses refer to number of deaths/total patients in that group.

Glycolytic subtyping correlates with histopathologic classification of gliomas. Previous studies have demonstrated a positive correlation between the WHO histologic grade of gliomas with FDG uptake as measured by PET imaging (17, 20). Thus, we wondered if grade 3 gliomas would be enriched in the high-glycolytic group. While grade 3 gliomas were enriched in male high-glycolytic gliomas (61%), they were not exclusively composed of them. In fact, 50% of the male low-glycolytic gliomas were also grade 3 (Figure 5 and Supplemental Figure 3). Conversely, grade 2 gliomas were enriched in 50% low-glycolytic gliomas and in 39% male high-glycolytic gliomas (Figure 5 and Supplemental Figure 3). Although similar trends were seen in the female groups, no statistically significant enrichment was present in either males or females. Survival analyses demonstrated that glycolytic stratification was more robust in grade 3 males. Grade 3 male high-glycolytic gliomas had a median OS of 25.5 months versus low-glycolytic gliomas with a median OS of 74.97 months (P < 0.0001, Figure 5). Glycolysis also stratified male grade 2 gliomas but with less significance; high-glycolytic males had a median OS of 62.91 months compared with low-glycolytic males with a median OS of 144.94 months (P = 0.0011; Figure 5). Females were not significantly stratified. Together, these findings suggested that glycolytic stratification of males could refine tumor grading and support the use of FDG-PET in conjunction with histology for patient stratification.

Figure 5

Glycolytic subtyping correlates with histopathologic classification of gliomas. (A) Visualization of glycolytic groups, metabolic subtypes derived from those groups, and relationship to the histologic classification and WHO grade of the tumor. Both male and female high-glycolytic groups are enriched for astrocytoma histology, where only male low-glycolytic groups are enriched for oligoastrocytomas and oligodendrogliomas. Survival analysis of (B) grade 2 and (C) grade 3 gliomas, (D) astrocytomas, and (E) oligodendrogliomas reveal more robust glycolytic stratification for grade 3 versus grade 2 males and approximately equivalent survival for male glycolytic astrocytomas and oligodendrogliomas. P values were calculated using the log-rank test. Numbers in parentheses refer to number of deaths/total patients in that group.

Although the effect of tumor histology on glucose uptake is not as well defined compared with the glioma grade, we investigated the impact of glycolytic classification on tumor histology. We determined that astrocytomas were significantly enriched in 63% of male high-glycolytic gliomas, but only 31% of male low-glycolytic gliomas (P < 0.0001, Figure 5 and Supplemental Figure 3). Females had a similar distribution. Conversely, oligodendrogliomas and oligoastrocytomas were significantly enriched in the male low-glycolytic group, but not in the female group. Oligodendrogliomas showed a more robust enrichment in the male low-glycolytic group, with only 22% enrichment in the male high-glycolytic group versus 41% in the low-glycolytic group (P < 0.01, Figure 5 and Supplemental Figure 3). Survival analyses paralleled these findings, demonstrating that male high-glycolytic astrocytomas had the poorest median OS of 36.33 months compared with male low-glycolytic astrocytomas with a median OS of 98.16 months (P < 0.0001, Figure 5). Oligodendrogliomas were also characterized by robust glycolysis-based stratification, with a median OS of 26.74 months for the high-glycolytic males versus 117.31 months for the low-glycolytic males (P < 0.0001, Figure 5). Although patients with astrocytomas typically have shorter OS than patients with oligodendrogliomas (21), our glycolytic stratification scheme suggests that males with glycolytic astrocytomas perform equally poorly compared with males with glycolytic oligodendrogliomas.

Glycolytic subtyping correlates with genomic classification of gliomas. Several genetic alterations that are key drivers of LGGs have the ability to modulate glucose metabolism. We hypothesized that genomic alterations known to modulate glycolysis would be enriched in the high-glycolytic group and modulate male-specific survival. We focused on a group of key genomic alterations that have been characterized in LGG, specifically TP53, ATRX, IDH1, IDH2, PTEN, EGFR, NF1, CIC, and FUBP1 mutations as well as 1p/19q codeletion (21, 22). Mutations in TP53 and PTEN tumor suppressors and the epidermal growth factor receptor (EGFR) oncogene are known to activate glycolysis (23–27). Mutations of isocitrate dehydrogenase 1 and 2 (IDH1 and IDH2) and the resultant production of 2-hydroxyglutarate (2-HG) can potentially inhibit glucose metabolism (28–30). Mutations in the neurofibromin 1 (NF1) tumor suppressor may have the ability to regulate glucose metabolism, at least in part through enhanced Akt/mTOR activity and increased expression of glucose transporters (31–33). The effects of mutations in the Capicua transcriptional repressor (CIC) gene on glucose metabolism are not well characterized; however, CIC mutations can cooperatively regulate 2-HG levels with IDH1 mutations in cell lines (34). Although the effects of 1p/19q codeletion on glucose metabolism have not been characterized in detail, oligodendroglial tumors with this codeletion are characterized by enhanced FDG uptake (35).

To identify mutations that were enriched in the high-glycolytic or low-glycolytic groups, we plotted the LGG patient samples as a function of patient sex, glycolytic classification, and genomic alterations. We also incorporated the genomic subtype classification of the tumors as previously described (21, 22). Genomic subtype 1 tumors are classified by the presence of both IDH mutations and 1p/19q codeletion, subtype 2 tumors are classified by IDH mutations without 1p/19q codeletion, and subtype 3 tumors are classified as IDH wild type. The graphical analysis disclosed multiple important findings. First, our unbiased glycolytic subtyping classification correlated with the genomic subtype of the tumors. Genomic subtype 1, characterized by the highest OS, was significantly enriched in both male and female low-glycolytic groups. Only 3% of genomic subtype 1 tumors were classified as male high-glycolytic compared with 42% of these tumors in the male low-glycolytic category (P < 0.0001, Figure 6 and Supplemental Figure 4). Females demonstrated a similar pattern, characterized by 16% in the high-glycolytic group and 41% in the low-glycolytic group (P < 0.0001). In contrast, genomic subtype 3 tumors that are characteristically the poorest prognostic group were significantly enriched in the male high-glycolytic group. A total of 53% of genomic subtype 3 tumors were categorized as high-glycolytic compared with 9% that were low-glycolytic. (P < 0.0001 Figure 6 and Supplemental Figure 4). Similarly, 42% of female genomic subtype 3 tumors were classified as glycolytic compared with 7% that were low-glycolytic (P < 0.0001).

Figure 6

Glycolytic subtyping correlates with genomic classification of gliomas. (A) Visualization of glycolytic groups and metabolic subtypes reveal 3 classes of genomic alterations: those enriched in the low-glycolytic groups, those enriched in the high-glycolytic groups, and those that are not significantly different among groups. Survival analysis of (B) samples with both TP53 and ATRX mutations, (C) samples that are both TP53 and ATRX wild type, (D) samples with either an IDH1 or IDH2 mutation, and (E) samples with wild-type IDH1 and IDH2 reveal more robust glycolytic stratification for wild-type TP53/ATRX gliomas. Glycolytic classification unexpectedly stratifies wild-type IDH females, but not males. P values were calculated using the log-rank test. Numbers in parentheses refer to number of deaths/total patients in that group.

Next, we conducted an analysis of individual genomic alterations comprising the high-glycolytic and low-glycolytic groups. Based upon the enrichment of these alterations in the high-glycolytic and low-glycolytic groups, we defined 3 classes of genomic alterations. The first class consisted of alterations that were significantly depleted in the glycolytic group (FUBP1, CIC, and IDH mutations and 1p/19q codeletion). The second class consisted of alterations that were significantly enriched in the glycolytic group (PTEN, EGFR, and NF1 mutations) and the third class consisted of mutations that were not significantly enriched in either group (TP53 and ATRX mutations) (Figure 6 and Supplemental Figure 4). The majority of the genomic alterations showed significant differences in both males and females (i.e., IDH, CIC, EGFR, and PTEN mutations, and 1p/19q codeletion), with only 2 mutations that were significantly enriched only in males (i.e., NF1 and FUBP1). With the exception of TP53, our findings not only support previously observed glycolytic effects from these genomic alterations but ascribe potentially novel glycolytic effects to alterations including CIC and FUBP1 mutations.

TP53 and ATRX mutations are seen more commonly in astrocytomas compared with wild-type TP53 and ATRX oligodendrogliomas (21). Prompted by our previous findings that male astrocytomas and oligodendrogliomas could be stratified by glycolysis, we performed a survival analysis incorporating these mutations. Although our glycolytic classification scheme stratified males with both mutant and wild-type ATRX and TP53 tumors, wild-type TP53/ATRX tumors showed the most robust stratification, with high-glycolytic wild-type males performing even more poorly than those with TP53/ATRX mutations. Male high-glycolytic TP53/ATRX–mutant tumors had a median OS of 62.91 months (6 deaths in 46 total patients) compared with the low-glycolytic group with a median OS of 105.12 months (19 deaths in 79 total patients, P = 0.0360, Figure 6). Male high-glycolytic wild-type TP53/ATRX tumors performed even worse, with a median OS of 24.38 months (17 deaths in 34 total patients) compared with the median OS of the low-glycolytic group of 134.17 months (19 deaths in 101 total patients, P < 0.0001, Figure 6). As expected, females in both categories were not stratified by glycolysis. Similar findings were made when TP53 and ATRX were analyzed separately (Supplemental Figure 5). Although the role of TP53 as a prognostic biomarker in gliomas remains controversial and no consistent relationship has been found between the presence of TP53 mutations and prognosis (36, 37), loss of ATRX function is associated with better prognosis in gliomas (38, 39). This suggests that ATRX mutations could, in part, be driving the better survival seen in male patients with high-glycolytic mutant tumors compared with the male patients with the high-glycolytic wild-type tumors.

Because of the known inhibitory effects of IDH-mutant activity and its effects on prolonging patient survival (28–30), we investigated the effects of glycolytic stratification as a function of IDH mutation status in males and females. When the OS of IDH-mutant patients was analyzed, glycolysis significantly stratified males but not females. Male high-glycolytic IDH-mutant patients had a median OS of 44.38 months compared with 117.31 months (P < 0.0045, Figure 6). Females did not stratify. However, an unexpected discovery was made when wild-type IDH patients were stratified with glycolysis. Although the OS of these patients was expectedly shorter than that seen in the IDH-mutant patient population, males did not stratify. Females not only stratified, but the stratification was opposite to the findings seen in males. Female wild-type IDH high-glycolytic patients performed significantly better than all other patient groups, with a median OS of 51.84 months compared with female low-glycolytic patients (median OS 16.82 months, P = 0.0016). Both male groups had a similar median OS of 21.29 and 19.45 months for high-glycolytic and low-glycolytic groups, respectively (Figure 6). Together, these findings suggest a previously uncharacterized phenomenon that could involve the metabolism-modulating activity of IDH mutations and the potential for sexually dimorphic preferences in tumor nutrient consumption, analogous to what has been documented with carbohydrate utilization in males and lipid utilization in females (5).

Metabolomic identification of sex differences in glycolysis. To validate our transcriptome-based survival data, we analyzed quantitative metabolome-level metadata obtained previously (16) as well as a new dataset from patients with grade 2 gliomas and correlated metabolite levels to survival (see Methods section). We analyzed a total of 26 male and 19 female grade 2 gliomas. We hypothesized that glycolytic metabolite levels would selectively stratify males, similar to the effects we observed with glycolytic transcripts. However, we needed to develop an alternative approach to identify sex differences using metabolite data. Metabolite levels represent a snapshot of metabolism at a single point in time and are a complex function of metabolite synthesis and consumption. Therefore, enhanced glycolytic flux can result in both enriched and depleted metabolites within the pathway as has been seen in glioblastoma (16). This suggested to us that our K-means and scalable Z-score methods would not be an optimal approach here, as these methods assume that poor prognostic patients would be characterized by either enrichment or depletion of all metabolites. Instead, we used a previously published algorithm to measure metabolite-specific thresholds that could selectively stratify males (40). Within the metabolomics datasets, we identified 7 glycolytic pathway metabolites (glucose, 3-phosphoglycerate, fructose-6-phosphate, dihydroxyacetone phosphate, phosphoenolpyruvate, pyruvate, and lactate). Levels of these metabolites were statistically similar between males and females, paralleling the similarities in overall transcript abundance seen in LGG (Supplemental Figure 6).

Out of the 7 metabolites, 3 selectively stratified male, but not female grade 2 glioma patients. First, we determined that males whose gliomas had lower levels of pyruvate did worse than males whose gliomas had higher pyruvate levels (Figure 7). No such effect was seen in female glioma patients. This was especially interesting, as pyruvate is a substrate for 3 proteins whose transcripts were identified as prognostic in our analysis. Lactate dehydrogenase A (LDHA) that converts pyruvate to lactate (and an important component of poor prognostic metabolic subtype 3) and the MCTs that transport lactate and/or pyruvate (SLC16A3 and SLC16A8) across the plasma membrane were part of this signature. LDHA overexpression in solid tumors including gliomas has been implicated in tumor progression, immune evasion, and reduced survival in solid tumors including gliomas (41–44). Intriguingly, there is evidence of a sexual dimorphism in LDHA expression in both serum and tissue samples in healthy animals (45, 46). SLC16A3 (MCT4) is widely expressed in glycolytic tissues including tumors. Less is known about SLC16A8 (MCT3) that has restricted expression in the eye (47–52). Together, these findings suggested that pyruvate metabolism and transport could represent a clinically measurable and actionable component of sex differences in glioma glycolysis.

Figure 7

Glycolytic metabolites modulate sex-specific survival in grade 2 gliomas. (A) Pyruvate (Pyr) and (C) lactate/pyruvate (Lac/Pyr) levels stratify survival in males, but not in females (panels B and D). Significance calculated with log-rank tests.

These findings prompted us to also investigate an additional metabolic parameter, the lactate/pyruvate (lac/pyr) ratio. The lac/pyr ratio is a normalized quantity that is a measure of the equilibrium between the lactate product and pyruvate substrate of LDH, a surrogate marker for the cellular oxidation-reduction state, and a relative measure of glycolytic activity relative to mitochondrial respiration that is enhanced in hypoxia (53, 54). Interestingly, we determined that males with elevated lac/pyr did more poorly (median OS = 50 months) than did males with low lac/pyr (undefined median OS, P = 0.0497). Female gliomas trended in the opposite direction, as tumors with higher lac/pyr levels did better, although this effect was not significant (P = 0.2367, Figure 7). This supported our hypothesis that pyruvate metabolism mediated through enhanced LDH activity could potentially drive sex differences in LGG. Other significant metabolites included fructose-6-phosphate, the metabolic product of GPI (a component of the 11-gene signature) and dihydroxyacetone phosphate that is upstream of GAPDH (Supplemental Figure 7). These findings validated our transcriptomic results, further suggesting that survival risk in males but not females might be stratified with measures of glycolysis such as FDG-PET.

Warburg observed that cancer cells, despite having adequate oxygen, preferentially generate energy through enhanced glucose uptake and conversion to lactate instead of further oxidation of pyruvate in the mitochondria (Warburg effect or aerobic glycolysis) (14). The relevance of aerobic glycolysis to glioma biology and patient outcome is evidenced by the prognostic significance of FDG-PET, where the uptake of FDG positively correlates with glioma grade and inversely correlates with survival (17, 55). Sex differences in glucose metabolism are evident immediately after fertilization and throughout life under normal circumstances. Here, we demonstrate that sex differences in glucose metabolism are also evident in glioma and have sex-specific effects on outcome. Within the male population, there is a subset of males with significantly altered levels of specific transcripts and metabolites within glycolysis who exhibit significantly decreased OS. Thus, while levels of glycolytic transcript expression are similar in male and female glioma patients overall, sex-specific glycolytic phenotypes are more completely described by both gene expression and metabolite profile. This enhanced glycolytic phenotyping is not readily identifiable through transcriptomic analyses alone and is likely to be informative regarding sex differences in glioma biology and patient outcome.

Our findings also propose that clinical imaging can be used to noninvasively identify prognostic sex differences in metabolism. Although our transcriptome and metabolome data suggest that FDG uptake in male and female gliomas will be similar, clinical FDG-PET imaging in glioma care may require different reference ranges and interpretations for males versus females. In addition to FDG-PET imaging, there has been significant interest in emerging metabolic imaging technologies that measure hyperpolarized (HP) [¹³C]pyruvate uptake and metabolism in tumors with magnetic resonance spectroscopic imaging (MRSI) (56, 57). The possibility of a prognostic sexual dimorphism in pyruvate/lactate transport and LDH activity in gliomas also suggests that in vivo quantification of lactate/pyruvate flux using HP imaging could also be used as a complement to FDG-PET to identify prognostic sex differences in glioma metabolism.

How sex exerts its effects on metabolism and glioma biology remains to be defined. Sex differences in metabolism in the immediate postfertilization period are driven by the differences in sex chromosome complement, sex-specific reprogramming of imprinted loci, as well as rapidly evolving sexual dimorphism in global epigenetics at the level of histone modification and DNA methylation. These processes highlight an interdependency between epigenetics and metabolism in which epigenetics determines sex-specific expression of metabolic enzymes and metabolites function as cofactors for epigenetic modifications (58).

The presence of sex differences in epigenetic modifications and metabolic pathways suggest that there may also be sex differences in the activity of signaling and metabolic pathways driven by oncogenes and tumor suppressors. Mutations in the PTEN tumor suppressor, for example, not only enhance proliferation and glucose metabolism through the upregulation of downstream Akt and mTOR signaling (24), but there is also ample evidence of sexual dimorphism in the expression and activity of the PI3K/Akt/mTOR pathway in brain liver, heart, skeletal muscle, and adipose tissue in both normal and noncancerous pathologic states (59–61). Our data showing sex-specific survival differences associated with glycolytic pathway gene expression highlight a previously uncharacterized finding that sex differences in metabolism may have a role in sex-specific survival. It further suggests that signaling cascades that regulate these metabolic pathways may be drivers of these sex differences that will need to be elucidated in future investigations.

Our findings add to the growing body of evidence that there are sex differences in the biology of gliomas. Although little is documented with LGG, some intriguing discoveries have been made with glioblastoma multiforme (GBM). In a mouse model of mesenchymal GBM, combined loss of Nf1 and Tp53 function in male astrocytes resulted in enhanced proliferation, greater induction of a stem cell–like population, and tumorigenesis relative to female transformed astrocytes. Possible mechanisms underlying these sex-specific effects include cell-intrinsic differences in RB1 activation, which were greater in females (62). A second rodent model of GBM in which male tumors exhibit higher proliferation, angiogenesis, and metabolic activity relative to females supports these findings (63). Clinically, there are also recognized phenotypic sex differences in GBM. While necrosis is significantly higher in male GBM, females are selectively stratified into prognostically significant higher- and lower-necrosis groups (64). Although it is currently unclear as to how the female phenotypic variability in GBM relates to the male molecular variability seen in LGG, many factors may be involved, such as the presence of specific mutational drivers (e.g., TP53 or MYC) that have been suggested to contribute to sex differences in GBM (64).

Alterations in TP53 function are at the core of cancer biology. The impact of TP53 loss of function has previously been described as sex dependent. Deletion of Tp53 in mice results in disproportionate loss of female embryos from neural tube defects, subsequently ascribed to differences in X chromosome dosage and Lyar function (65, 66). The presence of sexual dimorphism in response to loss of Tp53 activity has also been reported in neurofibromatosis 1–null (Nf1-null) mouse astrocytes. Here, combined loss of Tp53 and Nf1 function resulted in significantly enhanced growth rates, clonogenic potential, and in vivo tumorigenesis in male compared with female astrocytes (62). Despite the lack of differences in TP53 mutation enrichment in either the high-glycolytic or low-glycolytic groups, we found that male TP53/ATRX wild-type patients actually did worse than the patients with mutant tumors when glycolytic gene overexpression was considered. Because the TP53 tumor suppressor has numerous effects on cellular metabolism that can be modulated by mutations (23), the etiology for these survival effects are currently unclear. Although these survival differences could be attributable to the potentially beneficial effects of ATRX loss as described above, one additional possibility may involve interactions between TP53 and currently uncharacterized drivers of glycolysis. However, this will need to be interrogated in future studies.

Mutations of IDH were also observed to interact with the glycolytic phenotype. IDH mutations are found in gliomas as well as acute myelogenous leukemia (AML) (67). In general, IDH mutations are associated with enhanced prognosis in glioma patients (28). The mechanism by which IDH mutations lead to improved prognoses is unclear but involvement of the metabolite 2-HG and its potential inhibition of glucose metabolism has been proposed (29, 30). Our data support an interaction between IDH mutations and levels of 2-HG as determinants of survival. However, our data also show a previously uncharacterized, discordant effect of IDH status on survival in males and females when the level of glycolytic transcription is considered. Although increased glycolytic gene expression stratifies males with IDH mutations, as seen with the rest of the comparisons in this study, wild-type IDH patients stratify in the opposite direction where glycolysis identifies poor prognostic females but not males. This is a previously uncharacterized finding that must be investigated further and may even reflect the biology behind a potential sexual dimorphism in GBM pathogenesis.

Although there are inherent limitations to this data-mining analysis, as it is based upon previously acquired genomic, transcriptomic, and metabolomic data, many important questions arise that stress the importance of the sex of the patient and the metabolism of the patient’s tumor in both tumor classification and patient stratification. Large multicenter prospective trials are needed to further validate and establish the relevance of these findings. We propose that imaging studies of glioma glycolysis with FDG-PET should be reevaluated by the oncologic community to incorporate additional factors such as the sex of the patient, genomic alterations, gene expression, and biochemical/metabolic markers. Integration of these currently clinically available technologies through a new sex-specific lens may pave the way to new advancements in precision medicine.

Datasets. Level 3 RNA-Seq gene expression for TCGA LGG samples were obtained from the NCI Genomic Data Commons data portal and Broad GDAC Firehose data portal. The mutation information for the LGG samples was obtained from the GDAC firehose Oncotated Calls MAF files. Clinicopathologic data for these samples were downloaded from the cBioPortal for cancer genomics (http://www.cbioportal.org/). Neoplasm histologic type and neoplasm histologic grade were used to define the histology and grade of the LGG samples. Only tumor samples that represented primary tumors were used and all recurrent tumor samples were excluded from the analysis. In total, molecular data were available for 228 females and 285 males and OS data available for 227 females and 283 males.

Inferring 1p/19q codeletions of LGG samples. Because 1p/19q deletions for samples are not annotated in TCGA, we inferred the codeletion status of the LGG samples using SNP-based loss-of-heterozygosity (LOH) analysis based on the copy number variation data (CNV) obtained from the Broad GDAC Firehose database (68). In brief, the focal somatic CNV in LGG samples were determined using GISTIC 2.0 (69). The segment mean is the log₂ ratio of the tumor intensity to the normal intensity. Conversion to an absolute CN value can be done by applying 2^{segment mean} × 2. Similar to a previously published approach using the TCGA (21), we further inferred regions of LOH with an absolute CN value less than 1.8, and aggregated different focal CNV into the corresponding chromosome arm positions, and determined the 1p1/19q codeletion by assessing whether the 1p and 19q are over 80% deleted. One hundred sixty-eight out of 513 samples were determined to be 1p/19q codeleted, with an average of 94.5% of 1p and 86% of 19q being deleted with a standard deviation less than 1%. Comparison with the published TCGA analysis (with 293 samples analyzed in that publication; see ref. 21) showed that we were able to identify 86 additional samples with 1p/19q codeletion apart from their 83 samples. Two of the samples previously reported as codeleted in TCGA (TCGA-CS-5394-01 and TCGA-DU-5870-01) were excluded from this group as they are mostly 1p deleted, but only 60% deleted in 19q.

Glycolytic pathway gene expression analyses. Gene expression values from 36 genes that characterize hexose uptake (SLC2A1, SLC2A2, SLC2A3, SLC2A4, and SLC2A5), glycolysis (HK1, HK2, HK3, GCK, GPI, PFKM, PFKL, PFKP, ALDOA, ALDOB, ALDOC, GAPDH, GAPDHS, PGK1, PGK2, PGAM1, PGAM2, ENO1, ENO2, ENO3, PKM2, PKLR, LDHAL6B, LDHA, LDHB, and LDHC), and monocarboxylate (lactate and pyruvate) transport (SLC16A1, SLC16A3, SLC16A7, and SLC16A8) were analyzed.

To identify sex differences in expression level among glycolytic genes, RSEM expression for each gene was plotted individually for males and females. A 2-tailed t test was used to determine significant expression differences between males and females. For all subsequent analyses, the gene expression value was transformed into a Z score that was specific to the sex of the patient.

K-means clustering of LGG samples. The sex-specific Z scores of the 36 candidate glycolytic genes were input into a K-means clustering analysis using the Hartigan-Wong algorithm. To overcome potentially different cluster assignment using different starting points, 1,000 random starting points were used to assign female and male patients each into 2 clusters, while the cluster membership of male and female patients was ultimately defined by the random starting point leading to the optimal separation of the 2 clusters. MDS plots were used to illustrate that the smaller of the 2 clusters corresponds to the extreme groups in the first MDS coordinate. A heatmap of the glycolytic Z scores of each sex was generated with samples ordered by cluster membership and genes were clustered based on Euclidean distance and average linkage. The 2 clusters of each sex were associated with survival endpoints by the Kaplan-Meier method and the log-rank test was used to compare survival difference between clusters. The HRs with 95% CIs of patient groups were estimated. All analyses were conducted in R (version 3.3.1).

Survival-based algorithm for sex differences in glycolysis. To determine a Z score for glycolytic gene expression that could maximally stratify males, the Z-score threshold using all 36 genes simultaneously was scaled from 0 to 3 in 0.25-unit increments. The log-rank test was used to calculate the P value indicating statistical significance of survival difference and derive the HR with 95% CI. The Z-score threshold that maximized differences in male survival (Z = 1.75) was selected for transcript-specific analyses.

To determine specific glycolytic transcripts whose overexpression (Z > 1.75) resulted in decreased OS, each of the 36 transcripts was analyzed individually among males and females separately. Any patient sample with a sex-specific Z score for that specific transcript that was over the threshold was labeled as selected. Conversely, any sample whose Z score was below the threshold was labeled as unselected.

The OS, DFS, and survival status of the selected samples were extracted and compared with those from the unselected samples by fitting the data by the Kaplan-Meier method to determine the median survival times of both male groups. The log-rank test was also used to calculate P values to determine significance of the survival differences. Genes whose overexpression resulted in a significantly (P ≤ 0.05) decreased median OS in the selected samples versus the unselected samples were saved. The same procedures were applied for all female samples independently.

Genomic alteration analyses. Genomic alteration data were compiled for all LGG samples. All mutations for a given gene were used. IDH-mutant samples were defined as having either IDH1 or IDH2 mutations. To determine significant differences in mutation enrichment between groups, a 2-tailed Fisher exact test was performed using the total number of samples with a given mutation versus wild type and high-glycolytic versus low-glycolytic classification. Data were plotted as the fraction of samples with a mutation within either the high-glycolytic or low-glycolytic group.

Metabolomic analysis of patient survival. All glioma samples were obtained retrospectively from the H. Lee Moffitt Cancer Center (Tampa, Florida, USA). Quantitative metabolomic datasets of grade 2 gliomas (9 females and 8 males) performed at Metabolon Inc. were published previously (16). An additional 10 female and 18 male grade 2 glioma datasets that were also generated by Metabolon Inc. using the same methods published previously (16). For statistical analyses, missing values (usually due to detection limit) were imputed with the compound minimum value (16, 70–72). Glycolysis metabolites (glucose, 3-phosphoglycerate, fructose-6-phosphate, dihydroxyacetone phosphate, phosphoenolpyruvate, pyruvate, and lactate) were isolated from the datasets. The lac/pyr ratio was calculated by dividing the lactate quantity by the pyruvate quantity. Metabolite quantities were converted to sex-specific Z scores and datasets were merged. A 2-tailed Mann-Whitney-Wilcoxon test was performed to establish significance in metabolite levels between males and females. For survival analyses, a biomarker cutoff optimization algorithm was used to determine metabolite levels that can maximally stratify male gliomas (40). The OS and survival status of the patient samples were fitted to a Kaplan-Meier model. The log-rank test was used to determine the P value and significance of the differences. Optimal Z-score thresholds used for the significant metabolites were as follows: fructose-6-phosphate (Z = –0.6113), dihydroxyacetone phosphate (Z = 0.1054), pyruvate (Z = –0.01736), and lac/pyr ratio (Z = –0.2818).

Statistics. Statistics involved with transcriptome and metabolome profiling are detailed in the sections above. All Student’s t tests and Fisher exact tests were performed 2-tailed. All survival analyses were performed using the Kaplan-Meier method and log-rank test using GraphPad Prism. A P value less than 0.05 was considered significant for all tests.

Study approval. All transcriptome data were obtained from TCGA. For metabolome data, Institutional Review Board/Human Subjects approval from H. Lee Moffitt Cancer Center was obtained prior to the study.