Pscmtoolbox

LeakTestRecon[source]

Leak test for a closed whole-body/Recon model.

Closes all exchange (Ex_, EX_), demand (DM_), and sink (sink_) reactions of the workspace model modelClosed, then runs fastLeakTest to identify metabolites that can be produced without any input (leaks).

changeCobraSolver(‘tomlab_cplex’,’lp’); modelClosed = modelConsistent;

Test4HumanFctExtv5(model, test, optionSinks)[source]

Test a model for the ~460 human metabolic functions

Tests a model for the ~460 human metabolic functions (duplicates removed). For each function an FBA problem is solved and the optimal objective value is recorded.

USAGE:

[TestSolution, TestSolutionName, TestedRxns, PercTestedRxns] = Test4HumanFctExtv4 (model, test, optionSinks)

INPUTS:

  • model – Model structure (Recon1, with the desired in silico condition), with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .c - objective coefficients

  • test – Which test set to run: Recon1, IECori, IEC, or all (default). Choose IECori to test the IEC model or a model that contains lumen (‘u’) as a compartment, otherwise choose IEC; all checks Recon1 and IEC

  • optionSinks – If true, set sink reactions to 0 (default false, leave unchanged). Note that all lower bounds of exchanges and demands will be set to 0

OUTPUTS:

  • TestSolution – Array containing the optimal objective value for the different tests

  • TestSolutionName – Array containing the names for the different tests

  • TestedRxns – Reactions of the original model that were active in at least one test

  • PercTestedRxns – Percentage of the original model reactions that were tested

adjustWholeBodyRxnCoeff(model, listOrgan, listCoeff)[source]

Adjust the stoichiometric coefficients of the whole-body maintenance reaction

This function adjusts the coefficients of the whole-body biomass maintenance (WBM) reaction. The WBM reaction contains each organ present in the whole-body metabolic reconstruction. For each organ, the stoichiometric coefficients represent the fractional weight contribution of the respective organ to the whole body weight. These coefficients can be updated to reflect individual-specific body contributions; e.g. in obese individuals the ratio of muscle and adipose tissue differs from that of a normal-BMI individual.

USAGE:

[model] = adjustWholeBodyRxnCoeff (model, listOrgan, listCoeff)

INPUTS:

  • model – Whole-body metabolic model with fields:

    • .rxns - reaction identifiers

    • .S - stoichiometric matrix

    • .A - constraint matrix (stoichiometry plus coupling constraints); created from .S if absent

    • .mets - metabolite identifiers

  • listOrgan – List of organs whose stoichiometric coefficient should be updated

  • listCoeff – List of coefficients that replace the current ones in the WBM reaction (order must match the order of organs in listOrgan)

OUTPUT:

model – Whole-body metabolic model with adjusted stoichiometric coefficients in the whole-body maintenance reaction

annotateHH[source]

annotate Harvey and Harvetta

annotateModel(model, annotateRxns, annotateMets, modelID, modelName, modelAnnotation)[source]

Annotate a model with VMH reaction and metabolite identifiers

Annotates a model with VMH reaction and metabolite identifiers by reading the VMH reaction and metabolite tables and populating the corresponding annotation fields.

USAGE:

model = annotateModel (model, annotateRxns, annotateMets, modelID, modelName, modelAnnotation)

INPUTS:

  • model – Model structure with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .modelID - model identifier

    • .modelName - model name

    • .modelAnnotation - model annotation string

    • .rxnECNumbers - EC numbers for reactions

    • .rxnMetaNetXID - MetaNetX reaction identifiers

    • .rxnKEGGID - KEGG reaction identifiers

    • .rxnSBOTerms - SBO terms for reactions

    • .rxnSEEDID - SEED reaction identifiers

    • .metNames - metabolite names

    • .metFormulas - metabolite elemental formulas

    • .metCharges - metabolite charges

    • .metChEBIID - ChEBI identifiers

    • .metHMDBID - HMDB identifiers

    • .metInChIString - InChI strings

    • .metKEGGID - KEGG metabolite identifiers

    • .metSmiles - SMILES strings

    • .metMetaNetXID - MetaNetX metabolite identifiers

    • .metPubChemID - PubChem identifiers

    • .metBiGGID - BiGG metabolite identifiers

    • .metBioCycID - BioCyc metabolite identifiers

    • .metSEEDID - SEED metabolite identifiers

    • .metSBOTerms - SBO terms for metabolites

    • .metChemSpider - ChemSpider identifiers

    • .metInchiKey - InChIKey identifiers

  • annotateRxns – If 1 (default), annotate reactions

  • annotateMets – If 1 (default), annotate metabolites

