A Predictive Model for Selective Targeting of the Warburg Effect through GAPDH Inhibition with a Natural Product

Abstract

Targeted cancer therapies that use genetics are successful, but principles for selectively targeting tumor metabolism that is also dependent on the environment remain unknown. We now show that differences in rate-controlling enzymes during the Warburg effect (WE), the most prominent hallmark of cancer cell metabolism, can be used to predict a response to targeting glucose metabolism. We establish a natural product, koningic acid (KA), to be a selective inhibitor of GAPDH, an enzyme we characterize to have differential control properties over metabolism during the WE. With machine learning and integrated pharmacogenomics and metabolomics, we demonstrate that KA efficacy is not determined by the status of individual genes, but by the quantitative extent of the WE, leading to a therapeutic window in vivo. Thus, the basis of targeting the WE can be encoded by molecular principles that extend beyond the status of individual genes.

Keywords: Warburg effect; cancer metabolism; glucose metabolism; metabolic control analysis; metabolic flux analysis; metabolomics; natural product; pharmacogenomics; precision medicine; systems biology.

Figure 1. Thermodynamic and kinetic analysis of rate-control in glycolysis

(A) Glycolysis pathway with metabolites (black) and enzymes (blue). Model for flux control coefficients (FCCs) and reaction free energy (ΔG) values. f= flux. (B) Median FCCs vs. median ΔGs of glycolytic enzymes. HK, hexokinase; PFK, phosphofructokinase; GAPDH, glyceraldehyde-3-phosphate dehydrogenase. (C) Schematic showing flux inputs and outputs during Non-Warburg Effect and Warburg Effect conditions. (D) Relationship between FCC values of GAPDH, HK, and PFK and extent of the Warburg Effect defined as ratio of lactate production flux to oxidative phosphorylation flux. (E) Schematic demonstrating the high rate-limiting effect of GAPDH over cells carrying out the Warburg Effect compared to those undergoing oxidative phosphorylation. GA3P, glyceraldehyde-3-phosphate; 1,3-BPG, 1,3-bisphosphoglycerate.

Figure 2. Comparative metabolomics nominates a specific GAPDH inhibitor

(A–E) Clustered heatmaps with pathway and dose annotations of dose-dependent global metabolic responses to putative GAPDH inhibitors. (A) Arsenate (As₂O₄). (B) Arsenic trioxide (As₂O₃). (C) 3-bromopyruvate (3BP). (D) Iodoacetate (IA). (E) Koningic acid (KA). (F) Network-based pathway analysis for each compound. (G) Pathway analysis showing top 4 highest scoring pathways in response to KA. (H) KA docking analysis to GAPDH active site Cys152.

Figure 3. Expression of a fungal-derived KA-resistant GAPDH allele renders human cells completely resistant to KA and reverses their metabolic profile

(A) Schematic showing expression of a resistant allele of GAPDH by Trichoderma virens. (B) Immunoblotting of parental, EV or KAr-GAPDH expressing HEK293T cells. (C) Cell viability of HEK293T cells expressing KAr-GAPDH or EV (top left). Representative images of well (bottom left, top right) KAr-GAPDH or EV expressing cells treated with vehicle (0 μM) or KA (200μM). (D) Volcano plots showing metabolite profiles of HEK293T cells expressing EV compared to those expressing KAr-GAPDH after treatment with DMSO or 5μM KA. Log₂ fold change versus −log₁₀ (p-value). Dotted lines along x-axis represent ±log₂(2) fold change and dotted line along y-axis represents −log₁₀(0.05). Glycolysis metabolites shown as red points. All other metabolites are black points. (E) Glycolytic metabolite levels (F) Pentose Phosphate Pathway levels. (G) TCA cycle metabolite levels. G6P, Glucose-6-Phoshate; F 1,6-BP, Fructose 1,6-Bisphosphate; DHAP, Dihydroxyacetone Phosphate; 3PG, 3-Phosphoglycerate; PEP, Phosphoenolpyruvate; S7P, Sedoheptulose-7-Phosphate; E4P, Erythrose-4-Phosphate; R5P, Ribose-5-Phosphate. All data are represented as mean ±SEM from n = 3 biological replicates unless otherwise noted.

Figure 4. The cytotoxic response to KA treatment is heterogeneous

(A) Waterfall plot showing the difference in response of KA to 60 independent cell lines treated with vehicle (0.01% DMSO) or 10μM KA. Representative KA-resistant cell lines (Red, *) and KA-sensitive cell lines (Green, *). (B) Relative GAPDH activity in representative KA-sensitive and resistant cell lines in response to vehicle (DMSO) or KA. SK-MEL-28 and UACC-257 were treated with vehicle or 1μM KA; NCI-H522 and NCI-H226 were treated with vehicle or 0.4μM KA; BT-549 and MCF-7 were treated with 0.7μM KA (n=2). (C) Pearson correlation of KA IC₅₀ values for KA-sensitive and resistant cell lines versus percent of GAPDH activity. (D) Volcano plots showing metabolite profiles of breast cancer cell lines after treatment with DMSO or 90μM KA. Log₂ fold change versus −log₁₀ (p-value). Dotted lines along x-axis represent ±log₂(2) fold change and dotted line along y-axis represents −log₁₀(0.05). Glycolysis metabolites shown as red points. All other metabolites are black points. (E) Melanoma cell lines as in (D). (F) Non-small cell lung cancer cell lines as in (D). (G) Kinetic flux profiling for lactate labeling from ^13-C-glucose. (H) Change in lactate flux in response to KA based on fold changes vs. relative fluxes of the vehicle group from ¹³C-lactate enrichment from U-¹³C-glucose. (I) Glycolysis profiles for KA-sensitive and resistant breast, melanoma, and NSCLC cell lines treated with vehicle or their respective KA IC₅₀ concentrations. S7P, sedoheptulose-7-phosphate; F 1,6-BP, fructose 1,6-bisphosphate; G6P, glyceraldehyde-6-phosphate; DHAP, dihydroxyacetone phosphate; X5P, xylulose-5-phosphate; PEP, phosphoenolpyruvate; UDP-D-glucose, uridine diphosphate-D-glucose. All data are represented as mean ±SEM from n = 3 biological replicates unless otherwise noted.

