Test physiologically relevant ATP yields from different carbon sources for a metabolic model
Author(s): Ines Thiele, Ronan M. T. Fleming, LCSB, University of Luxembourg.
Reviewer(s):
INTRODUCTION
In this tutorial, we show how to compute the ATP yield from different carbon sources under aerobic or anaerobic conditions. The theoretical values for the corresponding ATP yields are also provided. The tutorial can be adapted for Recon 3 derived condition- and cell-type specific models to test whether these models are still able to produce physiologically relevant ATP yields.
EQUIPMENT SETUP
If necessary, initialize the cobra toolbox with
initCobraToolbox(false) % false, as we don't want to update
For solving linear programming problems in FBA analysis, certain solvers are required:
changeCobraSolver ('glpk', 'all', 1);
> Solver for LP problems has been set to glpk.
> Solver for MILP problems has been set to glpk.
> Solver glpk not supported for problems of type MIQP. Currently used: gurobi
> Solver glpk not supported for problems of type NLP. Currently used: matlab
> Solver glpk not supported for problems of type QP. Currently used: gurobi
This tutorial can be run with GLPK package as linear programming solver, which does not require additional installation and configuration. However, for the analysis of large models, such as Recon 3, it is not recommended to use GLPK, but rather industrial-strength solvers, such as the GUROBI package.
PROCEDURE
Before proceeding with the simulations, the path for the model needs to be set up:
modelFileName = 'Recon2.0model.mat';
modelDirectory = getDistributedModelFolder(modelFileName); %Look up the folder for the distributed Models.
modelFileName= [modelDirectory filesep modelFileName]; % Get the full path. Necessary to be sure, that the right model is loaded
model = readCbModel(modelFileName);
In this tutorial, the used model is the generic model of human metabolism, Recon 3 or Recon2.0 model. The metabolites structures and reactions are from the Virtual Metabolic Human database (VMH, http://vmh.life). Harmonization of abbreviation usage
First, we will harmonize different bracket types used in different model versions, e.g., different version of the human metabolic reconstruction.
model.rxns = regexprep(model.rxns, '\(', '\[');
model.rxns = regexprep(model.rxns, '\)', '\]');
model.mets = regexprep(model.mets, '\(', '\[');
model.mets = regexprep(model.mets, '\)', '\]');
Recon 3 uses ATPSm4mi instead of ATPS4m as an abbreviation for the ATP synthetase: model.rxns = regexprep(model.rxns, 'ATPS4mi', 'ATPS4m');
Similarly, the glucose exchange reaction has been updated:
if length(strmatch('EX_glc[e]', model.rxns))>0
model.rxns{find(ismember(model.rxns, 'EX_glc[e]'))} = 'EX_glc_D[e]';
Add ATP hydrolysis reaction to the model. If the reaction exist already, nothing will be added by the rxnIDexists variable will contain the index of the reaction that is present in the model. In this case, we will rename the reaction abbreviation to ensure that the tutorial works correctly.
[model, rxnIDexists] = addReaction(model, 'DM_atp_c_', 'h2o[c] + atp[c] -> adp[c] + h[c] + pi[c] ');
Warning: Reaction with the same name already exists in the model, updating the reaction
DM_atp_c_ h2o[c] + atp[c] -> adp[c] + h[c] + pi[c]
if length(rxnIDexists) > 0
model.rxns{rxnIDexists} = 'DM_atp_c_'; % rename reaction in case that it exists already
Close model
Now, we will set the lower bound ('model.lb') of all exchange and sink (siphon) reactions to ensure that only those metabolites that are supposed to be taken up are indded supplied to the model.
First, we will find all reactions based on their abbreviation ('model.rxns')
modelexchanges1 = strmatch('Ex_', modelClosed.rxns);
modelexchanges4 = strmatch('EX_', modelClosed.rxns);
modelexchanges2 = strmatch('DM_', modelClosed.rxns);
modelexchanges3 = strmatch('sink_', modelClosed.rxns);
Grab also the biomass reaction(s) based on the reaction abbreviation.
BM= (find(~cellfun(@isempty,strfind(lower(modelClosed.mets), 'bioma'))));
As these measures may not identify all exchange and sink reactions in a model, depending on the used nomencalture, we will also grab all reactions based on stoichiomettry. Here, we will identify all reactions that contain only one non-zero entry in the S matrix (column).
selExc = (find(full((sum(abs(modelClosed.S)==1, 1)==1) & (sum(modelClosed.S~=0) == 1))))';
We will now put all these identified reactions together into one variable 'modelexchanges' and set the lower bound for these reactions to 0.
modelexchanges = unique([modelexchanges1; modelexchanges2; modelexchanges3; modelexchanges4; selExc; BM]);
modelClosed.lb(find(ismember(modelClosed.rxns, modelClosed.rxns(modelexchanges))))=0;
Also, set all upper bounds t 1000 (representing infinity). This may be important if other constraints had been applied to the model, which may interfere with the newly set lower bound of lb=0 for all exchange reactions. Note that this may affect any constraints that had been applied, e.g., condition-specific constraints based on measured uptake or secretion rates.
modelClosed.ub(selExc) = 1000;
Define the ATP hydrolysis reactioDefinen DM_atp_c to be the objective reaction, for which we will maximize for in the following sections. modelClosed = changeObjective(modelClosed, 'DM_atp_c_');
Store the original closed model setup for consequent use in the variable modelClosedOri.
modelClosedOri = modelClosed;
Test for ATP yield from different carbon sources
Now, we re ready to thest for the different individual carbon sources under aerobic and anaerobic conditions for their ATP yield. Therefore, we will provide 1 mol/gdw/hr of a carbon source and maximize the flux through the DM_atp_c_. The results will be stored in the table 'Table_csources'. The table will also contain the theoretical ATP yield, as given by. The table also provides the information of how much flux is going throught he ATP syntheatse. Note that the computed flux distribution is not garantied to be unique, although we use the option 'zero', which approximates the sparsest possible flux distribution with an maximal ATP yield. Carbon source: Glucose (VMH ID: glc_D), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_glc_D[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_glc_D[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{2, 1} = strcat(modelName, ': ATP yield');
Table_csources{3, 1} = strcat(modelName, ': ATPS4m yield');
Table_csources{4, 1} = 'Theoretical';
Table_csources{1, k} = 'glc - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '31';
For this carbon source (glucose), we will also print all reactions that are non-zero in the sparse flux distribution and thus contribute to the maximale ATP yield.
ReactionsInSparseSolution = modelClosed.rxns(find(FBA.x));
We now initiate the next test and delete the variable 'FBA'.
Carbon source: Glucose (VMH ID: glc_D), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns,'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns,'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns,'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns,'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns,'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns,'EX_glc_D[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns,'EX_glc_D[e]'))) = -1;
FBA = optimizeCbModel(modelClosed,'max');
Table_csources{1, k} = 'glc - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns,'ATPS4m'))));
Table_csources{4, k} = '2';
Carbon source: Glutamine (VMH ID: gln_L), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_gln_L[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_gln_L[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'gln_L - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = 'NA';
Carbon source: Glutamine (VMH ID: gln_L), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_gln_L[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_gln_L[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'gln_L - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = 'NA';
Carbon source: Fructose (VMH ID: fru), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_fru[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_fru[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'fru - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3 ,k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '31';
Carbon source: Fructose (VMH ID: fru), Oxygen: No
modelClosed = changeObjective(modelClosed, 'DM_atp_c_');
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_fru[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_fru[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'fru - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '2';
Carbon source: Butyrate (VMH ID: but), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_but[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_but[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'but - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '21.5';
Carbon source: Butyrate (VMH ID: but), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_but[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_but[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'but - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Caproic acid (VMH ID: caproic), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_caproic[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_caproic[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'caproic - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '35.25';
Carbon source: Caproic acid (VMH ID: caproic), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_caproic[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_caproic[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'caproic - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Octanoate (VMH ID: octa), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_octa[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_octa[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'octa - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '49';
Carbon source: Octanoate (VMH ID: octa), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, modelClosed.rxns(modelexchanges)))) = 0;
modelClosed.c = zeros(length(modelClosed.rxns), 1);
modelClosed = changeObjective(modelClosed, 'DM_atp_c_');
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_octa[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_octa[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'octa - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Decanoate (VMH ID: dca), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, modelClosed.rxns(modelexchanges))))=0;
modelClosed.c = zeros(length(modelClosed.rxns),1);
modelClosed = changeObjective(modelClosed, 'DM_atp_c_');
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_dca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_dca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'dca - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '62.75';
Carbon source: Decanoate (VMH ID: dca), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_dca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_dca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'dca - anaerobic';
% fill in only when the LP problem was feasible
Table_csources(2, k) = num2cell(FBA.f);
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Laureate (VMH ID: ddca), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.c = zeros(length(modelClosed.rxns), 1);
modelClosed = changeObjective(modelClosed, 'DM_atp_c_');
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ddca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ddca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ddca - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '76.5';
Carbon source: Laureate (VMH ID: ddca), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ddca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ddca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ddca - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3,k) = num2cell(FBA.x(find(ismember(modelClosed.rxns,'ATPS4m'))));
Table_csources{4,k} = '0';
Carbon source: Tetradecanoate (VMH ID: ttdca), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ttdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ttdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ttdca - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '90.25';
Carbon source: Tetradecanoate (VMH ID: ttdca), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ttdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ttdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ttdca - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns,'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Hexadecanoate (VMH ID: hdca), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_hdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_hdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'hdca - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '104';
Carbon source: Hexadecanoate (VMH ID: hdca), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_hdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_hdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'hdca - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Octadecanoate (VMH ID: ocdca), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ocdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ocdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ocdca - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '117.75';
Carbon source: Octadecanoate (VMH ID: ocdca), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_ocdca[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_ocdca[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'ocdca - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Arachidate (VMH ID: arach), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_arach[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_arach[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'arach - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '131.5';
Carbon source: Arachidate (VMH ID: arach), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_arach[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_arach[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'arach - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Behenic acid (VMH ID: docosac), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_docosac[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_docosac[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'docosac - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '145.25';
Carbon source: Behenic acid (VMH ID: docosac), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_docosac[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_docosac[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'docosac - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Lignocerate (VMH ID: lgnc), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_lgnc[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_lgnc[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'lgnc - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3 ,k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '159';
Carbon source: Lignocerate (VMH ID: lgnc), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_lgnc[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_lgnc[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'lgnc - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
Carbon source: Cerotate (VMH ID: hexc), Oxygen: Yes
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = -1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_hexc[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_hexc[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'hexc - aerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3, k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '172.75';
Carbon source: Cerotate (VMH ID: hexc), Oxygen: No
modelClosed = modelClosedOri;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_o2[e]'))) = 0;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = -1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_h2o[e]'))) = 1000;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_co2[e]'))) = 1000;
modelClosed.lb(find(ismember(modelClosed.rxns, 'EX_hexc[e]'))) = -1;
modelClosed.ub(find(ismember(modelClosed.rxns, 'EX_hexc[e]'))) = -1;
FBA = optimizeCbModel(modelClosed, 'max');
Table_csources{1, k} = 'hexc - anaerobic';
Table_csources(2, k) = num2cell(FBA.f);
% fill in only when the LP problem was feasible
% set all flux values less than tol to 0
FBA.x(find(abs(FBA.x)<=tol)) = 0;
Table_csources(3 ,k) = num2cell(FBA.x(find(ismember(modelClosed.rxns, 'ATPS4m'))));
Table_csources{4, k} = '0';
The table contains all computed ATP yields.
Table_csources = Table_csources'
Warning: A value of class "com.mathworks.mde.cmdwin.XCmdWndView" was indexed with no subscripts specified. Currently the result of this operation is the indexed value itself, but in a future release, it will be an error.
Warning: A value of class "com.mathworks.mde.cmdwin.XCmdWndView" was indexed with no subscripts specified. Currently the result of this operation is the indexed value itself, but in a future release, it will be an error.
TIMING
This tutorial takes only a few minutes.
REFERENCES
[1] Brunk, E. et al. Recon 3D: A Three-Dimensional View of Human Metabolism and Disease. Submited