OPTIONAL INPUTS:

  • modelID – Identifier to assign to model.modelID

  • modelName – Name to assign to model.modelName

  • modelAnnotation – Annotation string to assign to model.modelAnnotation

OUTPUT:

model – Updated model structure with populated VMH annotation fields

annotateOrganAtlas[source]

annotate OrganAtlas

calcOrganFract(model, IndividualParameters)[source]

Extrapolate organ weight fractions from body-weight polynomials

This function extrapolates the organ weight fractions based on the polynomials given in http://www.ams.sunysb.edu/~hahn/psfile/pap_obesity.pdf (PMID 19267313), Table 3. Organs not given there are assumed to remain constant with weight, using the fractions from the reference man and reference woman.

USAGE:

[organs, OrganWeight, OrganWeightFract, IndividualParameters] = calcOrganFract (model, IndividualParameters)

INPUTS:

  • model – Model structure with fields:

    • .rxns - reaction identifiers

  • IndividualParameters – Structure of individual parameters with fields:

    • .sex - ‘male’ or ‘female’

    • .bodyWeight - individual body weight (kg)

    • .CardiacOutput - cardiac output (ml), used here as the blood volume

    • .BloodVolume - blood volume (ml)

    • .OrgansWeights - per-organ weights (working field)

    • .OrgansWeightsRefMan - per-organ weights of the reference man/woman

OUTPUTS:

  • organs – List of organs

  • OrganWeight – List of organ weights (same order as organs)

  • OrganWeightFract – List of organ weight fractions (same order as organs)

  • IndividualParameters – Updated structure of individual parameters

calculateBMR(sex, weight, height, age)[source]

Calculate the basal metabolic rate (BMR) from phenomenological equations

Calculates the basal metabolic rate using the phenomenological model proposed by Harris-Benedict and derivations thereof (see references). See also the corresponding Wikipedia entry: https://en.wikipedia.org/wiki/Harris%E2%80%93Benedict_equation

USAGE:

BMR = calculateBMR (sex, weight, height, age)

INPUTS:

  • sex – ‘male’ or ‘female’

  • weight – body weight in kg

  • height – body height in cm

  • age – age in years

OUTPUT:

BMR – 3 x 1 array of basal metabolic rate estimates, calculated from

  1. the original Harris-Benedict equations [1], [2]

  2. the Harris-Benedict equations revised by Roza and Shizgal in 1984 [3]

  3. the Harris-Benedict equations revised by Mifflin and St Jeor in 1990 [4]

References: [1] Harris JA, Benedict FG (1918). “A Biometric Study of Human Basal Metabolism”. Proceedings of the National Academy of Sciences of the United States of America. 4 (12): 370?3. doi:10.1073/pnas.4.12.370. PMC 1091498?Freely accessible. PMID 16576330. [2] A Biometric Study of Basal Metabolism in Man. J. Arthur Harris and Francis G. Benedict. Washington, DC: Carnegie Institution, 1919. [3] Roza AM, Shizgal HM (1984). “The Harris Benedict equation reevaluated: resting energy requirements and the body cell mass”. The American Journal of Clinical Nutrition. 40 (1): 168?82. PMID 6741850. [4] Mifflin MD, St Jeor ST, Hill LA, Scott BJ, Daugherty SA, Koh YO (1990). “A new predictive equation for resting energy expenditure in healthy individuals”. The American Journal of Clinical Nutrition. 51 (2): 241?7. PMID 2305711.

