Tomer Shlomi School of Computer Science, Tel-Aviv University, Tel-Aviv, Israel January, 2006

1. Constraint-Based Modeling of Metabolic Networks based on: “Genome-scale models of microbial cells:Evaluating the consequences of constraints”, Price, et. al (2004) Tomer Shlomi School of Computer Science, Tel-Aviv University, Tel-Aviv, Israel January, 2006

2. Outline • Metabolism and metabolic networks • Kinetic models vs. constraints-based modeling • Flux Balance Analysis • Exploring the solution space • Altering phenotypic potential: gene knockouts

3. Cellular Metabolism • The essence of life.. • Catabolism and anabolism • The metabolic core – production of energy – anaerobic and aerobic metabolism • Probably the best understood of all cellular networks: metabolic, PPI, regulatory, signaling • Tremendous importance in Medicine; antibiotics, metabolic disorders, liver disorders, heart disorders • Bioengineering; efficient production of biological products.

4. Metabolites and Biochemical Reactions • Metabolite: an organic substance, e.g. glucose, oxygen • Biochemical reaction: the process in which two or more molecules (reactants) interact, usually with the help of an enzyme, and produce a product Glucose + ATP Glucokinase Glucose-6-Phosphate + ADP

6. Kinetic Models • Dynamics of metabolic behavior over time • Metabolite concentrations • Enzyme concentrations • Enzyme activity rate – depends on enzyme concentrations and metabolite concentrations • Solved using a set of differential equations • Impossible to model large-scale networks • Requires specific enzyme rates data • Too complicated

7. Constraint Based Modeling • Provides a steady-state description of metabolic behavior • A single, constant flux rate for each reaction • Ignores metabolite concentrations • Independent of enzyme activity rates • Assume a set of constraints on reaction fluxes • Genome scale models Flux rate: μ-mol / (mg * h)

8. Constraint Based Modeling • Find a steady-state flux distribution through all biochemical reactions • Under the constraints: • Mass balance: metabolite production and consumption rates are equal • Thermodynamic: irreversibility of reactions • Enzymatic capacity: bounds on enzyme rates • Availability of nutrients

9. Metabolic Networks Biochemistry Cell Physiology Genome Annotation Inferred Reactions Network Reconstruction Analytical Methods Metabolic Network

10. Mathematical Representation • Stoichiometric matrix – network topology with stoichiometry of biochemical reactions Glucokinase Glucose + ATP Glucokinase Glucose-6-Phosphate + ADP Glucose -1 ATP -1 G-6-P +1 ADP +1 Mass balance S·v = 0 Subspace of R Thermodynamic vi > 0 Convex cone Capacity vi < vmax Bounded convex cone n

11. Growth Medium Constraints • Exchange reactions enable the uptake of nutrients from the media and the secretion of waste products Lower bound Upper bound Glucose 0 2.5 Oxygen 0 Inf CO2 -Inf 0 G-Ex O-Ex Co2-Ex Glucose 1 Oxygen 1 CO2 1

12. Determination of Likely Physiological States • How to identify plausible physiological states? • Optimization methods • Maximal biomass production rate • Minimal ATP production rate • Minimal nutrient uptake rate • Exploring the solution space • Extreme pathways • Elementary modes

13. Outline: Optimization Methods • Predicting the metabolic state of a wild-type strain • Flux Balance Analysis (FBA) • Predicting the metabolic state after a gene knockout • Minimization Of Metabolic Adjustment • Regulatory On/Off Minimization

14. Biomass Production Optimization • Metabolic demands of precursors and cofactors required for 1g of biomass of E. coli • Classes of macromolecules: Amino Acids, Carbohydrates Ribonucleotides, Deoxyribonucleotides Lipids, Phospholipids Sterol, Fatty acids • These precursors are removed from the metabolic network in the corresponding ratios • We define a growth reaction Z = 41.2570 VATP - 3.547VNADH+18.225VNADPH + ….

