Metabolomics

constrainRxns(model, specificData, param, mode, printLevel)[source]

Function used to apply constraints in the XomicsToModel function

USAGE:

[newModel, rxnsConstrained, rxnBoundsCorrected, newSpecificData] = constrainRxns (model, specificData, param, mode, printLevel)

INPUTS:

  • model – A COBRA model structure with fields:

    • .S - m x n stoichiometric matrix

    • .rxns - n x 1 cell array of reaction identifiers

    • .rxnNames - n x 1 cell array of reaction names

    • .lb - n x 1 vector with lower bounds

    • .ub - n x 1 vector with upper bounds

    • .SIntRxnBool - n x 1 boolean, true for internal reactions

    • .SinkRxnBool - n x 1 boolean, true for sink reactions

    • .DMRxnBool - n x 1 boolean, true for demand reactions

    • .c - n x 1 linear objective coefficient vector

    • .C - k x n coupling constraint matrix on reactions

    • .D - k x nEvars coupling constraint matrix on extra variables

    • .E - m x nEvars matrix coupling extra variables into the steady state

    • .evars - identifiers of the extra variables

    • .evarc - objective coefficients of the extra variables

    • .evarlb - lower bounds of the extra variables

    • .evarub - upper bounds of the extra variables

  • specificData – A structure containing context-specific data with fields:

    • .rxns2constrain - table of custom reaction constraints

    • .mediaData - table of fresh-media constraints

    • .essentialAA - exchange reactions for essential amino acids

    • .exoMet - table of exometabolomic fluxes, with fields:

      • .rxns - k x 1 cell array of reaction identifiers

      • .mean - k x 1 numeric array of mean measured reaction flux

      • .SD - k x 1 numeric array of standard deviation in measured reaction flux

      • .varID - identifiers of the measured extra variables

  • param – A structure containing the parameters for the function:

    • .TolMinBoundary - the reaction boundary’s minimum value

    • .TolMaxBoundary - the reaction boundary’s maximum value

    • .boundPrecisionLimit - precision limit of a flux estimate

    • .metabolomicWeights - String indicating the type of weights to be applied to penalise the difference between predicted and experimentally measured fluxes, where ‘SD’ weights = 1/(1+exoMet.SD^2), ‘mean’ weights = 1/(1+exoMet.mean^2), ‘RSD’ weights = 1/((exoMet.SD./exoMet.mean)^2)

    • .eta - lower bound on significant bound constraint perturbation (Default: feasTol*10)

    • .printLevel - verbose level controlling printed output

  • mode – String indicating the type of constraints to be applied: ‘customConstraints’ uses specificData.rxns2constrain, ‘mediaDataConstraints’ uses specificData.mediaData, ‘exometabolomicConstraints’ uses specificData.exoMet

  • printLevel – verbose level controlling printed output

OPTIONAL INPUTS:

  • param.weightLower_flx_default – scalar or n x 1 double User defined penalty applied to lower bounds of all fluxes. It overwrites the default penalty. It is overwritten for measured fluxes only by specificData.exoMet.penaltyLowerBoundPerturbation when available

  • param.weightUpper_flx_default – scalar or n x 1 double User defined penalty applied to upper bounds of all fluxes. It overwrites the default penalty. It is overwritten for measured fluxes only by specificData.exoMet.penaltyUpperBoundPerturbation when available

  • param.weightLower_e_default – scalar or number of extra variables x 1 double. User defined penalty applied to lower bounds of extra variables. It overwrites the default penalty. It is overwritten for measured extra variables only by specificData.exoMet.penaltyLowerBoundPerturbation when available

  • param.weightUpper_e_default – scalar or number of extra variables x 1 double. User defined penalty applied to upper bounds of extra variables. It overwrites the default penalty. It is overwritten for measured extra variables only by specificData.exoMet.penaltyUpperBoundPerturbation when available

  • param.weightExpFct – scalar. Factor by which to multiply experimental weights (weightExp_flx ,weightExp_e) in mode ‘allConstraints’