checkIEM_WBM(model, IEMRxns, BiomarkerRxns, minRxnsFluxHealthy, reverseDirObj, fractionKO, minBiomarker, fixIEMlb, LPSolver)[source]

Perform inborn-error-of-metabolism (IEM) simulations on a whole-body model

This function performs inborn error of metabolism simulations by deleting (or reducing) the flux through the reaction(s) affected by a gene defect and optimising the flux through a defined set of biomarker reactions.

USAGE:

[IEMSol] = checkIEM_WBM (model, IEMRxns, BiomarkerRxns, minRxnsFluxHealthy, reverseDirObj, fractionKO, minBiomarker, fixIEMlb, LPSolver)

INPUTS:

  • model – Whole-body metabolic reconstruction or Recon3D model, with fields:

    • .A - constraint matrix (stoichiometry plus coupling constraints)

    • .S - stoichiometric matrix

    • .C - additional coupling-constraint matrix

    • .rxns - reaction identifiers

    • .c - objective coefficients

    • .osense - objective sense (-1 maximise, +1 minimise)

    • .osenseStr - objective sense (‘max’ or ‘min’)

  • IEMRxns – Reaction(s) affected by the inborn error of metabolism

  • BiomarkerRxns – Biomarker reactions and their reported direction of change

  • minRxnsFluxHealthy – Minimum flux value(s) through the IEMRxns

  • reverseDirObj – The function maximises the objective flux by default; if set to 1, it also checks the minimisation problem

  • fractionKO – Fraction of knockout applied to the IEM reactions (default 1 for 100% knockout)

  • minBiomarker – Minimise through the biomarker reaction (default 0)

  • fixIEMlb – Fix the IEM reaction lb = ub (default 0, i.e. lb = 0 while ub = (1-fractionKO)*solution.v(find(model.c)))

  • LPSolver – LP solver to use (‘ILOGcomplex’ default, ‘tomlab_cplex’)

OUTPUTS:

IEMSol – Predicted biomarker fluxes and comparison with the reported biomarkers

EXAMPLE:

Exampe of preparation of a set of inputs to checkIEM_WBM

R = {‘_2OXOADPTm’;’_2AMADPTm’;’_r0879’}; RxnsAll2 = ‘’; for i = 1: length(R)

RxnsAll = model.rxns(find(~cellfun(@isempty,strfind(model.rxns,R{i})))); RxnsAll2 =[RxnsAll2;RxnsAll];

end IEMRxns = unique(RxnsAll2); RxnMic = model.rxns(find(~cellfun(@isempty,strfind(model.rxns,’Micro_’)))) ; IEMRxns = setdiff(IEMRxns,RxnMic);

if ~strcmp(modelName,’Recon3D’)

% add demand reactions to blood compartment for those biomarkers reported for blood % biomarker based on https://www.omim.org/entry/204750 model = addDemandReaction(model, ‘L2aadp[bc]’);

BiomarkerRxns = {

‘DM_L2aadp[bc]’ ‘Increased (blood)’ ‘EX_2oxoadp[u]’ ‘Increased (urine)’ ‘EX_adpoh[u]’ ‘Increased (urine)’ };

else
BiomarkerRxns = {

‘EX_2oxoadp[u]’ ‘Increased (urine)’ ‘EX_adpoh[u]’ ‘Increased (urine)’ };

end

Then call the checkIEM_WBM function

[IEMSol_2OAA] = checkIEM_WBM(model,IEMRxns, BiomarkerRxns,minRxnsFluxHealthy);

Example IEMSol returned from the above