15. Biomass Composition Issues • Varies across different organisms • Depends on the growth medium • Depends on the growth rate • The optimum does not change much with changes in composition within a class of macromolecules • The optimum does change if the relative composition of the major macromolecules changes

16. growth Flux Balance Analysis (FBA) • Successfully predicts: • Growth rates • Nutrient uptake rates • Byproduct secretion rates • Solved using Linear Programming (LP) • Finds flux distribution with maximal growth rate Max vgro, - maximize growth s.t S∙v = 0, - mass balance constraints vmin  v  vmax - capacity constraints Fell, et al (1986), Varma and Palsson (1993)

17. FBA Example (1)

18. FBA Example (2)

19. FBA Example (2)

20. Linear Programming Basics (1)

21. Linear Programming Basics (2)

22. Linear Programming Basics (3)

23. Linear Programming: Types of Solutions (1)

24. Linear Programming: Types of Solutions (2)

25. Linear Programming Algorithms • Simplex • Used in practice • Does not guarantee polynomial running time • Interior point • Worse case running time is polynomial growth

26. Phenotype Predictions: Evolving Growth Rate

28. Exploring the Convex Solution Space

29. growth growth growth Alternative Optima • The optimal FBA solution is not unique One solution Optimal solutions Near-optimal solutions • Basic solutions enumeration – MILP (Lee, et. al, 2000) • Flux variability analysis (Mahadevan, et. al. 2003) • Hit and run sampling (Almaas, et. al, 2004) • Uniform random sampling (Wiback, et. al, 2004)

30. What Do Multiple Solutions Represent ? • Some of the solutions probably do not represent biologically meaningful metabolic behaviors as there are missing constraints • Previous studies tackled this problem by: • Incorporating additional constraints: regulatory constraints (Covert, et. al., 2004) • Looking for reactions for which new constraints may significantly reduce the solution space (Wiback, et. al., 2004) FBA solution space Meaningful solutions

31. Interpretations of Metabolic Space • Effect of exogenous factors – the metabolic space corresponds to growth in a medium under various external conditions that are beyond the model’s scope such as stress or temperature • Heterogeneity within a population - the metabolic space represents heterogenous metabolic behaviors by individuals within a cell population (Mahadevan, et. al., 2003, Price, et. al., 2004) • Alternative evolutionary paths – the metabolic space represents different metabolic states attainable through different evolutionary paths (Mahadevan, et. al., 2003, Fong, et. al., 2004) • The three interpretations are obviously not mutually exclusive

32. Alternative Optima: Basic Solutions Enumeration • Lee, et. al, 2000 • Basic solutions – metabolic states with minimal number of non-zero fluxes • Different solutions differ in at least a single zero flux • Use Mixed Integer Linear Programming • Formulate optimization as to identify new solutions that are different from the previous ones • Applicable only to small scale models growth

33. Alternative Optima: Flux Variability Analysis • Mahadevan, et. al. 2003 • Find metabolic states with extreme values of fluxes • Use linear programming to minimize and maximize the flux through each reaction while satisfying all constraints Max / Min vi, - maximize growth s.t S∙v = 0, - mass balance constraints vmin  v  vmax - capacity constraints Vgro = Vopt - set maximal growth rate

34. Alternative Optima: Hit and Run Sampling • Almaas, et. al, 2004 • Based on a random walk inside the solution space polytope • Choose an arbitrary solution • Iteratively make a step in a random direction • Bounce off the walls of the polytope in random directions

35. Alternative Optima: Uniform Random Sampling • Wiback, et. al, 2004 • The problem of uniform sampling a high-dimensional polytope is NP-Hard • Find a tight parallelepiped object that binds the polytope • Randomly sample solutions from the parallelepiped • Can be used to estimate the volume of the polytope

36. Topological Methods • Not biased by a statement of an objective • Network based pathways: • Extreme Pathways (Schilling, et. al., 1999) • Elementary Flux Modes (Schuster, el. al., 1999) • Decomposing flux distribution into extreme pathways • Extreme pathways defining phenotypic phase planes • Uniform random sampling