OUTPUTS:

  • newModel – A COBRA model with constraints applied

  • rxnsConstrained – List of the reactions constrained

  • rxnBoundsCorrected – list of reactions whose restrictions should have been modified to comply with TolMinBoundary, TolMaxBoundary and boundPrecisionLimit

  • newSpecificData – A new structure containing arguments for the XomicsToModel function

March 2026: extended to deal with extra variables - Tania

exoMetLOOCV(model, measuredFluxes, param)[source]

Leave-one-out cross-validation of exometabolomic flux fitting

USAGE:

[V, compareFluxes] = exoMetLOOCV (model, measuredFluxes, param)

INPUTS:

  • model – COBRA model with fields:

    • .S - m x n stoichiometric matrix

    • .rxns - n x 1 cell array of reaction identifiers

    • .osenseStr - objective sense passed to entropic FBA (‘min’ or ‘max’)

    • .cf - forward reaction linear weight for entropic FBA

    • .cr - reverse reaction linear weight for entropic FBA

    • .g - flux entropy weight for entropic FBA

    • .u0 - reference chemical potential for entropic FBA

    • .f - concentration entropy weight for entropic FBA

  • measuredFluxes – table with the fluxes obtained from exometabolomics experiments, with variables:

    • .rxns - k x 1 cell array of reaction identifiers

    • .mean - k x 1 double of measured mean flux

    • .SD - k x 1 double standard deviation of the measured flux

OPTIONAL INPUT:

param – a structure containing the parameters for the function:

  • .alpha - alpha in (alpha/2)(v-h)’*H*(v-h) (Default: 10000)

  • .approach - flux-fitting approach passed to the solver (Default: ‘QEFBA’)

  • .relaxBounds - True to relax bounds on reactions whose fluxes are fitted to exometabolomic data

  • .metabolomicWeights - String indicating the type of weights to be applied to penalise the difference between predicted and experimentally measured fluxes, where ‘SD’ weights = 1/(1+exoMet.SD^2), ‘mean’ weights = 1/(1+exoMet.mean^2), ‘RSD’ weights = 1/((exoMet.SD./exoMet.mean)^2)

OUTPUTS:

  • V – table with the fluxes predicted by leave one out cross validation, with variables:

    • .rxns - n x 1 cell array of reaction identifiers, one for each in input model

    • .v - n x 1 double of predicted mean flux

    • .lb - n x 1 double lower bound on reaction flux

    • .ub - n x 1 double upper bound on reaction flux

  • compareFluxes – table mapping the predicted fluxes V onto the measured fluxes for comparison

fitExperimentalFlux(model, vExp, weightLower, weightUpper, weightExpFlux, param)[source]

fit a vector of change in concentration, over a time interval, to the range of a stoichiometric matrix, minimising the weighed norm of the perturbations to the given lower and upper bounds on the corresponding exchange reactions of a model, to ensure that a steady state solution exists.

USAGE:

[v, p, q, dv, obj] = fitExperimentalFlux (model, vExp, weightLower, weightUpper, weightExpFlux, param)

INPUTS:

  • model –

    COBRA model with the following required fields (others can be supplied):

    • .S - m x n stoichiometric matrix

    • .mets - m x 1 metabolite identifiers

    • .b - m x 1 accumulation (change in concentration with time)

    • .c - n x 1 linear objective coefficients

    • .osense - objective sense

    • .lb - n x 1 lower bounds

    • .ub - n x 1 upper bounds

    vExp: n x 1 experimental flux vector (NaN if no info)

  • weightLower – n x 1 positive weight penalty on relaxation of lower bounds

  • weightUpper – n x 1 positive weight penalty on relaxation of upper bounds

  • weightExpFlux – n x 1 positive weight penalty on deviation from experimental flux

growthMediaToModel(model, specificData, param, coreRxnAbbr, modelGenerationReport)[source]

Apply growth-media constraints from context-specific media data to a COBRA model