{‘IEM Rxns All obj - Healthy’} {‘65403.5393’ } {0×0 double } {‘IEM Rxns All obj - Disease’} {‘65403.5393’ } {0×0 double } {‘WB obj - Healthy’ } {‘NA’ } {0×0 double } {‘WB obj - Disease’ } {‘1’ } {‘1’ } {‘Healthy:DM_L2aadp[bc]’ } {‘-5.5324e-08’} {‘Disease - Reported:Increas…’} {‘Disease:DM_L2aadp[bc]’ } {‘40.2418’ } {‘Disease - Reported:Increas…’} {‘Healthy:EX_2oxoadp[u]’ } {‘-2.3603e-10’} {‘Disease - Reported:Increas…’} {‘Disease:EX_2oxoadp[u]’ } {‘3.7368’ } {‘Disease - Reported:Increas…’} {‘Healthy:EX_adpoh[u]’ } {‘-1.6985e-10’} {‘Disease - Reported:Increas…’} {‘Disease:EX_adpoh[u]’ } {‘3.7368’ } {‘Disease - Reported:Increas…’}

compareBounds2Models(model1, model2)[source]

Compare the reaction bounds between two models

Finds the reactions shared by two models and compares their lower and upper bounds, reporting reactions missing from either model and reactions whose bounds match or mismatch.

USAGE:

[R1_missing, R2_missing, R12_bounds, R12_bounds_mismatch] = compareBounds2Models (model1, model2)

INPUTS:

  • model1 – COBRA model structure with fields:

    • .rxns - reaction identifiers

    • .lb - lower bounds

    • .ub - upper bounds

  • model2 – COBRA model structure with fields:

    • .rxns - reaction identifiers

    • .lb - lower bounds

    • .ub - upper bounds

OUTPUTS:

  • R1_missing – Reactions absent from model1 (present in model2)

  • R2_missing – Reactions absent from model2 (present in model1)

  • R12_bounds – Cell array of shared reactions whose lower/upper bound is identical in both models, with the bound values

  • R12_bounds_mismatch – Cell array of shared reactions whose lower/upper bound differs between the two models, with the bound values

compareMaleFemale(male, female)[source]

Compare basic features of the male and female whole-body metabolic models

This function compares basic features of the male and female whole-body metabolic models, reporting shared and gender-unique reactions, per-organ reaction counts, gallbladder subsystem enrichment, and biofluid exchanges.

USAGE:

[ResultsMaleFemale] = compareMaleFemale (male, female)

INPUTS:

  • male – Male whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .subSystems - subsystem annotations

  • female – Female whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .subSystems - subsystem annotations

OUTPUT:

ResultsMaleFemale – Structure of differences and commonalities between the two models, with fields:

  • .MaleOnly - reactions unique to the male model

  • .FemaleOnly - reactions unique to the female model

  • .BothGender - reactions present in both models

  • .OrgansNumRxnMale - per-organ reaction counts (male): total, male-only, fraction

  • .OrgansNumRxnFemale - per-organ reaction counts (female): total, female-only, fraction

  • .maleOrgans - list of organs in the male model

  • .femaleOrgans - list of organs in the female model

  • .MaleGallSSEnrich - subsystems enriched among male-only gallbladder reactions

  • .FemaleGallSSEnrich - subsystems enriched among female-only gallbladder reactions

  • .MaleOnlyBiofluid - male-only biofluid exchange (_EX_) reactions

  • .FemaleOnlyBiofluid - female-only biofluid exchange (_EX_) reactions

convertATPflux2StepNumer(ATP_hydrolysis_flux, sex, weight, height)[source]

Convert an ATP hydrolysis flux into distance walked and step number

This function converts an ATP hydrolysis flux (e.g. through Muscle_DM_atp_c_) into distance walked and step number. See the body for the assumptions and calculation details.

USAGE:

[Energy_kJ, Energy_kcal, Meter, StepNumber] = convertATPflux2StepNumber (ATP_hydrolysis_flux, sex, weight, height)

INPUTS:

  • ATP_hydrolysis_flux – Flux value through the Muscle_DM_atp_c_ reaction (mmol/person/day)

  • sex – ‘male’ or ‘female’

  • weight – body weight in kg

  • height – body height in cm

OUTPUTS:

  • Energy_kJ – Energy (kJ) corresponding to the ATP hydrolysis flux

  • Energy_kcal – Energy (kcal) corresponding to the ATP hydrolysis flux

  • Meter – Corresponding distance walked (m) that can be achieved

  • StepNumber – Corresponding number of steps that can be achieved

