Probabilistic model checking of cancer metabolism

Abstract

Cancer cell metabolism is often deregulated as a result of adaption to meeting energy and biosynthesis demands of rapid growth or direct mutation of key metabolic enzymes. Better understanding of such deregulation can provide new insights on targetable vulnerabilities, but is complicated by the difficulty in probing cell metabolism at different levels of resolution and under different experimental conditions. We construct computational models of glucose and glutamine metabolism with focus on the effect of IDH1/2-mutations in cancer using a combination of experimental metabolic flux data and patient-derived gene expression data. Our models demonstrate the potential of computational exploration to reveal biologic behavior: they show that an exogenously-mutated IDH1 experimental model utilizes glutamine as an alternative carbon source for lactate production under hypoxia, but does not fully-recapitulate the patient phenotype under normoxia. We also demonstrate the utility of using gene expression data as a proxy for relative differences in metabolic activity. We use the approach of probabilistic model checking and the freely-available Probabilistic Symbolic Model Checker to construct and reason about model behavior.

Figure 1

Modeled pathway and model creation workflow. (A) Glycolysis and TCA cycle reactions modeled, including neomorphic activity of IDH Mutant. For glycolysis the following metabolites and abbreviations were used: Gluc, glucose; G6P, glucose-6-phosphate; Pyr, pyruvate; Lac, lactate; ACoA, acetyl-CoA; PPP, pentose phosphate pathway; Lip, lipids. For TCA cycle the following metabolites and abbreviations were used: Gln, glutamine; Cit, citrate; aKG, alpha-ketoglutarate; Mal, malate; Oac, oxaloacetate; 2HG, 2-hydroxyglutarate. (B) Workflow of PRISM model creation: specification of reactions; instantiation of rates with Grassian et al. metabolic flux rates of isogenic colorectal carcinoma cell lines (with the parental cell line have IDH wild type, and a derived cell line with induction of IDH mutation expression) under normoxia and hypoxia; use of TCGA mRNA expression from patient derived glioma (IDH WT and IDH Mutant) and normal brain to estimate relative rates of metabolic activity for the pairs; specification of properties in probabilistic logic; instantiation of nutrient concentration; model checking of resulting models.

Figure 2

Metabolic phenotypes based on flux rates measured on isogenic cell lines. (A–H) Rates from Grassian et al. were used to instantiate the models of behavior under normoxia (A–D) and hypoxia (E–H). The model checking problem was expressed as characterizing the transient quantity of molecular species at a given instant in time. Glucose and glutamine are initiated at maximum quantity 5 (fmol/cell) and all other species at 0.

Figure 3

Effects of varying nutrient concentration on metabolites under normoxia and hypoxia conditions. (A) Effect of varying inital glucose concentration (from 1 to 5) with initial glutamine set to 5 (fmol/cell). (B) Effect of varying initial glutamine concentration (from 1 to 5) with initial glucose set to 5 (fmol/cell). The PRISM models were instantiated with corresponding phenotype and condition rates as in the experimental models. All other molecular species were initialized with 0. The model checking problem was expressed as characterizing the transient quantity of molecular species at a given instant in time.

Figure 4

Probability of lactate saturation as a function of initial concentration of nutrients. The PRISM models were instantiated with corresponding phenotype and conditions as in the experimental models. Glucose and glutamine initial concentrations are varied over 1–5 (fmol/cell) and all other species are set to 0. The model checking problem was expressed as computing the probability of reaching a state in which lactate is saturated at a given point in time.

Figure 5

Metabolic phenotypes of patient derived samples. (A) Patient-derived mRNA expression data from The Cancer Genome Atlas of genes encoding for TCA and Glycolysis enzymes were processed and aggregated to estimate enzyme levels for three cohorts: Glioblastoma (GBM) samples with IDH-WT phenotype, Lower Grade Glioma (LGG) Samples with IDH-Mutant phenotype, and Normal Brain samples; displayed are across-samples-normalized expression heatmaps. (B) Estimated median expression of the enzymes for each cohort. (C) Using the experimental model flux rates from Grassian et al.  as the predicted GBM IDH WT rates, the flux rates for the other two cohorts were estimated using the ratios of the expression between GBM IDH-WT and LGG IDH-Mutant, and between GBM IDH-WT and Normal Brain. The PRISM models were instantiated with corresponding rates. Glucose and glutamine are initiated at maximum quantity 5 (fmol/cell) and all other species at 0.