Figure 5. A multi-omics analysis reveals that only the extent of the WE predicts KA response

(A) Schematic workflow of spearman rank correlation calculations from metabolite consumption and excretion rates, cell size, doubling time, mutation and copy number alteration (CNA) frequencies, and gene and protein expression patterns with KA IC₅₀ response across NCI-60 cell line collection. (B) Cell size correlation with KA IC₅₀ response across NCI-60 cell line panel. (C) Doubling time correlation as in (B). (D) Quantile-Quantile (Q-Q) plot comparing the quantiles from the correlation between mRNA expression levels and KA IC₅₀ values in NCI-60 cell line screen with random values. (E) Glycolysis gene expression correlation as in (B). (F) Q-Q plot for protein expression levels as in (D). (G) Glycolysis protein expression correlation as in (B). (H) Correlations of uptake and secretion rates of metabolites. (I) Glucose uptake correlation as in (B). (J) Lactate secretion correlation as in (B). (K) Q-Q plot for metabolic flux as in (D). (L) Receiver-operating characteristic (ROC) curves of standard cancer therapies with specific targets and/or biomarkers. Area under the curve (AUC) assesses biomarker accuracy.

Figure 6. KA is bioavailable and induces dynamic changes in glycolysis in vivo

(A) Schematic of timeline and treatment regimen of E2 treated female nu/nu mice injected with BT-474 cells orthotopically. These mice were treated with either saline or 1mg/kg, 2.5mg/kg, 5mg/kg, or 10mg/kg to identify the maximum tolerated dose (MTD) for 24 hours. (B) Plasma from MTD analysis for doses of 0–10mg/kg KA. (C) Levels of upper glycolytic intermediates from 0–24 hours (n=2 per group). (D) Levels of lower glycolytic intermediates as in (C). (E) Volcano plot showing metabolic profile at 16 hour time-point. Log₂ fold change versus −log₁₀ (p-value). Dotted lines along x-axis represent ±log₂(2) fold change and dotted line along y-axis represents −log₁₀(0.05). Glycolysis and related metabolites shown as red points. All other metabolites are black points (n=2). (F) 10 minute time-point as in (E). (G) Schematic showing the dynamic response KA treatment. G6P, glucose-6-phosphate; F 1,6-BP, fructose 1,6-bisphosphate; DHAP, dihydroxyacetone phosphate; GA3P, glyceraldehyde-3-phosphate; 1,3-BPG, 1,3-bisphosphoglycerate; 3PG, 3-phosphoglycerate; PEP, phosphoenolpyruvate. All data are represented as mean ±SEM from biological replicates.

Figure 7. A therapeutic window for precise targeting of the WE with KA can be achieved in vivo

(A) Schematic of timeline and treatment of E2 treated female nu/nu mice injected with KA-sensitive BT-474 cells and female nu/nu mice injected with KA-resistant MDA-MB-231 cells orthotopically and treated with either saline or 1mg/kg KA. (B) Average tumor volume (mm³) over 14 days in mice injected with BT-474 or MDA-MB-231 cells treated with vehicle or 1mg/kg KA (n=10 per group). (C) Tumor volume of each individual mouse from day 0 and day 14 in the BT-474 and MDA-MB-231 tumor models (n=10). (D) Representative BT-474 and MDA-MB-231 tumors from vehicle treated and 1mg/kg KA treated E2 treated female nu/nu mice on day 14. (E) Representative hematoxylin and eosin (H&E) and immunohistochemical (IHC) staining of serial sections from nu/nu mice injected with BT-474 cells treated with vehicle or 1mg/kg KA. (F) Volcano plot showing metabolic profile from skeletal muscle treated with vehicle or 1mg/kg KA from BT-474 tumor model. Log₂ fold change versus −log₁₀ (p-value). Dotted lines along x-axis represent ±log₂(2) fold change and dotted line along y-axis represents −log₁₀(0.05). Significantly different metabolites shown as red points. All other metabolites are black points (n=3). (G) Complete blood count (CBC) for BT-474 mice treated with vehicle or 1mg/kg KA after 14 days (n=9 for vehicle group; n=10 for treatment group). (H) Relative levels of creatinine from plasma (n=10 per group). (I) Relative levels of bilirubin from plasma (n=10 per group). (J) Pyruvate, alanine, and pyruvate to alanine ratio in liver (n=10 per group). Data are represented as mean ±SEM, biological replicates. *p<0.05, **p<0.01, ***p<0.001, ****p<0.0001 as determined by two-way ANOVA.