createModelNewCompartment(model, OldComp, NewComp, NewCompName, LB, UB, RemoveExch)[source]

Convert a two-compartment metabolic model into a three-compartment model

This function converts a two-compartment metabolic model into a three-compartment metabolic model by replicating an existing compartment under a new name.

USAGE:

[modelComp] = createModelNewCompartment (model, OldComp, NewComp, NewCompName, LB, UB, RemoveExch)

INPUTS:

  • model – Model structure with fields:

    • .rxns - reaction identifiers

    • .rxnNames - reaction names

  • OldComp – Name of the current (to be replicated) compartment

  • NewComp – Abbreviation for the new compartment, e.g. ‘lu’

  • NewCompName – Name of the new compartment, e.g. ‘lumen’

  • LB – Lower bound on the new exchange compartment (default -1000)

  • UB – Upper bound on the new exchange compartment (default 1000)

  • RemoveExch – If set to 1, remove the old compartment (default 0)

OUTPUT:

modelComp – Model with the added (third) compartment

determineFluxValuesOnBoundary(model, solution)[source]

Determine the fraction of reactions whose flux values lie on the bounds

This function determines the fraction of reactions in the flux distribution whose flux values are on the lower or upper bounds.

USAGE:

[OfConstraint, OfAll] = determineFluxValuesOnBoundary (model, solution)

INPUTS:

  • model – Model structure with fields:

    • .lb - lower bounds

    • .ub - upper bounds

    • .rxns - reaction identifiers

  • solution – Solution structure with fields:

    • .full - full flux vector (primal solution)

OUTPUTS:

  • OfConstraint – Fraction of flux values on the lower/upper bounds of all constrained reactions (assuming +/- 1,000,000 as infinity)

  • OfAll – Fraction of flux values on the lower/upper bounds of all reactions in the model

findMetinOrgan(WBModel, metabolite)[source]

Find all organs in which a metabolite participates

USAGE:

[OrganListLong, OrganListOnly] = findMetinOrgan (WBModel, metabolite)

INPUTS:

  • WBModel – Whole-body metabolic model, with fields:

    • .mets - metabolite identifiers

  • metabolite – Abbreviation of the metabolite to be looked up

OUTPUTS:

  • OrganListLong – Metabolite identifiers (with organ prefix) in which the metabolite participates

  • OrganListOnly – Unique list of organs in which the metabolite occurs

getBasicHarveyStats(male, female)[source]

Compile basic statistics on the male and female whole-body metabolic models

USAGE:

[TableHHStats] = getBasicHarveyStats (male, female)

INPUTS:

  • male – Male whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .genes - gene identifiers

    • .subSystems - subsystem annotations

  • female – Female whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .genes - gene identifiers

    • .subSystems - subsystem annotations

OUTPUT:

TableHHStats – Cell array (table) of basic statistics for the two models

getListOfUniqueIsozymes(model)[source]

Get a preliminary list of unique isozymes for a given model

This function gets a preliminary list of unique isozymes for a given model or reconstruction by inspecting the gene-protein-reaction rules.

USAGE:

listIsozymes = getListOfUniqueIsozymes (model)

INPUT:

model – Model or reconstruction structure, with fields:

  • .genes - gene identifiers

OUTPUT:

listIsozymes – Preliminary list of unique isozymes for the model

getOrgansFromHarvey(modelWBM, runTestsOnly, OrganCompendium, printLevel)[source]

Extract the individual organs from a whole-body metabolic model

This function cuts the organs from the whole-body metabolic model. The different biofluid compartments are retained, but all constraints on the exchange and transport reactions are overwritten. Once the organs are extracted, the function runs a sanity check on each organ; this step can also be done independently by setting runTestsOnly to 1, in which case the OrganCompendium must be provided as input. The function also loads Recon3* so that the test results can be compared with the organ test results. The organ compendium for each sex is saved as OrganAtlas_Harvetta.mat and OrganAtlas_Harvey.mat along with the test results.

USAGE:

[OrganCompendium, TableCSources] = getOrgansFromHarvey (modelWBM, runTestsOnly, OrganCompendium, printLevel)

