Optenvelope

addEnv(origModel, biomass, desiredProduct, varargin)[source]

addEnv adds a production envelope to the current figure

Adds a production envelope (with optional knockouts) to the active figure. The algorithm is able to knock out genes as well as reactions to produce the production envelope.

USAGE:

line = addEnv (origModel, biomass, desiredProduct, varargin)

INPUTS:

  • origModel – COBRA model structure with fields:

    • .S - Stoichiometric matrix

    • .rxns - Reaction identifiers

    • .mets - Metabolite identifiers

    • .lb - Lower bounds

    • .ub - Upper bounds

    • .b - Right hand side values for metabolite constraints

    • .c - Objective coefficients

  • biomass – Reaction name of biomass [char]

  • desiredProduct – Reaction name of desired product [char]

OPTIONAL INPUTS:

  • KnockOuts – List of knockouts (reaction or gene IDs) for the production envelope [cell array] (default: {})

  • colour – Colour of the plotted line (any valid MATLAB colour) (default: ‘r’)

  • prodMol – Molar mass of the target product for a yield plot [double] (default: [])

  • subUptake – Uptake of substrate for a yield plot [double] (default: 10)

  • molarSum – Molar mass of substrate for a yield plot [double] (default: 180)

OUTPUTS:

line – Line object returned by plot for the maximum edge of the envelope

Note

Sometimes the last point of the envelope drops to zero (possibly a rounding error); this function connects the last points of the lines so the graph forms a continuous line. This algorithm only adds the graph; it does not change labels.

Example

line = addEnv(model, ‘BIOMASS_Ecoli’, ‘EX_ac_e’, {‘GHMT2r’, ‘GND’}, ‘m’)

milpOEReinserts(model, data, K, minP, numKO, toDel, timeLimit, printLevel)[source]

milpOEReinserts finds a set of knockouts using a MILP

Creates and solves a MILP to find a given number of knockouts that achieve the best possible production envelope.

USAGE:

[knockouts] = milpOEReinserts (model, data, K, minP, numKO, toDel, timeLimit, printLevel)

INPUTS:

  • model – COBRA model structure in irreversible form with fields:

    • .rxns - Reaction identifiers

    • .S - Stoichiometric matrix

    • .lb - Lower bounds

    • .ub - Upper bounds

    • .c - Objective coefficients

    • .C_chemical - objective coefficient vector marking the desired product, set internally by this function

  • data – Struct with the field:

    • .mainActive - List of active reactions for the main envelope [cell array]

  • K – List of reactions that cannot be selected for knockout (reaction IDs) [double array]

  • minP – Struct with the fields:

    • .bioID - Reaction ID of biomass [double]

    • .proID - Reaction ID of desired product [double]

  • numKO – Number of knockouts to achieve [double]

OPTIONAL INPUTS:

  • toDel – Numeric flag for what to delete: 0 = reactions, 1 = genes, 2 = enzymes (default: 0)

  • timeLimit – Time limit for the MILP solver in seconds [double] (default: inf)

  • printLevel – Print level for the MILP solver [double] (default: 0)

OUTPUTS:

knockouts – List of reactions that, when removed, give the optimal envelope [cell array]

Note

This function is not designed for stand-alone use. It should be used through optEnvelope.m with the numKO parameter set.

minActiveRxns(model, matchRev, K, minP, toDel, timeLimit, midPoints, printLevel)[source]

minActiveRxns finds the minimum set of active reactions for an envelope

Finds a minimum-size set of reactions that must be active at a specific point on the production envelope to obtain the best possible production envelope.

USAGE:

[data] = minActiveRxns (model, matchRev, K, minP, toDel, timeLimit, midPoints, printLevel)

INPUTS:

  • model – COBRA model structure in irreversible form with fields:

    • .S - Stoichiometric matrix

    • .rxns - Reaction identifiers

    • .b - Right hand side values for metabolite constraints

    • .lb - Lower bounds

    • .ub - Upper bounds

    • .genes - Gene identifiers

    • .rxnGeneMat - Reaction-gene incidence matrix

  • matchRev – Matching of forward and backward reactions of a reversible reaction [double array]

  • K – List of reactions that cannot be selected for knockout (reaction IDs) [double array]

  • minP – Struct with the fields:

    • .bioID - Reaction ID of biomass [double]

    • .bioMin - 1% of the max biomass in the wild-type [double]

    • .proID - Reaction ID of desired product [double]

    • .proMin - 1% of the max desired product in the wild-type [double]

    • .proMax - Max desired product in the wild-type [double]

OPTIONAL INPUTS:

  • toDel – Numeric flag for what to delete: 0 = reactions, 1 = genes, 2 = enzymes (default: 0)

  • timeLimit – Time limit for the MILP solver in seconds [double] (default: inf)

  • midPoints – Number of mid points to calculate active reactions for [double] (default: 0)

  • printLevel – Print level for the MILP solver [double] (default: 0)

OUTPUTS:

data – Struct with the fields:

  • .pro - Minimal product for the best envelope

  • .bio - Maximum biomass for the best envelope

  • .results - Results of the MILP for each point

  • .mainModel - Model for the main envelope

  • .mainActive - List of active reactions for the main envelope

  • .models - Models for all mid envelopes

  • .active - Lists of active reactions for mid envelopes

Note

This function is designed to be used with optEnvelope and was not designed as a stand-alone function.

optEnvelope(model, desiredProduct, varargin)[source]

optEnvelope finds a minimal set of knockouts for a production envelope

Uses a MILP to find the minimum set of active reactions and then finds the smallest set of reactions in the pool of inactive reactions that offers the same production envelope. The algorithm provides multiple ways to reinsert reactions: sequential, MILP and GA (under construction).

USAGE:

[main, mid] = optEnvelope (model, desiredProduct, varargin)

INPUTS:

  • model – COBRA model structure with fields:

    • .S - Stoichiometric matrix

    • .rxns - Reaction identifiers

    • .mets - Metabolite identifiers

    • .lb - Lower bounds

    • .ub - Upper bounds

    • .b - Right hand side values for metabolite constraints

    • .c - Objective coefficients

    • .grRules - Readable gene-protein-reaction rules

    • .metNames - Metabolite names

    • .metFormulas - Elemental formulas

  • desiredProduct – Reaction name of desired product [char]

OPTIONAL INPUTS:

varargin – Parameters given as parameter name / value pairs:

  • protectedRxns - Additional reactions to ignore (must be in irreversible form) [cell array] (default: {})

  • numTries - Iterations for finding the best set of deletions [double] (default: [])

  • numKO - Number of reactions to remove for the final result (triggers MILP reinsertion) [double] (default: [])

  • prodMol - Molar mass of product for a yield plot (g/mol) [double] (default: [])

  • midPoints - Number of points to check along the edge for the best envelope [double] (default: 0)

  • timeLimit - Time limit for the solver in seconds (also limits numTries) [double] (default: inf)

  • printLevel - Print level for the solver [double] (default: 0)

  • drawEnvelope - Whether the algorithm should draw envelopes [logical] (default: true)

  • delGenes - Delete genes (unfinished) [logical] (default: false)

  • delEnzymes - Delete enzymes (unfinished) [logical] (default: false)

  • GAon - Use genetic algorithm (unfinished) [logical] (default: false)

OUTPUTS:

  • main – Struct with information about the reactions to remove for the optimal envelope and its most probable point (fields .knockouts and .peak)

  • mid – Struct with information about the reactions to remove for the midpoint envelopes and their most probable points (fields .midKnockoutsTable and .peak)

Example

[mainKnockouts, midKnockouts] = optEnvelope(model, ‘EX_ac_e’, ‘timeLimit’, 600, ‘midPoints’, 15)

Note

A figure (desired product versus biomass) including plots for the wild-type and the opt envelope is presented after running optEnvelope. Mid envelopes currently work only for sequential (default) reinsertions.

seperateTransposeJoinOE(A, Ay, B, C, ySize, yMax, aSizeCol, highNum, zSize)[source]

seperateTransposeJoinOE transposes and joins arrays for a MILP

Helper for milpOEReinserts.m that transposes and joins the constraint arrays needed to build the dual MILP used for reaction reinsertion.

USAGE:

[A_w, Ay_w, B_w, C_w, lb_w, ub_w, wSize] = seperateTransposeJoinOE (A, Ay, B, C, ySize, yMax, aSizeCol, highNum, zSize)

INPUTS:

  • A – Constraint matrix to transpose and join [matrix]

  • Ay – Constraint matrix associated with the y (knockout) variables [matrix]

  • B – Right hand side vector of the constraints [double array]

  • C – Objective coefficient vector [double array]

  • ySize – Number of y (knockout) variables [double]

  • yMax – Upper bound for the y variables [double]

  • aSizeCol – Number of columns of A [double]

  • highNum – Large finite number used as a big-M bound [double]

  • zSize – Number of z (auxiliary) variables [double]

OUTPUTS:

  • A_w – Transposed and joined constraint matrix [matrix]

  • Ay_w – Joined matrix for the y variables [matrix]

  • B_w – Joined right hand side vector [double array]

  • C_w – Joined objective coefficient vector [double array]

  • lb_w – Lower bounds for the joined problem [double array]

  • ub_w – Upper bounds for the joined problem [double array]

  • wSize – Number of w (dual) variables [double]

Note

This function is not designed for stand-alone use. It is used by milpOEReinserts.m.

sequentialOEReinserts(modelOld, data, K, toDel, minP, midPoints, numTries, timeLimit)[source]

sequentialOEReinserts reinserts reactions one by one to minimise knockouts

Goes through the inactive reactions sequentially and reinserts them one by one to obtain the best possible set of knockouts while retaining the optimal envelope. With the numTries parameter this can be repeated several times by randomizing the list of inactive reactions.

USAGE:

[mainKnockouts, finalMidKnockouts] = sequentialOEReinserts (modelOld, data, K, toDel, minP, midPoints, numTries, timeLimit)

INPUTS:

  • modelOld – COBRA model structure in irreversible form with fields:

    • .rxns - Reaction identifiers

    • .lb - Lower bounds

    • .ub - Upper bounds

  • data – Struct with the fields:

    • .mainActive - List of active reactions for the main envelope [cell array]

    • .mainModel - Model for the main envelope [struct]

    • .active - Lists of active reactions for mid envelopes [cell array]

    • .models - Models for all mid envelopes [struct]

  • K – List of reactions that cannot be selected for knockout (reaction IDs) [double array]

  • toDel – Numeric flag for what to delete: 0 = reactions, 1 = genes, 2 = enzymes

  • minP – Struct with the fields:

    • .bioID - Reaction ID of biomass [double]

    • .proID - Reaction ID of desired product [double]

  • midPoints – Number of middle points to calculate mid envelopes for [double]

  • numTries – Number of iterations for randomizing the indexes of inactive reactions [double]

  • timeLimit – Time limit for the solver in seconds [double]

OUTPUTS:

  • mainKnockouts – List of reactions that, when removed, give the optimal envelope [cell array]

  • finalMidKnockouts – List of reactions that, when removed, give the optimal envelopes for the middle envelopes [cell array]

Example

[mainKnockouts, finalMidKnockouts] = sequentialOEReinserts(model, data, K, 0, minP, 10, 10, 60)

Note

This function is not designed for stand-alone use. It should be used through optEnvelope.m.