Metabotools

analyzeSingleGeneDeletion(ResultsAllCellLines, path, samples, cutoff, OverViewResults)

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)

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)

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)

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)

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)

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)

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)

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

  1. Control TP1,

  2. Control TP2,

  3. Condition 1 TP1,

  4. Condition 1 TP 2

input_B: Matrix, 4 columns

  1. Control2 TP1,

  2. Control2 TP2,

  3. Condition 2 TP1,

  4. 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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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)

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:

    1. The uptake rates have been established based on the medium composition and concentrations (i.e., setMediumConstraints), and

    1. 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)

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)

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