USAGE:

[model, specificData, coreRxnAbbr, modelGenerationReport] = growthMediaToModel (model, specificData, param, coreRxnAbbr, modelGenerationReport)

INPUTS:

  • model – COBRA model with fields:

    • .rxns - n x 1 reaction identifiers

    • .lb - n x 1 lower flux bounds

    • .ub - n x 1 upper flux bounds

  • specificData – A structure containing context-specific data with fields:

    • .mediaData - table of fresh-media constraints (with .mets and .rxns columns)

    • .exoMet - table of exometabolomic fluxes

  • param – A structure containing the parameters for the function:

    • .metabolomicsBeforeExtraction - Logical, whether metabolomics is added before extraction

    • .debug - Logical, whether to save progress for debugging

    • .workingDirectory - directory where debug files are written

    • .printLevel - verbose level controlling printed output

    • .TolMinBoundary - the reaction boundary’s minimum value

    • .TolMaxBoundary - the reaction boundary’s maximum value

    • .relaxOptions - options passed to relaxedFBA if the problem becomes infeasible

  • coreRxnAbbr – core reaction identifiers

  • modelGenerationReport – report struct accumulated during model generation

OUTPUTS:

  • model – the input COBRA model with media constraints applied (fields .rxns, .rxnNames, .lb, .ub updated)

  • specificData – the input specificData with its media and exometabolomic fields updated

  • coreRxnAbbr – updated core reaction identifiers

  • modelGenerationReport – updated model-generation report

mediaConstraints(model, uptakeRates, uptakeInChi, uptakeNames)[source]

Find the highest uptake rate for each metabolite in the uptake data and set it as a lower bound on the exchange reaction for that metabolite. If the maximum uptake rate is higher than 0 (the metabolite is only secreted by the cells) then the lower bound is set to 0

USAGE:

[model, metsUnmapped, rxnsUpdated] = mediaConstraints (model, uptakeRates, uptakeInChi, uptakeNames)

INPUTS:

  • model – COBRA model

  • uptakeRates – numeric uptake/secretion rates (mmol/gDW/h) for each measured metabolite

  • uptakeInChi – InChI codes identifying the measured metabolites

  • uptakeNames – chemical names of the measured metabolites

OUTPUTS:

  • model – COBRA model with constrained uptake rates

  • metsUnmapped – metabolites that could not be mapped to the model

  • rxnsUpdated – exchange reactions whose bounds were updated

metabolomicsTomodel(model, specificData, param, coreRxnAbbr, modelGenerationReport)[source]

Integrates metabolomics data to COBRA models obtained in the cell culture media either by metabolomics experiments or by the content of the culture medium

USAGE:

[model, specificData, coreRxnAbbr, modifiedFluxes, modelGenerationReport] = metabolomicsTomodel (model, specificData, param, coreRxnAbbr, modelGenerationReport)

INPUT:

  • model – A generic COBRA model with fields:

    • .S - m x n stoichiometric matrix

    • .mets - m x 1 metabolite identifiers

    • .rxns - n x 1 reaction identifiers

    • .lb - n x 1 lower flux bounds

    • .ub - n x 1 upper flux bounds

    • .SIntRxnBool - n x 1 boolean, true for internal reactions

    • .lb_preconstrainRxns - n x 1 lower bounds saved before constraining reactions

    • .ub_preconstrainRxns - n x 1 upper bounds saved before constraining reactions

  • specificData – A structure containing the context-specific data, with fields:

    • .exoMet - table of exometabolomic fluxes

    • .rxns2constrain - table of custom reaction constraints

    • .essentialAA - exchange reactions for essential amino acids

    • .mediaData - table of fresh-media constraints

  • param – a structure containing the parameters for the function:

    • .printLevel - verbose level controlling printed output

    • .curationOverOmics - Logical, use curated data with priority over omics data

    • .addSinksexoMet - Logical, add sink reactions for metabolites in exometabolomic data

    • .metabolomicsBeforeExtraction - Logical, whether metabolomics is added before extraction

    • .debug - Logical, whether to save progress for debugging

    • .workingDirectory - directory where debug files are written

    • .TolMaxBoundary - the reaction boundary’s maximum value

    • .excludeMetsFromRelax - metabolites excluded from bound relaxation

    • .boundsToRelaxExoMet - which bounds may be relaxed to fit exometabolomic data

    • .relaxOptions - options passed to relaxedFBA if the problem becomes infeasible

  • coreRxnAbbr – Set of core reactions

  • modelGenerationReport – A struct array where the data will be saved, with fields:

    • .sinksAddedFromexoMet - reactions for which sinks were added from exometabolomic data