INPUTS:

  • modelWBM – Whole-body metabolic model, with fields:

    • .sex - ‘male’ or ‘female’

    • .S - stoichiometric matrix

    • .A - constraint matrix (stoichiometry plus coupling constraints)

    • .mets - metabolite identifiers

  • runTestsOnly – If 1, only run the organ tests (requires OrganCompendium as input); default 0

  • OrganCompendium – Previously extracted organ compendium (used when runTestsOnly is 1), with fields:

    • .sex - ‘male’ or ‘female’

    • .Recon3DHarvey - Recon3D reference model used to compare organ tests

  • printLevel – Verbosity level (default 0)

OUTPUTS:

  • OrganCompendium – Structure containing the individual organs and the basic tests that each organ passed

  • TableCSources – Overview table of ATP yield per carbon source under aerobic and anaerobic conditions for each organ in the model

getRxnsFromGene(model, gene, causal)[source]

Get all reaction(s) associated with a particular gene

This function gets all reaction(s) associated with a particular gene by screening through the grRules provided in the model structure.

USAGE:

[Rxns, grRules] = getRxnsFromGene (model, gene, causal)

INPUTS:

  • model – Model structure with fields:

    • .grRules - readable gene-protein-reaction rules

  • gene – Gene of interest

  • causal – If 1, get only genes whose loss would cause loss of function of the associated reactions; otherwise get all associated reactions (default)

OUTPUTS:

  • Rxns – List of reaction(s) associated with the input gene

  • grRules – List of grRules containing the input gene (same order as Rxns)

getStatsOrganComp(female, male, OrganCompendium_female, OrganCompendium_male, violinPlots)[source]

Compile general statistics on the male and female organ compendia

This function compiles general statistics on the male and female organ compendia derived from the male and female whole-body metabolic models. A more comprehensive comparison output is provided in the file Results_StatsOrganComp.mat, which is created at the end of this function.

USAGE:

[TableProp_female, TableProp_male, TableGRM, TableMetsNum_female, TableMetsNum_male, TableGenes_femaleNum] = getStatsOrganComp (female, male, OrganCompendium_female, OrganCompendium_male, violinPlots)

INPUTS:

  • female – Female whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .genes - gene identifiers

    • .subSystems - subsystem annotations

  • male – Male whole-body metabolic model, with fields:

    • .rxns - reaction identifiers

    • .mets - metabolite identifiers

    • .genes - gene identifiers

    • .subSystems - subsystem annotations

  • OrganCompendium_female – Structure of the female organs (from getOrgansFromHarvey)

  • OrganCompendium_male – Structure of the male organs (from getOrgansFromHarvey)

  • violinPlots – Plot violin plots (does not work below MATLAB 2016; default 0)

OUTPUTS:

  • TableProp_female – Table of organ-specific property information (female)

  • TableProp_male – Table of organ-specific property information (male)

  • TableGRM – Table comparing gene, reaction, and metabolite counts across organs

  • TableMetsNum_female – Metabolite-by-organ presence matrix (female)

  • TableMetsNum_male – Metabolite-by-organ presence matrix (male)

  • TableGenes_femaleNum – Gene-by-organ presence matrix (female)

linearRegression(x, y)[source]

Calculate the linear regression y = a0 + a1*x

USAGE:

[a0, a1, Rsqr, Residuals] = linearRegression (x, y)

INPUTS:

  • x – Values for the explanatory variable x in y = a0 + a1*x

  • y – Values for the dependent variable y in y = a0 + a1*x

OUTPUTS:

  • a0 – Intercept

  • a1 – Slope of the line

  • Rsqr – Square of the correlation coefficient

  • Residuals – Regression residuals (an objective measure of the goodness of fit of the linear regression equation)

optimizeWBModel(model, param)[source]

Optimise whole body metabolic model

Solves LP problems of the form

\[\begin{split}max/min ~& c^T v \\ s.t. ~& S v = dxdt ~~~~~~~~~~~:y \\ ~& C v \leq d~~~~~~~~:y \\ ~& lb \leq v \leq ub~~~~:w\end{split}\]
USAGE:

solution = optimizeWBModel (model, param)

INPUT:

model – Whole-body model structure with fields:

  • .S - m x n stoichiometric matrix

  • .c - n x 1 linear objective coefficients

  • .lb - n x 1 lower bounds

  • .ub - n x 1 upper bounds

  • .C - k x n left-hand side of the coupling constraints C*v <= d

  • .E - extra-variable constraint matrix

  • .osense - objective sense (-1 maximise, +1 minimise)

  • .osenseStr - objective sense (‘max’ or ‘min’)

OPTIONAL INPUT:

param – Additional parameters as a struct. Fields that are not COBRA solver parameters (see getCobraSolverParamsOptionsForType) for this problem type are passed on to the solver in a solver-specific manner. Recognised fields include:

  • .minNorm - if nonzero, minimise the Euclidean norm of the solution to the LP problem (0 by default; ~1e-6 is usually enough for regularisation)

  • .secondsTimeLimit - solver time limit in seconds

  • .verify - if true, verify the returned solution

  • .printLevel - verbosity level (-1 silent, 0 warnings only (default), 1-4 increasing solver detail)

  • .solveWBMmethod - method used to solve the WBM problem

  • .ScaleFlag - solver scaling flag

  • .scaind - solver scaling index

  • .emphasis_numerical - numerical emphasis flag

  • .MSK_DPAR_OPTIMIZER_MAX_TIME - MOSEK optimizer time limit

  • .MSK_IPAR_WRITE_DATA_PARAM - MOSEK data-writing parameter

  • .MSK_IPAR_LOG_INTPNT - MOSEK interior-point log level

  • .MSK_IPAR_LOG_PRESOLVE - MOSEK presolve log level

  • .MSK_IPAR_INTPNT_SCALING - MOSEK interior-point scaling

  • .MSK_IPAR_SIM_SCALING - MOSEK simplex scaling

OUTPUT:

solution – solution object:

  • f - Objective value

  • v - Reaction rates (Optimal primal variable, legacy FBAsolution.x)

  • y - Dual for the molecular species

  • w - Reduced costs of the reactions

  • s - Slacks of the molecular species

  • stat - Solver status in standardized form: * -1 - No solution reported (timelimit, numerical problem etc) * 0 - Infeasible * 1 - Optimal solution * 2 - Unbounded solution

  • origStat - Original status returned by the specific solver

  • ctrs_y - the duals for the constraints from C

  • ctrs_slack - Slacks of the additional constraints

organEssentiality(model, LPSolver)[source]

Compute organ essentiality in a whole-body model

This function computes the organ essentiality in a whole-body model by setting each organ’s value to zero in the whole-body objective reaction, setting all organ-specific reaction bounds to zero (lower and upper), and then computing whether a non-zero flux through this objective is still possible.

USAGE:

[ResultsOrganEss] = organEssentiality (model, LPSolver)

INPUTS:

  • model – Whole-body metabolic model, with fields:

    • .sex - ‘male’ or ‘female’

    • .A - constraint matrix (stoichiometry plus coupling constraints)

  • LPSolver – LP solver to use (‘ILOGcomplex’, ‘tomlab_cplex’ default)

OUTPUT:

ResultsOrganEss – Maximal flux value for the whole-body reaction for each organ: col 1 organ name, col 2 max flux, col 3 min flux (if minimisation is activated; default inactive), col 4 solver status (1 feasible, 5 feasible with numerical difficulties, 3 infeasible)

performIEMAnalysis(model, geneMarkerList, compartment, urine, minRxnsFluxHealthy, causal, reverseDirObj, fractionKO, minBiomarker, fixIEMlb, LPSolver)[source]

Perform inborn-error-of-metabolism (IEM) analysis from a list of genes

This function performs the IEM analysis from a list of genes, testing for the defined biomarker metabolites in one or more biofluid compartments.

USAGE:

[IEMSolutions, IEMTable, missingMetAll] = performIEMAnalysis (model, geneMarkerList, compartment, urine, minRxnsFluxHealthy, causal, reverseDirObj, fractionKO, minBiomarker, fixIEMlb, LPSolver)

INPUTS:

  • model – Whole-body metabolic model, with fields:

    • .modelID - model identifier

  • geneMarkerList – Cell array of gene marker lists and the biomarkers to test for

  • compartment – List of biofluid compartments in which the biomarkers appear

  • urine – Also test for urine excretion of the biomarker metabolite (default true)

  • minRxnsFluxHealthy – Minimum flux value(s) through the IEMRxns (default 0.75)

  • causal – If 1, use only genes whose loss causes loss of function of the associated reactions; otherwise use all associated reactions (default)

  • reverseDirObj – The function maximises the objective flux by default; if set to 1, it also checks the minimisation problem

  • fractionKO – Fraction of knockout applied to the IEM reactions (default 1)

  • minBiomarker – Minimise through the biomarker reaction (default 0)

  • fixIEMlb – Fix the IEM reaction lb = ub (default 0)

  • LPSolver – LP solver to use (‘ILOGcomplex’ default, ‘tomlab_cplex’)

OUTPUTS:

  • IEMSolutions – Structure of predictions for each gene (‘NA’ where a metabolite does not occur in a biofluid)

  • IEMTable – Cell array of predictions for each gene (same content as IEMSolutions)

  • missingMetAll – Metabolites not appearing in a biofluid

performSanityChecksonRecon(model, resultsFileName, ExtraCellCompIn, ExtraCellCompOut, runSingleGeneDeletion, resultsPath, param)[source]

Perform quality-control and quality-assurance checks on a reconstruction

This function performs various quality control and quality assurance tests on a metabolic reconstruction.

USAGE:

[TableChecks, Table_csources, CSourcesTestedRxns, TestSolutionNameOpenSinks, TestSolutionNameClosedSinks, TableSingleDel] = performSanityChecksonRecon (model, resultsFileName, ExtraCellCompIn, ExtraCellCompOut, runSingleGeneDeletion, resultsPath, param)

INPUTS:

  • model – Model structure, with fields:

    • .rxns - reaction identifiers

    • .rxnGeneMat - reaction-gene incidence matrix

    • .rev - reversibility indicator per reaction

    • .lb - lower bounds

    • .genes - gene identifiers

  • resultsFileName – File name of the generated output file

  • ExtraCellCompIn – Extracellular uptake compartment ([e] by default)

  • ExtraCellCompOut – Extracellular secretion compartment ([e] by default)

  • runSingleGeneDeletion – If 1, run single gene deletion (default 0)

  • resultsPath – Path where the results should be saved (default: location of the ‘MethodSection3.mlx’ folder)

  • param – Testing parameters, with fields:

    • .checkDuplicates - if 1, check for duplicated reactions using the ‘FR’ method (default 0)

    • .testRev - if 1, test whether .rev is consistent with .lb and .ub (default 0)

    • .epsilon - numerical tolerance for the duplicate-reaction check (default feasTol*100)

    • .modeFlag - mode flag for the duplicate-reaction check (default 0)

OUTPUTS:

  • TableChecks – Table overview of the performed tests and their outcomes

  • Table_csources – Table of ATP-yield test results from various carbon sources under aerobic and anaerobic conditions

  • CSourcesTestedRxns – List of reactions active when testing ATP yield from the various carbon sources

  • TestSolutionNameOpenSinks – Results when testing for 460 metabolic functions with all sinks open

  • TestSolutionNameClosedSinks – Results when testing for 460 metabolic functions with all sinks closed

  • TableSingleDel – Single gene deletion results (col 1 gene, col 2 growth ratio)

runIEM_HH[source]

Predict biomarker metabolites in whole-body models for inborn errors of metabolism (IEMs).

Loops over 57 IEMs and predicts known biomarker metabolites in the urine, blood, and CSF compartments of a whole-body model. Select the model via the sex variable (‘male’, ‘female’, or ‘Recon3D’). Metabolite and IEM abbreviations are defined at www.vmh.life.