Metabotools¶
- analyzeSingleGeneDeletion(ResultsAllCellLines, path, samples, cutoff, OverViewResults)[source]¶
This function conducts a single gene deletion analysis and provides statistcs about the frequency of essential genes, and partial effects and no effects across the unique set of genes and across the models in ResultsAllCellLines defined by samples.
- USAGE
[genes, ResultsAllCellLines, OverViewResults] = analyzeSingleGeneDeletion (ResultsAllCellLines, path, samples, cutoff, OverViewResults)
- INPUTS
ResultsAllCellLines – Structure containing the models (the two models can be the same or pruned and starting model) of the format
ResultsAllCellLines.samples.modelMin
ResultsAllCellLines.samples.modelPruned
path – Path where picture of heatmap is to be saved
samples – Name of samples
cutoff – Tolerance for small deviation from the optimal value whic his used to distinguish partial effects
- OPTIONAL INPUTS
OverViewResults – Overview of results of the model building and analysis (generated by the function setQuantConstraints)
- OUTPUTS
genes – Table that lists unique gene, category, number of ‘no effect’, number of ‘KO’, number of ‘partial effect’, and number of models where the gene is absent ‘gene absent’
ResultsAllCellLines – Updated structure ResultsAllCellLines.sample.singleGeneDeletion
- OPTIONAL OUTPUTS:
OverViewResults: Updated overview
- calculateLODs(theo_mass, lod_ngmL)[source]¶
This function converts detection limits of the unit ng/mL to mM using the theoretical mass (g/mol) for the metabolites
- USAGE
[lod_mM] = calculateLODs (theo_mass, lod_ngmL)
- INPUTS
theo_mass – Vector that specifies the theoretical mass (g/mol) of each metabolite
lod_ngmL – Instrumental limit of detection (ng/mL)
- OUTPUT
lod_mM – Detection limits in mM
- calculateQuantitativeDiffs(data_RXNS, slope_Ratio, ex_RXNS, lod_mM, cond1_uptake, cond2_uptake, cond1_secretion, cond2_secretion)[source]¶
This function provides sets of exchange reactions with higher uptake and secretion in condition 1 and condition 2.
- USAGE
[cond1_upt_higher, cond2_upt_higher, cond2_secr_higher, cond1_secr_higher, cond1_uptake_LODs, cond2_uptake_LODs, cond1_secretion_LODs, cond2_secretion_LODs] = calculateQuantitativeDiffs (data_RXNS, slope_Ratio, ex_RXNS, lod_mM, cond1_uptake, cond2_uptake, cond1_secretion, cond2_secretion)
- INPUTS
data_RXNS – Exchange reactions same order as Input_A and Input_B
slope_Ratio – For length of inputA/B, relative difference
ex_RXNS – Exchange reactions in the same order as LODs(lod_mM)
lod_mM – Detection limits in mM
cond1_uptake – List of exchanges that specify consumed metabolites in condition 1
cond2_uptake – List of exchanges that specify consumed metabolites in condition 2
cond1_secretion – List of exchanges that specify released metabolites in condition 1
cond2_secretion – List of exchanges that specify released metabolites in condition 2
- OUTPUTS
cond1_upt_higher – Exchange reactions and relative differences with higher uptake in condition 1 compared to condition 2
cond2_upt_higher – Exchange reactions and relative differences with higher uptake in condition 2 compared to condition 1
cond2_secr_higher – Exchange reactions and relative differences with higher secretion in condition 2 compared to condition 1
cond1_secr_higher – Exchange reactions and relative differences with higher secretion in condition 1 compared to condition 2
cond1_uptake_LODs – Detection limits for metabolites with higher uptake in condition 1
cond2_uptake_LODs – Detection limits for metabolites with higher uptake in condition 2
cond1_secretion_LODs – Detection limits for metabolites with higher secretion in condition 1
cond2_secretion_LODs – Detection limits for metabolites with higher secretion in condition 2
- checkEffectRxnKO(samples_to_test, fill, Genes_to_test, samples, ResultsAllCellLines)[source]¶
This function checks the effect of constraining reactions associated with a single or set of genes on the ability of the model to satisfy an objective.
- USAGE
[FBA_Rxns_KO, ListResults] = checkEffectRxnKO (samples_to_test, fill, Genes_to_test, samples, ResultsAllCellLines)
- INPUTS
ResultsAllCellLines – uses modelMin
samples – Name of samples
samples_to_test – Name of samples that should be tested (can be samples if all should be tested)
fill – Identifier if the rxns is not in the model (e.g.,100, num(‘NAN’))
Genes_to_test – Set of genes to be tested
- OUTPUTS
FBA_Rxns_KO – FBA results for constraining one reaction at a time to zero.
ListResults – Reactions associated with Genes_to_test, same order as FBA_Rxns_KO.
- checkExchangeProfiles(samples, path, nmets)[source]¶
The Function generates a summary of the number of uptake exchanges and secretion exchanges per samples.
- USAGE
[mapped_exchanges, minMax, mapped_uptake, mapped_secretion] = checkExchangeProfiles (samples, path, nmets)
- INPUTS
samples
path – path to output of make exchangeprofiles
nmets – number of metabolites in data set (max number of uptake or secreted metabolites)
- OUTPUTS
mapped_exchanges – table listing for each sample the number of uptake and secretion that were mapped to the model
minMax – lists minimal and maximal number of uptakes, secretions, and total number of exchanges
mapped_uptake – table summarizing for each sample the uptake exchange reactions
mapped_secretion – table summarizing for each sample the secretion exchange reactions
- computeFluxSplits(model, mets, V, coeffSign)[source]¶
Compute relative contributions of fluxes (V) to the net production (P) and consumption (C) of a set of metabolites (mets)
- USAGE
[P, C, vP, vC, s] = computeFluxSplits (model, mets, V, 0)
- INPUTS
model – a COBRA model structure with required fields
.S: the m x n stoichiometric matrix
.mets: an m x 1 cell array of metabolite identifiers
mets – a list of metabolite identifiers from model.mets
V – an n x k matrix of k flux distributions
- OPTIONAL INPUT
coeffSign – only use the sign of the stoichiometric coefficient (i.e. +1 or -1, Default = false)
- OUTPUTS
P – an n x k matrix of relative contributions to the production of mets
C – an n x k matrix of relative contributions to the consumption of mets
vP – an n x k matrix of net producing fluxes
vC – an n x k matrix of net consuming fluxes
s – net stoichiometry of reactions containing metabolites
Example: [P,C,vP,vC,s] = computeFluxSplits(model,mets,V,1);
- conc2Rate(metConc, cellConc, t, cellWeight)[source]¶
Converts metabolite concentration and (viable) cell concentration into uptake rate. CellConc consumed MetConc in t.
- USAGE
[flux] = conc2Rate (metConc, cellConc, t, cellWeight)
- INPUTS
metConc – Change in metabolite concentration (mM)
cellConc – Cell concentration (cells per 1 ml)
t – Time in hours
cellWeight – gDW per cell
- OUTPUT
flux – mmol/gDW/hr
- defineUptakeSecretionProfiles(input_A, input_B, data_RXNS, tol, essAA_excl, exclude_upt, exclude_secr, add_secr, add_upt)[source]¶
This function calculated the slope ratios and gives out uptake&secretion profiles for condition 1 and condition 2
- USAGE
[cond1_uptake, cond2_uptake, cond1_secretion, cond2_secretion, slope_Ratio] = defineUptakeSecretionProfiles (input_A, input_B, data_RXNS, tol, essAA_excl, exclude_upt, exclude_secr, add_secr, add_upt)
- INPUTS
input_A – Matrix, 4 columns
Control TP1,
Control TP2,
Condition 1 TP1,
Condition 1 TP 2
input_B – Matrix, 4 columns
Control2 TP1,
Control2 TP2,
Condition 2 TP1,
Condition 2 TP 2
data_RXNS – Exchange reactions same order as input_A and input_B
tol – minimal threshold to call uptake or secretion (Default 0.05)=5%
- OPTIONAL INPUTS
essAA_excl – If essential amino acids should be excluded from the secretion profile 1(yes), or 0(no) (Default = 0); essAAs = {EX_his_L(e); EX_ile_L(e); EX_leu_L(e); EX_lys_L(e); EX_met_L(e); EX_phe_L(e); EX_thr_L(e); EX_trp_L(e); EX_val_L(e)}
exclude_upt – Exclude uncertain metabolites from uptake (e.g., metabolites from GlutaMax, e.g., EX_gln_L(e), EX_cys_L(e), EX_ala_L(e))
exclude_secr – Exclude uncertain metabolites from secretion (e.g., metabolites from GlutaMax, e.g., EX_gln_L(e), EX_cys_L(e), EX_ala_L(e))
add_secr – Due to mising data points automatic analysis might do wrong assignment of metabolites to secretion
add_upt – Due to mising data points automatic analysis might do wrong assignment of metabolites to uptake
- OUTPUTS
cond1_uptake – List of exchanges that specify consumed metabolites in condition 1
cond2_uptake – List of exchanges that specify consumed metabolites in condition 2
cond1_secretion – List of exchanges that specify released metabolites in condition 1
cond2_secretion – List of exchanges that specify released metabolites in condition 2
slope_Ratio – For length of inputA/B, relative difference
- extractConditionSpecificModel(model, threshold)[source]¶
The function prunes a subnetwork based on a user-defined threshold. The subnetwork does not contain blocked reactions. Please note that Recon has blocked reactions which will always be removed. Thus, not all reactions are removed as a consequence of the data integration. Depends on fastFVA, fluxVariability analysis
- USAGE
[modelPruned] = extractConditionSpecificModel (model, threshold)
- INPUTS
model – Global metabolic model (Recon) - constrained
threshold – Fluxes below the threshold will be considered zero and respective reactions as blocked, e.g., 10e-6.
- OUTPUTS
modelPruned – submodel without blocked reactions
- findMinCardModel(model, Ex_Rxns)[source]¶
Function to find the minimal cardinality model.
- USAGE
[modelMin, AddedExchange] = findMinCardModel (model, Ex_Rxns)
- INPUTS
model – Metabolic model
Ex_Rxns – Vector of exchange reactions for FVA
- OUTPUTS
modelMin – Updated model with new reaction bounds
AddedExchange – Vector of exchanged reactions
- findOptExchRxns(model, Ex_Rxns, parFVA)[source]¶
- USAGE
[OptExchRxns] = findOptExchRxns (model, Ex_Rxns, parFVA)
- INPUTS
model – Metabolic model
Ex_Rxns – Vector of exchange reactions for FVA
fastFVA – use FastFVA (default=0)
- OUTPUTS
OptExchRxns
- generateCompactExchModel(model, minGrowth, biomassRxn, prune, fastFVA)[source]¶
This function identifies a subnetwork with the least number of possible exchange reactions given the model and the applied constraints. It returns the resulting pruned model.
- USAGE
[modelMin, modelPruned, Ex_Rxns] = generateCompactExchModel (model, minGrowth, biomassRxn, prune, fastFVA)
- INPUTS
model – model structure
minGrowth – minimal Growth rate to be set on biomass reaction
biomassRxn – biomass reaction name (default: ‘biomass_reaction2’)
prune – optional: to prune the model based on exchange reactions (default: 1)
fastFVA – optional: to use fastFVA instead of fluxvariability for computing FVA results (default: 0)
medium – (default: {})
- OUTPUTS
modelUpdated – same as input model but constraints on blocked reactions are set to be 0
modelPruned – pruned model, where all blocked reactions are removed (attention this seems to cause issues with GPRs)
Ex_Rxns – List of exchange reactions in pruned model
- illustrate_ppp(ResultsAllCellLines, mets, path, samples, label, fonts, tol)[source]¶
This function generates and saves heatmaps for the results of the function performPPP for all sample models.
- USAGE
illustrate_ppp (ResultsAllCellLines, mets, path, samples, label, fonts, tol)
- INPUTS
ResultsAllCellLines – Result structure
mets – Metabolites that were tested in the phase plane analysis
step_size – Step size of each metabololite tested
path – Path where output is saved
samples – Names of conditions
label – Defining label of X-axis, y-axis and z-axis, e.g., {Glucose uptake (fmol/cell/hr); Oxygen uptake (fmol/cell/hr); Growth rate (hr-1)} The z-axis is the growth rate which is color coded
fonts – Font size for labels on heatmap
tol – Fluxes are considered zero if they are below the tol
Individual pdf files showing the result of the PPP are being saved automatically for each condition.
- integrateGeneExpressionData(model, dataGenes)[source]¶
This function sets constrains based on sets of absent genes. It does not test for the functionality of the model, to maintain full control over the model generation process. Alternatively, the COBRA function createTissueSpecificModel can be used, and which can produce a functional model. See the supplemental tutorial of DOI:10.1371/journal.pone.0049978 (Aurich and Thiele, PloS One, 2012) for description on potentially necessary curation when using this crude method of data integration.
- USAGE
[modelGE] = integrateGeneExpressionData (model, dataGenes)
- INPUTS
model – Metabolic model (e.g., Recon)
dataGenes – Vector of absent genes. Follow the supplemental tutorial of DOI: 10.1371/journal.pone.0049978 (Aurich and Thiele, PloS One, 2012) for the generation of P/A calls, e.g., Absent_genes = [535;1548];
- OUTPUT
modelGE – Model, where all genes in DataGenes have been constrained to zero. Follow the supplemental tutorial of DOI:10.1371/journal.pone.0049978 (Aurich and Thiele, PloS One, 2012) for description of potentially necessary curation.
- make3Dplot(PHs, maximum_contributing_flux, fonts, path, diff_view)[source]¶
The function allows illustration on different phenotypes predicted in MetTB onto the 3D illustration plot. The phenotypes are defined as colors, the location of the points is based on the ATP flux predicted by the function predictFluxSplits.
- USAGE
make3Dplot (PHs, maximum_contributing_flux, fonts, path, diff_view)
- INPUTS
PHs – Phenotypes = [samples PHENOTYPE], is generated from different functions (predictFluxSplits,…)
maximum_contributing_flux – Output of the function predictFluxSplits, option ATP.
fonts – Font size
path – Path to save output
diff_view – Next to the standard view, print graphic from 3 different(2D) viewpoints (default = 0/off)
PDFs are automatically saved to the user-defined path
- makeSummaryModels(ResultsAllCellLines, samples, model, mk_union, mk_intersect, mk_reactionDiff)[source]¶
This function generates the union and the intersect model from the modelPruned in the ResultsAllCellLines structure.
- USAGE
[unionModel, intersectModel, diffRxns, diffExRxns] = makeSummaryModels (ResultsAllCellLines, samples, model, mk_union, mk_intersect, mk_reactionDiff)
- INPUTS
ResultsAllCellLines – structure containing samples and models for the samples, e.g., ResultsAllCellLines.UACC_257.modelPruned
samples – conditions or cell lines, e.g., UACC_257
- OPTIONAL INPUTS
mk_union – make union model, yes=1, no=0 (Default = 1)
mk_intersect – make intersect model, yes=1, no=0 (Default = 1)
mk_reactionDiff – make reactionDiff, yes=1, no=0 (can only be 1 if union and intersect are 1 or []) (Default = 1).
- OUTPUTS
unionModel – model containing all reactions appearing at least once in the models in the ResultsAllCellLines.sample.modelPruned
intersectModel – model containing all reactions shared by all models in the ResultsAllCellLines.sample.modelPruned
diffRxns – all differential reactions that distinguish unionModel and intersectModel
diffExRxns – all differential exchange reactions that distinguish unionModel and intersectModel
- mkTableOfAddedExchanges(ResultsAllCellLines, samples, Ex_added_all_unique)[source]¶
The function generates a table of the added exchanges defined by the function generateCompactExchModel called by the function setQuantConstraints.
- USAGE
[Added_all] = mkTableOfAddedExchanges (ResultsAllCellLines, samples, Ex_added_all_unique)
- INPUTS
ResultsAllCellLines
samples
Ex_added_all_unique – Output of the function statisticsAddedExchanges
- OUTPUT
Added_all – Table overview of all added Exchanges for each sample
- networkTopology(model)[source]¶
Analysis of the metabolite connectivity of a metabolic model
- USAGE
[MetConn, RxnLength] = NetworkTopology (model)
- INPUTS
model – Model structure
- OUTPUTS
MetConn – Vector of metabolite connectivity (how many reactions a metabolite participates in (in same order as model.mets)
RxnLength – Vector of reaction participation, i.e., how many metabolites per reaction (in same order as mode.rxns)
- performPPP(ResultsAllCellLines, mets, step_size, samples, step_num, direct)[source]¶
This function performs the a phase plane analysis. The analysis starts from zero and proceed step_num*step_size in the specified direction and for two exchanges. the results of the analysis will be saved into ResultsAllCellLines.
- USAGE
[ResultsAllCellLines] = performPPP (ResultsAllCellLines, mets, step_size, samples, step_num, direct)
- INPUTS
ResultsAllCellLines
samples – Conditions
mets – Set of two metabolites to test e.g. mets ={EX_glc(e), EX_gln(e); EX_o2(e), EX_o2(e)}
step_size – Step size for metabolites specified in mets e.g., step_size = 100
direct – Direction for metabolites specified in mets: uptake (-1) or secretion (1), e.g., direct = [-1,-1,-1,-1];
step_num – Number of steps (1000/20 = 50), step_num = [5,5;5,5]
- OUTPUT
ResultsAllCellLines – The matrix of growth rates are added to the ResultsAllCellLines structure alond with the values or the bounds using the same of the manipulated exchanges
- performSampling(model, warmupn, fileName, nFiles, pointsPerFile, stepsPerPoint, fileBaseNo, maxTime, path)[source]¶
The function performs the sampling analysis for one model, combines multiple COBRA toolbox functions
- USAGE
performSampling (model, warmupn, fileName, nFiles, pointsPerFile, stepsPerPoint, fileBaseNo, maxTime, path)
- INPUTS
model – Model to be sampled
warmupn – Number of warm-up points
fileName – Name for output files
nFiles – Number of files saved
pointsPerFilePoints – Points saved per file
stepsPerPoint – Steps skipped between two saved points, in order to increase mixing
fileBaseNo – Counter for the numbering of the output files, e.g., 0 to start with 1;
maxTime – Maximal running time after which the analysis should be terminated, e.g., 3600000;
path – Path to file
Automatically saved to path using specified file name
- predictFluxSplits(model, obj, met2test, samples, ResultsAllCellLines, dir, transportRxns, ATPprod, carbon_source, eucNorm)[source]¶
This function performs the flux splits analysis for the metabolites of interest, meaning it predicts the fraction of metabolite produced (or consumed) based all reactions producing (or consuming) the metabolite
- USAGE
[BMall, ResultsAllCellLines, metRsall, maximum_contributing_rxn, maximum_contributing_flux, ATPyield] = predictFluxSplits (model, obj, met2test, samples, ResultsAllCellLines, dir, transportRxns, ATPprod, carbon_source, eucNorm)
- INPUTS
model – Generic model, e.g., modelMedium
obj – objective function, e.g., biomass or ATPM
met2test – e.g., atp[c]. Mind that the metabolites are produced in multiple compartments.
samples – Name of conditions as in ResultsAllCellLines
ResultsAllCellLines – Structure containing the pruned submodels of the samples
- OPTIONAL INPUTS
dir – Production = 1 (default = production), else consumption = 0
eucNorm – Default: 1e-6
transportRxns – Vector of reactions that do not really produce or consume the metabolite, e.g., reactions that transport ATP from one compartment to the other. This input is optional only to allow the initial prediction of all producing reactions to define the exclude reaction set.
carbon_source – Reference uptake for calculation of ATP yield, e.g.,{EX_glc(e)}.
ATPprod
- OUTPUTS
BMall – Matrix of flux vectors used for calculations
ResultsAllCellLines – Structure containing results of run analysis
metRsall – Matrix of flux (producing or consuming) a defined metabolite
maximum_contributing_rxn – Reactions with highest flux (producing or consuming) a defined metabolite across analyzed samples (not necessarily >50%), if multiple reactions have the same contribution, all will be reported seperated by a back slash.
maximum_contributing_flux – Matrix containing:
highest flux (column 1),
sum of flux (producing or consuming) a defined metabolite (column 2),
percentage (column 3),
contribution of glycolysis (column 4),
contribution of ETC (column 5),
combined contribution of glycolysis and ETC (column 6),
contribution of TCA (column 7),
combined contribution of glycolysis, ETC, and TCA (column 8).
ATPyield – ATP yield calculated from the sum of ATP production divided by the predicted uptake flux of the metabolite specified as carbon_source. No extra constraints are applied, thus not only production flux from the specified carbon source is considered.
Example
% if met2test is not atp to reduce number of useless outputs [BMall, ResultsAllCellLines, metRsall] = predictFluxSplits(model, obj, met2test,samples,ResultsAllCellLines, dir, eucNorm, transportRxns, ATPprod, carbon_source)
- prepIntegrationQuant(model, metData, exchanges, samples, test_max, test_min, outputPath, tol, variation)[source]¶
This function generates individual uptake and secretion profiles from a data matrix (fluxes) with samples as columns and metabolites as rows. Negative values are interpreted as uptake and positive values are interpreted as secretion. The function tests based on the input model if uptake and secretion of all ‘Exchanges’ is possible. Subsequently, it removes from each sample the metabolite uptakes or secretions that cannot be realized by the model due to missing production or degradation pathways, or blocked reactions. If only secretion is not possible, only secretion is eliminated from the sample profile whereas uptake will still be mapped. The individual uptake and secretion profile for each sample is saved to the location specified in outputPath using the unique sample name.
- USAGE
prepIntegrationQuant (model, metData, exchanges, samples, test_max, test_min, outputPath, tol, variation)
- INPUTS
model – Prepared global model (e.g., model_for_CORE from prepModel)
exchanges – Vector containing exchange reactions
metData – Fluxes of uptake (negative) and secretion(positive) flux values. The columns are the samples and the rows are the metabolites (Unit matching remaining model constraints!).
samples – Vector of sample names (no dublicate names)
test_max – Maximal uptake/secretion set while testing if model can perform uptake and secretion of a metabolite (e.g., 500)
test_min – Minimal uptake/secretion set while testing if model can perform uptake and secretion of a metabolite (e.g., 0.00001)
outputPath – Path where output files should be saved (e.g. ‘Y:StudiesData_integrationCOREusingRecon1models')
- OPTIONAL INPUTS
tol – All fluxes below this value are considered to be zero, (e.g., 1e-6)
variation – Lower and upper bound are established with this value as error range in % (e.g.,20)
For every sample a file is automatically saved and contains the following variables:
Secretion_not_possible: Vector of metabolite (exchange reactions) that cannot be secreted by the model
Uptake_not_possible: Vector of metabolite (exchange reactions) that cannot be uptaken by the model
FBA_all_secreted: FBA results of test metabolite secretion
FBA_all_secreted_names: Name of exchange test metabolite secretion
FBA_all_uptake: FBA results of test metabolite uptake
FBA_all_uptake_names: Name of exchange test metabolite uptake
uptake: Vector of exchange reactions that are associated with uptake in the cell line (no additional exchanges, since these reactions will not be closed)
uptake_value: Matrix of flux values, constitute the lower (column 2) and upper (column 3) limits for the model uptake
secretion: Vector of exchange reactions that are associated with secretion in the cell line (no additional exchanges, since these reactions will not be closed)
secr_value: Matrix of flux values, constitute the lower (column 2) and upper (column 3) limits for the model secretion
- pruneModel(model, minGrowth, biomassRxn)[source]¶
This function prunes a model to its most compact subnetwork given some model constraints and a minimal growth rate by identifying and removing all blocked reactions.
- USAGE
[modelUpdated, modelPruned, Ex_Rxns] = pruneModel (model, minGrowth, biomassRxn)
- INPUTS
model – model structure
minGrowth – minimal Growth rate to be set on biomass reaction
biomassRxn – biomass reaction name (default: ‘biomass_reaction2’)
- OUTPUTS
modelUpdated – same as input model but constraints on blocked reactions are set to be 0
modelPruned – pruned model, where all blocked reactions are removed (attention this seems to cause issues with GPRs)
Ex_Rxns – List of exchange reactions in pruned model
- setConstraintsOnBiomassReaction(model, of, dT, tolerance)[source]¶
This function sets constrains biomass objective function in the model accoring to the condition specific doubling time.
- USAGE
[modelBM] = setConstraintsOnBiomassReaction (model, of, dT, tolerance)
- INPUTS
model – Metabolic model (e.g., Recon)
of – Objective funtion, e.g., biomass_reaction
dT – Doubling time
tolerance – Upper and lower limit of the growth rate are adjusted according to this tolerance value, e.g., 20 (%).
- OUTPUTS
modelBM – Model constrained with condition-specific growth rates
- setMediumConstraints(model, set_inf, current_inf, medium_composition, met_Conc_mM, cellConc, t, cellWeight, mediumCompounds, mediumCompounds_lb, customizedConstraints, customizedConstraints_ub, customizedConstraints_lb, close_exchanges)[source]¶
Calculates and sets constraints according to fresh medium composition in mM. Is based on the function Conc2Rate. Returns a model with constraints.
- USAGE
[modelMedium, basisMedium] = setMediumConstraints (model, set_inf, current_inf, medium_composition, met_Conc_mM, cellConc, t, cellWeight, mediumCompounds, mediumCompounds_lb, customizedConstraints, customizedConstraints_ub, customizedConstraints_lb, close_exchanges)
- INPUTS
model – Metabolic model (Recon)
current_inf – Models can have differently defined infinite constraints, e.g., 500
set_inf – New value for infinite constraints, e.g., 1000
medium_composition – Vector of exchange reactions of metabolites in the cell medium, e.g., RPMI medium
met_Conc_mM – Vector of the same length of ‘Exchanges’, providing the concentration (in mM) of each metabolite in the respective row
cellConc – Cell concentration (cells per 1 ml)
t – Time in hours
cellWeight – Cellular dry weight
mediumCompounds – Composition of basic medium, which are usually not defined in the composition of the medium and need to be defined in addition.
mediumCompounds_lb – Lower bound applied for all mediumCompounds (the same for all)
- OPTIONAL INPUTS
customizedConstraints – If additional constraints should be set apart from the basic medium and medium composition, e.g. EX_o2(e)
customizedConstraints_lb – Vecor of lower bounds to be set (mmol/gDW/hr), must be of same length as customizedConstraints
customizedConstraints_ub – Vecor of upper bounds to be set (mmol/gDW/hr), must be of same length as customizedConstraints
close_exchanges – If exchange reactions except those specified should be closed (Default = 1) 1= close exchanges, 0=no
- OUTPUTS
modelMedium – Model with set constraints, all exchanges not specified in mediumCompounds, ions or customizedConstraints are set to secretion only
basisMedium – Vector of adjusted constraints, for reference such that these constraints are not overwritten at a later stage
Please note that if metabolites of medium_composition and mediumCompounds overlap, the constraints of the Medium_composition will be set to 0 in the output model. The function depends on the functions conc2Rate, changeRxnBounds.
- setQualitativeConstraints(model, cond_uptake, cond_uptake_LODs, cond_secretion, cond_secretion_LODs, cellConc, t, cellWeight, ambiguous_metabolites, basisMedium)[source]¶
This function sets qualitative constraints (by enforcing minimal uptake or secretion based on individual detection limits), e.g., based on the uptake and secretion profile of metabolites measured through mass-spectrometry. Uptake is only possible if the lower bound has been set to a value >0 using setMediumConstraints The minimal allowable uptake and secretion flux is defined by enforcing uptake at or above the limit of detection (mass-spectrometry). If these values are not available, a very small value, e.g., 1.0E-06 can be used. Note that this value has to be below the concentrations defined in the medium. Otherwise the model will be infeasible. Alternatively to the LODs a small value can be used (e.g., 0.00001). Note that the value has to be above the threshold later called zero.
- USAGE
[modelLOD] = setQualitativeConstraints (model, cond_uptake, cond_uptake_LODs, cond_secretion, cond_secretion_LODs, cellConc, t, cellWeight, ambiguous_metabolites, basisMedium)
- INPUTS
model – Metabolic model (Recon), with set constraints (output of setMediumConstraints)
cond_uptake – Vector of exchanges of metabolites consumed by the cells in the experiment
cond_uptake_LODs – Vector of detection limits (LOD in Mm) for the compounds and in the experiment
cond_secretion – Vector of metabolite exchanges consumed by the cells in the experiment
cond_secretion_LODs – Vector of detection limits (LOD in Mm) for the compounds and in the experiment
cellConc – Cell concentration (cells per 1 ml)
t – Time in hours
cellWeight – gDW per cell
ambiguous_metabolites – Since all exchanges, except the ones specified in uptake and secretion are closed in modelLOD, this input variable allows to specify metabolite exchanges that should remain open. Thus, if these exchanges are open in the starting model they can be uptaken or secreted by the modelLOD, e.g., ambiguous_metabolites This can for example be the case if it is suspected that an uptake or secretion might have taken place in concentrations below the detection limit.
basisMedium – Vector defining the metabolite exchanges of the basic medium, i.e., ions, mediumCompounds
- OUTPUT
modelLOD – Model that is constrained qualitatively to the condition-specific uptake and secretion profile
- setQuantConstraints(model, samples, tol, minGrowth, obj, no_secretion, no_uptake, medium, addExtraExch, addExtraExch_value, path)[source]¶
This function takes a model and quantitative extracellular metabolomic data and returns a model in which the data is integrated as constraints described in “MetaboTools: A Comprehensive Toolbox for Analysis of Genome-Scale Metabolic Models”, 2016, Front. Physiol., and the supplemental data “Tutorial II: Workflow for the integration of quantitative extracellular metabolomic data into the network context.” (Data Sheet 3.pdf)
- USAGE
[ResultsAllCellLines, OverViewResults] = setQuantConstraints (model, samples, tol, minGrowth, obj, no_secretion, no_uptake, medium, addExtraExch, addExtraExch_value, path, epsilon)
- INPUTS
model – Global metabolic model (Recon)
samples – Vector specifying the samples used (there must be an output file of function …. for each sample)
tol – Cutoff value for small numbers (e.g., -1e-8). All number smaller than tol will be treated as zero
minGrowth – Will be the lower bound of the objective function (e.g., 0.008). Forces the output model(s) to be able to produce a minimal objective value
obj – Objective function, e.g. biomass_reaction2
no_secretion – Define metabolites that should not be secreted (e.g., {EX_o2(e)})
no_uptake – Define metabolites that should not be consumed (e.g., {EX_o2s(e), EX_h2o2(e)})
medium – Define if certain exchanges should be excluded from minimization of exchanges (e.g., {}, if no medium except the exometabolomic data has been defined)
addExtraExch – After adding secretions, models are still not growing, this variable allows one to recover exchanges with a defined small value
addExtraExch_value – e.g. 1 as arbitrary small flux value / the resulting ub = 1, lb = -1.
path – Location of the .mat files for samples.
- OUTPUTS
ResultsAllCellLines – Structure that contains pruned and unpruned model, Vector of the Exchange_reactions, the exchange reactions added by minExCard, minFLux and maxFlux of the added reactions, maximal objective value, and the results of the gene deletion
OverViewResults – Overview of model statistics, e.g., number of reactions, metabolites, genes, number of essential genes, min and max objective values for easy comparison between sets of models
Depends on changeRxnBounds, fluxVariability, optimizeCbModel, generateCompactExchModel as well its dependent functions pruneModel, findMinCardModel, findOptExchRxns
- setSemiQuantConstraints(modelA, modelB, cond1_upt_higher, cond2_upt_higher, cond2_secr_higher, cond1_secr_higher)[source]¶
This function establishes constraints for two models (modelA and modelB) based on the relative difference of signal intensities obtained from two samples. The prerequisite for the adjustment is that constraints have been applied on which the relative differences can be applied:
The uptake rates have been established based on the medium composition and concentrations (i.e., setMediumConstraints), and
qualitative constraints have been applied based on know detection limits or small values (i.e., setQualitativeConstraints).
- USAGE
[modelA_QUANT, modelB_QUANT] = setSemiQuantConstraints (modelA, modelB, cond1_upt_higher, cond2_upt_higher, cond2_secr_higher, cond1_secr_higher)
- INPUTS
modelA – model constrained according to condition 1
modelB – model constrained according to condition 2
cond1_upt_higher – Exchange reactions and relative differences with higher uptake in condition 1 compared to condition 2
cond2_upt_higher – Exchange reactions and relative differences with higher uptake in condition 2 compared to condition 1
cond2_secr_higher – Exchange reactions and relative differences with higher secretion in condition 2 compared to condition 1
cond1_secr_higher – Exchange reactions and relative differences with higher secretion in condition 1 compared to condition 2
- OUTPUTS
modelA_QUANT – model constrained according to condition 1
modelB_QUANT – model constrained according to condition 2
- statisticsAddedExchanges(ResultsAllCellLines, samples)[source]¶
This function performs statistics on the added exchange reactions using the structure ResultsAllCellLines generated by the function setQuantConstraints and across all samples.
- USAGE
[Ex_added_all_unique] = statisticsAddedExchanges (ResultsAllCellLines, samples)
- INPUTS
ResultsAllCellLines – Last output ResultsAllCellLines that contains results of all samples (last in sampleResult.mat) e.g. in ‘Y:StudiesData_integrationTutorial_paperToolboxUACC_257_2Result.mat’
samples
- OUTPUTS
Ex_added_all_unique – a table of exchange reactions which the models contain in addition to the mapped(metabolomics-based)-exchanges, ordered by the frequency of appearance.
Column 1 = exchange reactions,
Column 2 = frequency of use,
Column 3 = secretion,
Column 4 = uptake.
- summarizeSamplingResults(modelA, modelB, path, nFiles, pointsPerFile, starting_model, dataGenes, show_rxns, fonts, hist_per_page, bin, fileNameA, fileNameB)[source]¶
The function summarizes the sampling results of the two models modelA and modelB. Subsequently, it returnes medians from the two models’ sampling, FVA results for the respective reaction, simple test results (stats). Finally, histograms are produced, one for each reactions containing the distributions of the two models. The set of analyzed reactions can be limited (show_rxns). Reactions associated with genes of special interest, e.g. differentially expressed genes, can be marked to facilitate the analysis, and simplify the identification of interesting histograms.
- USAGE
[stats, statsR] = summarizeSamplingResults (modelA, modelB, path, nFiles, pointsPerFile, starting_model, dataGenes, show_rxns, fonts, hist_per_page, bin, fileNameA, fileNameB)
- INPUTS
modelA – Sampled modelA (condition 1)
modelB – Sampled modelB (condition 2)
path – Path to sampling output files
nFiles – Number of files saved, e.g., 20;
pointsPerFilePoints – Points saved per file, e.g., 5000;
starting_model – Original metabolic model (Recon)
dataGenes – Gene set, whose associated reactions should be emphasized by color, e.g., alternatively spliced or differentially expressed genes
- OPTIONAL INPUTS
show_rxns – Vector of reactions that should be displayed, e.g., certain pathways or non-loop reactions
fonts – Font size (default = 9)
hist_per_page – Defines the number of histogramms that are printed per page, e.g., 4, 9, or 25 (default = 9)
bin – Binning for histogramm (default = 30)
fileNameA – Name of files the sampling points are stored in
fileNameB – Name of files the sampling points are stored in
- OUTPUTS
stats – Statistics from the sampling (Columns: Median modelA, Median modelB, minFlux_A, maxFlux_A, minFlux_B, maxFlux_B )
statsR – Reaction vector for stats output
Depends on COBRA functions: loadsamples, findRxnsFromGenes, flux variability analysis / fastFVA