OUTPUTS:

  • model – A Context-specific COBRA model with the metabolomic data.

  • specificData – The exometabolomic data is updated according to the relaxations

  • coreRxnAbbr – A Context-specific COBRA model with the metabolomic data.

  • modifiedFluxes – New set of core reactions.

  • modelGenerationReport – an updated version of the struct array

metsToKeep(dataOrg, metRep, calRelErr, param)[source]

removes the metabolomics measurements from the dataOrg table that have low measurement quality based on the relative SD of technical replicate sample (RSDqc provided in the metRep variable) and/or the average relative error (measured vs actual concentration) of the calibration line samples (RE provided in the calRelErr variable)

USAGE:

[cleanedData, metToKeepSummary] = metsToKeep (dataOrg, metRep, calRelErr, param)

INPUTS:

  • dataOrg – A table with original metabolomics data in a long format with information about measured samples, compounds, retention times (RT), area(s), concentrations, etc. the following column is required for the analysis:

    • .compound - compound (metabolite names) identical

      to the names used in the metRep and calRE variables

  • metRep – A table with quality information based on the dataOrg in a long format (output from mzQuality tool). The following columns are required for the analysis:

    • .compound - compound (metabolite names) identical

      to the names used in the dataOrg and calRE variables

    • .RSDqc_* - one or more columns specifying the

      relative standard deviation of the repeated measurement of a (quality) sample (per batch) per metabolite

  • calRelErr – A table with quality information about the relative error (RE) of the concentration estimation of the calibration line samples based on the dataOrg in a long format with information about measured samples, compounds, retention times (RT), known concentrations, estimated concentrations, and relative error etc. The following columns are required for the analysis:

    • .compound - compound (metabolite names) identical

      to the names used in the dataOrg and calRE variables

    • .RE* - a columns specifying the relative error of

      the concentration estimation of the calibration line samples

  • param – a structure with fields:

    • .tresholdRSDqc - the treshold value for the relative SD of the repeated

      measurement of sample (default = 25 (%); based on the based on the “Guidelines and considerations for the use of system suitability and quality control samples in mass spectrometry assays applied in untargeted clinical metabolomic studies”)

    • .tresholdCalRE - the treshold value for the relative error of the

      concentration estimation of the calibration line samples (default = 25 (%)

OUTPUTS:

  • cleanedData – table in the same format as dataOrg without the metabolites of low quality (above set tresholds)

  • metToKeepSummary –

    table in the same format as metRep variable with an

    added columns:

    • .keep - specifies whether a metabolite is to be

      kept (1) or removed (0) from the further analysis

    • .REcheck - specifies whether a metabolite passed

      (1) or failed (0) the check based on the relative error of the concentration estimation of the calibration line sample

    • .RSDqcCheck - specifies whether a metabolite passed

      (1) or failed (0) the check based on the relative standard deviation of the repeated measurement of a (QC) sample

    • .sumRE - shows the average relative error of the

      concentration estimation of the calibration line sample per metabolite

EXAMPLE:

NOTE:

Author(s): Agnieszka Wegrzyn (2021)

plotExperimentalvsPredictedExchange(measuredFluxes, predictedFluxes, otherFluxes, extraFluxes, param)[source]

Plot a comparison of experimental and predicted exchange reaction rates

USAGE:

plotExperimentalvsPredictedExchange (measuredFluxes, predictedFluxes, otherFluxes, extraFluxes, param)

INPUTS:

  • measuredFluxes – table with the fluxes obtained from exometabolomics experiments, with variables:

    • .rxns - k x 1 cell array of reaction identifiers

    • .mean - k x 1 double of measured mean flux

    • .SD - k x 1 double standard deviation of the measured flux

    • .metNames - k x 1 metabolite names, used as labels when present

    • .labels - k x 1 cell array of labels to display on y axis

    • .Properties - table metadata; .Properties.Description is used in the plot legend

    • .zerobool - computed boolean, measured flux is statistically zero (mean +/- SD spans 0)

    • .posbool - computed boolean, measured flux is positive (secretion)

    • .negbool - computed boolean, measured flux is negative (uptake)

    • .uptakeProbability - computed probability that the reaction is an uptake

    • .toRankOrder - computed key used to rank reactions on the plot

  • predictedFluxes – table with the predicted fluxes, with variables:

    • .rxns - k x 1 cell array of reaction identifiers

    • .v - k x 1 double of predicted mean flux

    • .lb - k x 1 double lower bound on reaction flux

    • .ub - k x 1 double upper bound on reaction flux

    • .labels - k x 1 cell array of labels to display on y axis

    • .Properties - table metadata; .Properties.Description is used in the plot legend

    • .zerobool - computed boolean, predicted flux is zero

    • .posbool - computed boolean, predicted flux is positive (secretion)

    • .negbool - computed boolean, predicted flux is negative (uptake)

    • .uptakeProbability - computed probability that the reaction is an uptake

OPTIONAL INPUTS:

  • otherFluxes – table with the same variables as predictedFluxes or measuredFluxes:

    • .rxns - k x 1 cell array of reaction identifiers

    • .Properties - table metadata; .Properties.Description is used in the plot legend

    • .mean - k x 1 double of measured mean flux (if provided)

    • .SD - k x 1 double standard deviation (if provided)

    • .v - k x 1 double of predicted mean flux (if provided)

    • .lb - k x 1 double lower bound on reaction flux (if provided)

    • .ub - k x 1 double upper bound on reaction flux (if provided)

  • extraFluxes – an additional table with the same variables as otherFluxes

  • param – parameters structure with the following fields:

    • .measuredFluxes - {(1), 0} set to zero to not display experimental fluxes

    • .labelType - determines which y axis labels to display: ‘metabolitePlatform’ = platform left and metabolite name right, ‘metabolite’ = metabolite only right (default), ‘platform’ = platform only left

    • .expOrder - {(1), 0} 1 = order reactions by experimental flux, 0 = order by predicted flux

    • .saveFigures - set to 1 to save the figure as .fig and .png to the current directory

    • .plotTitle - title of the plot

    • .measuredBounds - Logical, plot error bars on the measured fluxes

    • .predictedBounds - Logical, plot bounds on the predicted fluxes

    • .otherBounds - Logical, plot bounds on otherFluxes

    • .extraBounds - Logical, plot bounds on extraFluxes

plotPredictedAndExperimentalReactionRates[source]

Plot experimental uptake and secretion rates as well as error bars on logarithmic scale options to be specified and loaded before: exoMet = glcValidationData; fullReport = comparisonData_C1.fullReport; condition = ‘Complex I inhibition’; objective = ‘unWeightedTCBMfluxConc’; comparisonObjective=[]; labelType = ‘metabolite’; %’platform’ saveFigures = 1; comparison = comparisonData_glc.fullReport; % second predicted flux to be added for comparison comparison_label = ‘Control’;

plotPredictedExchange(predictedFluxes, otherFluxes, extraFluxes, param)[source]

Plot predicted exchange reaction rates, optionally compared with other flux sets

USAGE:

rxnsOut = plotPredictedExchange (predictedFluxes, otherFluxes, extraFluxes, param)

INPUTS:

predictedFluxes – table with the predicted fluxes, with variables:

  • .rxns - k x 1 cell array of reaction identifiers

  • .v - k x 1 double of predicted mean flux

  • .lb - k x 1 double lower bound on reaction flux

  • .ub - k x 1 double upper bound on reaction flux

  • .metNames - k x 1 metabolite names, used as labels when present

  • .labels - k x 1 cell array of labels to display on y axis

  • .Properties - table metadata; .Properties.Description is used in the plot legend

  • .zerobool - computed boolean, predicted flux is zero

  • .posbool - computed boolean, predicted flux is positive (secretion)

  • .negbool - computed boolean, predicted flux is negative (uptake)

  • .uptakeProbability - computed probability that the reaction is an uptake

OPTIONAL INPUTS:

  • otherFluxes – table with the same variables as predictedFluxes:

    • .rxns - k x 1 cell array of reaction identifiers

    • .Properties - table metadata; .Properties.Description is used in the plot legend

    • .mean - k x 1 double of measured mean flux (if provided)

    • .SD - k x 1 double standard deviation (if provided)

    • .v - k x 1 double of predicted mean flux (if provided)

    • .lb - k x 1 double lower bound on reaction flux (if provided)

    • .ub - k x 1 double upper bound on reaction flux (if provided)

  • extraFluxes – an additional table with the same variables as otherFluxes

  • param – parameters structure with the following fields:

    • .rxns - subset of reactions to plot

    • .labelType - determines which y axis labels to display: ‘metabolitePlatform’ = platform left and metabolite name right, ‘metabolite’ = metabolite only right (default), ‘platform’ = platform only left

    • .expOrder - {(1), 0} 1 = order reactions by experimental flux, 0 = order by predicted flux

    • .saveFigures - set to 1 to save the figure as .fig and .png to the current directory

    • .plotTitle - title of the plot

    • .exchangesOnly - Logical, restrict the plot to exchange reactions

    • .measuredBounds - Logical, plot error bars on the measured fluxes

    • .predictedBounds - Logical, plot bounds on the predicted fluxes

    • .otherBounds - Logical, plot bounds on otherFluxes

    • .extraBounds - Logical, plot bounds on extraFluxes

OUTPUT:

rxnsOut – the reaction identifiers that were plotted

summaryConcentrations(cleanedData, param)[source]

function to summarise metabolite concentrations per group/perturbation (calculate mean and SD), function can detect and remove outliers before data is summarised (param.removeOutliers). Furthermore, function can map metabolite names into model IDs using either user specified metabolite lookup table (param.metLookupTable) or/and model structure (param.model).

USAGE:

[summaryConc, outliers] = summaryConcentrations (cleanedData, param)

INPUTS:

cleanedData – metabolomics data in long format (output from mzQuality; required columns: sample, type, compound, concentration, cellLine or group, perturbation), with fields:

  • .Properties - table metadata (e.g. .Properties.VariableNames)

  • .group - group/cell-line label for each measurement

  • .perturbation - perturbation label for each measurement

  • .compound - measured compound (metabolite) name

  • .sample - sample identifier for each measurement

OPTIONAL INPUTS:

param – parameters structure with the following fields:

  • .model - COBRA model used to map metabolite names to model metabolites

  • .metLookupTable - table of metabolite names and model identifiers to map onto (required columns: ChemicalName, AlternativeName, Mets)

  • .excludeGroup - groups to be excluded from the analysis (default: none)

  • .group - groups to be analysed (default: all unique groups in cleanedData)

  • .printLevel - set greater than zero to print plots with median concentration values (default 0)

  • .perturbation - perturbations to be summarised separately (default: all unique perturbations)

  • .removeOutliers - Logical, find and remove outliers from the dataset (default 0)

OUTPUTS:

  • summaryConc – table of summarised metabolite concentrations (mean and SD) per group/perturbation

  • outliers – the outlier measurements detected and removed