Flux Balance Analysis

From bio-physics-wiki

Jump to: navigation, search

Flux Balance Analysis is an approach to find a unique solution to the underdetermined system of equations arising from the Flux Balance Equations. The first step in a FBA is the reconstruction of the metabolic pathways from genome annotation data and literature. Ideally this results in a complete list of reactions occuring in the organism. To illustrate the concept consider the following hypothetical list of reactions


In general one has to be cautious, since there is not always a one to one correspondence between one gene and one reaction. Some enzymes

  • are formed by oligomerisation of gene products of several genes (e.g. Hemoglobin is a tetramer formed by two alpha and two beta globins).
  • catalyse more than one reaction. This behaviour is often called substrate promiscuity.
  • catalyse the same reaction. These enzymes are usually isofomrmes of each other, meaning that only a few amino acids in the protein sequence are changed. Thus these enzymes are also called Isomers. Isomers usually differ in their kinetics (they have different Michaelis Menten konstants $K_m$).

Bearing this in mind we can draw a wiring diagram for the considered network.


The reconstructed metabolic network can be represented in silico by a Stoichiometric Matrix.


The resulting Flux Balance Equations are usually underdetermined and there is no unique solution to the problem, but rather a whole subspace of solutions, the so called Nullspace of the Stoichiometric Matrix. If one splits the reversible reactions into two irreversible reactions (one forward, one backward - a physiological example would be the shared pathway of glycolysis and gluconeogenesis that could be separated, in two reactions), the steady state flux space spanned by the basis vectors form a convex cone. The intersection of the Nullspace with the semipositive orthant of the flux space forms a polyhedral called flux cone.

illustration of the flux cone

In FBA one can find a unique flux solution by assuming that biomass production is maximised. If the biomass composition of an organism is known, a pseudo-reaction can be formulated that produces biomass. Basically all sinks in the organism are collected to form this pseudo reaction. The flux through this biomass reaction also determines the growth rate of the modelled organism.

Given the Stoichiometric Matrix including the biomass composition vector and the objective function $\mathbf{c}$ the flux distribution that maximises the biomass production (flux through pseudo-reaction) can be found by solving the constraint based optimisation problem \begin{equation} Z=max \overset{!}{=}f(\mathbf{v})= \mathbf{c}^T \mathbf{v} \end{equation} subject to the constraints \begin{equation} \mathbf{S} \cdot \mathbf{v}= \mathbf{0} \tag{1} \end{equation} \begin{equation} v_i^{min} \leq v_i \leq v_i^{max} \end{equation} Linear Programming problems of this kind, can be solved with the so called Simplex Algorithm. Software packages like Matlab and Mathematica also support algorithms to solve these problems. Alternatively Solvers are freely available (e.g. GLPK) with GNU-Licence and commercial Solvers (e.g. LINDO). GLPK is also used by the Cobra Toolbox, a software package that allows to perform FBA.

The validity of the assumption of maximised biomass production is based on optimality considerations, see Optimality and Evolution.


To illustrate the Method, FBA is demonstrated at the central carbon metabolism of E. coli for aerobic and anaerobic conditions using Mathematica.


1. Download the Stoichiometric Matrix into a directory of your choice ('yourdirectory') stoichiometric matrix e coli core metabolism.
2. Import the Stoichiometric Matrix into Mathematica with the command:

s = Flatten[Import["yourdirectory\\Stoichiometric_matrix_e_coli_core_metabolism.xls"], 1];

3. The command 'Dimensions[s]' gives you the dimension of the Stoichiometric Matrix which should be $\{72,95\}$
4. In our problem the objective function $\mathbf{c}$ of the optimisation problem \begin{equation} Z=max \overset{!}{=}f(\mathbf{v})= \mathbf{c}^T \mathbf{v} \end{equation}

is simply the zero row vector except having a one at the position of the biomass reaction. Thus $Z$ picks simply the flux $v_{biomass}$ through the biomass pseudo reaction. In our Stoichiometric Matrix the pseudo-reaction is at row 13. Generate the vector $\mathbf{c}$ with the command:
objective = Join[Join[Table[0, {12}], {1}], Table[0, {82}]];

5. In Mathematica Linear Programming problems can be solved with the command

LinearProgramming[$c,m,b,\{ \{l_1,u_1\} ,\{ l_2,u_2\},...\}$]
minimizes c.x subject to the constraints specified by m and b and $l_i<=x_i<=u_i$.
(minimization of negative biomass flux, is the same as maximisation of the biomass flux, thus we have to take the negative objective function!)

6. To use the LinearProgramming feature we need to specify lower and upper bounds ($l_i,u_i$) for the fluxes ($v_i$). Download the constraints for aerobic and anaerobic conditions to 'yourdirectory' and import them into Mathematica with:

fluxconstaerob = Flatten[Import["yourdirectory\\fluxconstaerob.xls"], 1];
fluxconstanaerob = Flatten[Import["yourdirectory\\fluxconstanaerob.xls"], 1];
Under anaerobic conditions the exchange flux $EX_o2(e)$ is constrained to be zero (other bounds are unchanged).
The vector $\mathbf{b}$ is simply the zero vector since we are interested in the steady state reaction rates. Mathematica command:
constraints = Join[Table[{0, 0}, {72}]];

7. Solve the optimisation problems with the command:

fluxaerob = LinearProgramming[-objective, s, constraints, fluxconstaerob] // MatrixForm
fluxanaerob = LinearProgramming[-objective, s, constraints, fluxconstanaerob] // MatrixForm

8. If you wish, you can export the results as an xls-file with the command

Export["yourdirectory\\anaerobicflux.xls", fluxanaerob, "XLS"]
Export["yourdirectory\\aerobicflux.xls", fluxaerob, "XLS"]

Figure shows a flux map with thick edges indicating flux, thin edges indicate no flux. The left map shows the fluxes for aerobic conditions, the map on the right hand side shows fluxes for anaerobic conditions. Glycolysis (Glyc), pentose phosphate pathway (PPP), TCA cycle (TCA), oxidative phosphorylation (OxP), anaplerotic reactions (Ana), and fermentation pathways (Ferm).