37. Extreme Pathways andElementary Flux Modes • Unique set of vectors that spans a solution space • Consists of minimum number of reactions • Extreme Pathways are systematically independent (convex basis vectors)

38. Extreme Pathways andElementary Flux Modes • Inherent redundancy in metabolic networks (Price, et. al., 2002) • Robustness to gene deletion and changes in gene expression (Stelling, et. al., 2002) • Enzyme subsets (correlated reaction sets) in yeast (Papin, et. al., 2002) • Design strains (Carlson, et. al., 2002) • Assign functions to genes (Forster, et. al, 2002)

39. Altering Phenotypic Potential: Gene Knockouts

40. w v Altering Phenotypic Potential: Gene Knockouts • Minimization Of Metabolic Adjustment (MOMA) (Segre et. al, 2002) • The flux distribution after a knockout is close to the wild-type’s state under the Euclidian norm • Regulatory On/Off Minimization (ROOM) (Shlomi et. al, 2005) • Minimize the number of Boolean flux changes from the wild-type’s state

41. Altering Phenotypic Potential • Explaining gene dispensability (Papp, el. al., 2004) • Only 32% of yeast genes contribute to biomass production in rich media • Considered one arbitrary optimal growth solution • OptKnock – Identify gene deletions that generate desired phenotype (Burgard, et. al., 2003) • OptStrain – Identify strains which can generate desired phenotypes by adding/deleting genes (Pharkya, el., al., 2004)

42. Modeling Gene Knockouts • Gene knockout • Enzyme knockout • Reaction knockout

43. growth generations minutes Cellular Adaptation to Genetic and Environmental Perturbations • Transient changes in expression levels in hundreds of genes (Gasch 2000, Ideker 2001) • Convergence to expression steady-state close to the wild-type (Gasch 2000, Daran 2004, Braun 2004) • Drop in growth rates followed by a gradual increase (Fong 2004)

44. w v Regulatory On/Off Minimization (ROOM) • Predicts the metabolic steady-state following the adaptation to the knockout • Assumes the organism adapts by minimizing the set of regulatory changes Boolean Regulatory Change Boolean Flux Change • Finds flux distribution with minimal number of Boolean flux changes

45. ROOM: Implementation • Solved using Mixed Integer Linear Programming (MILP) • Boolean variable yi yi = 1 Flux vi change from wild-type Min yi - minimize changes s.t v – y ( vmax - w)  w - distance constraints v – y ( vmin - w)  w - distance constraints • S∙v = 0, - mass balance constraints • vj = 0, jG - knockout constraints • MILP is NP-Hard • Relax Boolean constraints - solve using LP • Relax strict constraint of proximity to wild-type

46. Example Network

47. ROOM’s Implicit Growth Rate Maximization • ROOM implicitly attempts to maintain the maximal possible growth rate of the wild-type organism • A change in growth requires numerous changes in fluxes M1 M2 Growth Reaction . . Biomass Mn

48. Intracellular Flux Measurements • Intracellular fluxes measurements in E. coli central carbon metabolism • Obtainedusing NMR spectroscopy in C labeling experiments • 5 knockouts: pyk, pgi, zwf, gnd, ppc in Glycolysis and Pentose Phosphate pathways • Glucose limited and Ammonia limited medias • FBA wild-type predictions above 90% accuracy 13 Emmerling, M. et al. (2002), Hua, Q. et al. (2003), Jiao, Z et al. (2003), Peng, et. al (2004)

49. Knockout Flux Predictions • ROOM flux predictions are significantly more accurate than MOMA and FBA in 5 out of 9 experiments • ROOM steady-state growth rate predictions are significantly more accurate than MOMA

50. ROOM vs. MOMA • ROOM predicts metabolic steady-state after adaptation • Provides accurate flux predictions • Preserved flux linearity • Finds alternative pathways • Predicts steady-state growth rates • MOMA predicts transient metabolic states following the knockout • Provides more accurate transient growth rates