E.coli Core Model for Beginners (PART 2)

(please run PART 1 of this tutorial first)

3. Flux Balance Analysis

Flux balance analysis (FBA) is used to calculate the flow of metabolites through a metabolic network making it possible to predict an organism's growth-rate or the production-rate of a bioproduct. Combining the stoichiometric matrix and the objective function can create a system of linear equations that can be used to calculate the fluxes through all the reactions in the network. In flux balance analysis, these equations are solved using linear programming algorithms that can quickly identify optimal solutions to large systems of equations.
Once the external conditions have been set, which include 1) defining the allowed carbon sources, 2) defining the oxygen uptake level, and 3) setting the objective function, then the simulation conditions are setup to perform FBA. This is accomplished through the use of the "optimizeCbModel(model,osenseStr)", a COBRA toolbox function where the first argument is the model name and the second argument determines if the optimization algorithm maximizes ('max') or minimizes ('min') the objective function. Below is an example for an aerobic environment with glucose as the carbon source optimizing for maximum growth-rate. [Timing: Seconds]
model = e_coli_core; % Starting with the original model
model = changeRxnBounds(model,'EX_glc(e)',-10,'l'); % Set maximum glucose uptake
model = changeRxnBounds(model,'EX_o2(e)',-30,'l'); % Set maximum oxygen uptake
model = changeObjective(model,'Biomass_Ecoli_core_w_GAM'); % Set the objective function
FBAsolution = optimizeCbModel(model,'max') % FBA analysis
FBAsolution =
full: [95×1 double] obj: 0.8739 rcost: [95×1 double] dual: [72×1 double] solver: 'gurobi' algorithm: 'default' stat: 1 origStat: 'OPTIMAL' time: 0.7348 basis: [1×1 struct] x: [95×1 double] f: 0.8739 y: [72×1 double] w: [95×1 double] v: [95×1 double]
“FBAsolution” is a Matlab structure that contains the following outputs. “FBAsolution.f “ is the value of objective function as calculated by FBA, thus if the biomass reaction is the objective function then “FBAsolution.f" corresponds to the growth-rate of the cell. In the example above, it can be seen that the growth-rate "FBAsolution.f" is listed as 0.8739 . “FBAsolution.x” is a vector listing the calculated fluxes flowing through the network. “FBAsolution.y” and “FBAsolution.w” contain vectors representing the shadow prices and reduced costs for each metabolite or reaction, respectively.
The flux values found in the structure "FBAsolution.x" can be printed out using the "printFluxVector(model,fluxData,nonZeroFlag,excFlag)" where the second argument is a vector of the flux values, the nonZeroFlag only prints nonzero rows (Default = false), and excFlag only prints exchange reaction fluxes (Default = false). Examples of printing non-zero fluxes and exchange reaction only fluxes are shown below. [Timing: Seconds]
printFluxVector(model,FBAsolution.x,true) % only prints nonzero rows
ACONTa 6.00725 ACONTb 6.00725 AKGDH 5.06438 ATPM 8.39 ATPS4r 45.514 Biomass_Ecoli_core_w_GAM 0.873922 CO2t -22.8098 CS 6.00725 CYTBD 43.599 ENO 14.7161 EX_co2(e) 22.8098 EX_glc(e) -10 EX_h(e) 17.5309 EX_h2o(e) 29.1758 EX_nh4(e) -4.76532 EX_o2(e) -21.7995 EX_pi(e) -3.2149 FBA 7.47738 FRD7 994.936 FUM 5.06438 G6PDH2r 4.95998 GAPD 16.0235 GLCpts 10 GLNS 0.223462 GLUDy -4.54186 GND 4.95998 H2Ot -29.1758 ICDHyr 6.00725 MDH 5.06438 NADH16 38.5346 NH4t 4.76532 O2t 21.7995 PDH 9.28253 PFK 7.47738 PGI 4.86086 PGK -16.0235 PGL 4.95998 PGM -14.7161 PIt2r 3.2149 PPC 2.50431 PYK 1.75818 RPE 2.67848 RPI -2.2815 SUCDi 1000 SUCOAS -5.06438 TALA 1.49698 TKT1 1.49698 TKT2 1.1815 TPI 7.47738
printFluxVector(model,FBAsolution.x,true,true) % only print exchange reaction fluxes
Biomass_Ecoli_core_w_GAM 0.873922 EX_co2(e) 22.8098 EX_glc(e) -10 EX_h(e) 17.5309 EX_h2o(e) 29.1758 EX_nh4(e) -4.76532 EX_o2(e) -21.7995 EX_pi(e) -3.2149
Printing all the zero and nonzero fluxes can be achieved using "printFluxVector(model,FBAsolution.x)."
These fluxes can also be overlayed on a map of the model as shown below, [Timing: Seconds]
map=readCbMap('ecoli_core_map');
options.zeroFluxWidth = 0.1;
options.rxnDirMultiplier = 10;
drawFlux(map, model, FBAsolution.x, options); % Draw the flux values on the map "target.svg"
Document Written
This overlayed map will be written to a file named "target.svg" that should be located in your working directory. Figure 7 is a screenshot of that map.
Figure 7. Screenshot of the network map of the E.coli core model with EX_glc(e) -10 and EX_o2(e) -30 .
As a cautionary note, the default condition for the E.coli core model sets the carbon source as glucose with an uptake rate of -10 , the oxygen uptake is -1000 which implies an aerobic environment with the objectve function defined as 'Biomass_Ecoli_core_w_GAM'. It is a good practice to define the conditions of your simulation explicity to avoid unexpected results and long troubleshooting times.

4. The Subsystems of the E.coli Core Model

Now with these basic Matlab and COBRA toolbox skills behind us, it is time to start exploring the subsytems that make up the E.coli core model. We will start by looking at the "energy production and management" section of the model that is referred to as the "oxidative phosphorylation" subsystem in this core model. This subsystem is located in the upper right corner of the E.coli core map as shown below in Figure 8.