Design


GDLS(model, targetRxns, varargin)

GDLS (Genetic Design through Local Search) attempts to find genetic designs with greater in silico production of desired metabolites.

USAGE:

[gdlsSolution, bilevelMILPProblem, gdlsSolutionStructs] = GDLS(model, varargin)

INPUTS:

model: Cobra model structure targetRxn: Reaction(s) to be maximized (Cell array of strings)

OPTIONAL INPUTS:
varargin: parameters entered using either a structure or list of

parameter, parameter value. List of optional parameters:

  • nbhdsz - Neighborhood size (default: 1)

  • M - Number of search paths (default: 1)

  • maxKO - Maximum number of knockouts (default: 50)

  • koCost - Cost for knocking out a reaction, gene set, or gene. A different cost can be set for each knockout (default: 1 for each knockout)

  • selectedRxns - List of reactions/geneSets that can be knocked out

  • koType - What to knockout: reactions, gene sets, or genes {(‘rxns’), ‘geneSets’, ‘genes’}

  • iterationLimit - Maximum number of iterations (default: 70)

  • timeLimit - Maximum run time in seconds (default: 252000)

  • minGrowth - Minimum growth rate

OUTPUTS:

gdlsSolution: GDLS solution structure (similar to OptKnock sol struct) bilevelMILPProblem: Problem structure used in computation gdlsSolutionStructs:

GetOptGeneSol(model, targetRxn, substrateRxn, generxnList, population, x, scores, isGeneList, saveFile, outputFolder)

Saves the solution from optGene and optGeneR in same format as OptKnock

USAGE:

[optGeneSol] = GetOptGeneSol(model, targetRxn, substrateRxn, generxnList, population, x, scores, isGeneList)

INPUTS:

model: model structure targetRxn: target reactions substrateRxn: substrate reactions generxnList: List of genes or rxns which can be knocked out population: population matrix x: the best solution scores: an array of scores isGeneList: boolean saveFile: boolean. Default = false; outputFolder: char. Default = pwd;

OUTPUT:

optGeneSol: Solution in the desired format

OptKnock(model, selectedRxnList, options, constrOpt, prevSolutions, verbFlag, solutionFileNameTmp)

Runs OptKnock in the most general form

USAGE:

OptKnock(model, selectedRxnList, options, constrOpt, prevSolutions, verbFlag, solutionFileNameTmp)

INPUTS:
model: Structure containing all necessary variables to described a

stoichiometric model

  • rxns - Rxns in the model

  • mets - Metabolites in the model

  • S - Stoichiometric matrix (sparse)

  • b - RHS of Sv = b (usually zeros)

  • c - Objective coefficients

  • lb - Lower bounds for fluxes

  • ub - Upper bounds for fluxes

  • rev - Reversibility of fluxes

selectedRxnList: n x 1 cell array where each cell is a reaction in model.rxns that can be knocked-out in OptKnock

OPTIONAL INPUTS:

options: OptKnock options

  • targetRxn - Target flux to be maximized

  • numDel - # of deletions allowed (Default: 5)

  • numDelSense - Direction of # of deletions constraint (G/E/L) (Default: L)

  • vMax - Max flux (Default: 1000)

  • solveOptKnock - Solve problem within Matlab (Default: true)

  • createGams - Create GAMS input file

  • gamsFile - GAMS input file name

constrOpt: Explicitly constrained reaction options

  • rxnList - Reaction list

  • values - Values for constrained reactions

  • sense - Constraint senses for constrained reactions (G/E/L)

prevSolutions: Previous solutions verbFlag: Verbose flag solutionFileNameTmp: File name for storing temporary solutions

OUTPUTS:
optKnockSol: optKnock solution structure

This is the same as the solution Structure from solveCobraMILP plus the following fields (if a solution exists): .rxnList - Reaction KO list .fluxes - the flux distribution

bilevelMILPproblem: optKnock problem structure

OptKnock uses bounds of -vMax to vMax or 0 to vMax for reversible and irreversible reactions. If you wish to constrain a reaction, use constrOpt.

analyzeGCdesign(modelRed, selectedRxns, target, deletions, maxKOs, objFunction, delPenalty, intermediateSlns)

Analyzes results with replacement knockouts should get closer to local maxima. Must have num KOs > 1.

USAGE:

[improvedRxns, intermediateSlns] = analyzeGCdesign(modelRed, selectedRxns, target, deletions, maxKOs, objFunction, delPenalty, intermediateSlns)

INPUTS:

modelRed: reduced model selectedRxns: selected reaction list from the reduced model target: exchange rxn to optimize deletions: initial set of KO rxns (must have at least 1 rxn)

OPTIONAL INPUTS:

maxKOs: maximum number of rxn KOs to allow (Default = 10) objFunction: pick an objective function to use (Default = 1):

  1. obj = maxRate (yield)

  2. obj = growth*maxRate (SSP)

  3. obj = maxRate*(delPenalty^numDels) (yield with KO penalty)

  4. obj = growth*maxRate*(delPenalty^numDels) (SSP with KO penalty)

  5. obj = maxRate*(slope^(-1)) (GC_yield)

  6. obj = growth*maxRate*(slope^(-1)) (GC_SSP)

  7. obj = maxRate*(delPenalty^numDels)*(slope^(-1)) (GC_yield with KO penalty)

  8. obj = growth*maxRate*(delPenalty^numDels)*(slope^(-1)) (GC_SSP with KO penalty)

delPenalty: penalty on extra rxn deletions (Default = .99) intermediateSlns: Previous set of solutions (Default = deletions)

OUTPUTS:

improvedRxns: the KO rxns for an improved strain intermediateSlns: all the sets of best KO rxns that are picked before the

final set is reached

analyzeOptKnock(model, deletions, target, biomassRxn, geneDelFlag)

Determines whether an optknock solution is growth coupled or not and what the maximum growth and production rates are

USAGE:

[type, maxGrowth, maxProd, minProd] = analyzeOptKnock(model, deletions, target, biomassRxn, geneDelFlag)

INPUTS:

model: COBRA model structure deletions: list of reaction or gene deletions (empty if wild type) target: the exchange reaction for the OptKnock target metabolite

OPTIONAL INPUTS:
biomassRxn: the biomass reaction name (Default = whatever is defined in

the model)

geneDelFlag: perform gene and not reaction deletions (Default = false)

OUTPUTS:

type: the type of OptKnock solution (growth coupled or not) maxGrowth: the maximum growth rate of the knockout strain maxProd: the maximum production rate of the target compound at the

maximum growth rate

minProd: the minimum production rate of the target compound at the

maximum growth rate

analyzeRxns(product, listProducts, listRxns)

Determines which knockout reactions occur most often when a specified product is produced

USAGE:

[allRxns, rxnCount] = analyzeRxns(product, listProducts, listRxns)

INPUTS:

product: the product to investigate listProducts: the list of all products produced in a RandKnock listRxns: the list of all rxns knocked out in a RandKnock

OUTPUTS:

allRxns: all of the rxns knocked out in strains producing the product rxnCount: the number of times each rxn was knocked out

augmentBOF(model, targetRxn, epsilon, printLevel)

Adjusts the objective function to eliminate “non-unique” optknock solutions by favoring the lowest production rate at a given predicted max growth-rate.

USAGE:

[model, rxn_name] = augmentBOF(model, targetRxn, epsilon)

INPUTS:
model: Structure containing all necessary variables to described a

stoichiometric model

targetRxn: objective of the optimization

OPTIONAL INPUT:

epsilon: degree of augmentation considering the biochemical objective printLevel: determine output level. (default: 0)

OUTPUTS:

model: Augmented model structure rxn_name: reaction that carries the augmented value

This funtion uses the outermembrane transport reaction, therefore:

  1. the model must have an extracellular compartment and a periplasm (e.g. iAF1260)

  2. there should not be an extracellular reaction acting on the metabolite besides the exchange reaction and the OM transport reaction

createBilevelMILPproblem(model, cLinear, cInteger, selRxns, selRxnMatch, constrOpt, measOpt, options, selPrevSol)

Creates the necessary matrices and vectors to solve a bilevel MILP with designated inner and outer problem objective functions

USAGE:

bilevelMILPproblem = createBilevelMILPProblem(model, cLinear, cInteger, selRxns, selRxnMatch, constrOpt, measOpt, options, selPrevSol);

INPUTS:

model: Model in irreversible format cLinear: Objective for linear part of the MILP problem (i.e. for fluxes) cInteger: Objective for integer part of the MILP problem selRxns: Reactions that participate in the integer part (e.g. ones

that can be deleted) (in the form [0 0 1 0 0 1 0 1 1 0 1])

selRxnMatch: Matching of the forward and reverse parts constrOpt: Constraint options

  • rxnInd

  • values

  • sense

measOpt: Measured flux options

  • rxnSel

  • values

  • weights

options: General options

  • vMax - Maximal/minimal value for cont variables

  • numDel - Number of deletions

  • numDelSense - # of ko <=/=/>= K (L/E/G)

selPrevSol: Previous solutions (optional)

OUTPUTS:

bilevelMILPproblem: struct with:

  • A - LHS matrix

  • b - RHS

  • c - Objective

  • csense - Constraint types

  • lb - Lower bounds

  • ub - Upper bounds

  • vartype - Variable types

  • contSolInd - Allows selecting the continuous solution (i.e. fluxes)

  • intsSolInd - Allows selecting the integer solution (i.e. what reactions are used)

doubleProductionEnvelope(model, deletions, prod1, prod2, biomassRxn, geneDelFlag, nPts)

Plots maximum growth rate as a function of the output of two specified products

USAGE:

[x1, x2, y] = doubleProductionEnvelope(model, deletions, prod1, prod2, biomassRxn, geneDelFlag, nPts)

INPUTS:

model: COBRA model structure deletions: The reactions or genes to knockout of the model prod1: One of the two products to investigate prod2: The other product to investigate

OPTIONAL INPUTS:
biomassRxn: The biomass objective function rxn name

(Default = ‘biomass_SC4_bal’)

geneDelFlag: Perform gene and not reaction deletions

(Default = false)

nPts: Number of points to plot for each product

(Default = 20)

OUTPUTS:

x1: The range of rates plotted for prod1 x2: The range of rates plotted for prod2 y: The plotted growth rates at each (x1, x2)

multiProductionEnvelope(model, deletions, biomassRxn, geneDelFlag, nPts, plotAllFlag, plotTools)

Calculates the byproduct secretion envelopes for every product (excreted metabolites with 1 or more Carbons)

USAGE:

[biomassValues, targetValues] = multiProductionEnvelope(model, deletions, biomassRxn, geneDelFlag, nPts, plotTools)

INPUT:

model: COBRA model structure

OPTIONAL INPUT:
deletions: List of reaction or gene deletions (empty if wild type)

(Default = {})

biomassRxn: Biomass rxn name (Default = whatever is defined in model) geneDelFlag: Perform gene and not reaction deletions (Default = false) nPts: Number of points in the plot (Default = 20) plotTools: boolean (default = false) - add tools for editing the figure and its properties plotAllFlag: plot all envelopes, even ones that are not growth

coupled (Default = false)

plotTools: boolean (default = false) - add tools for editing the figure and its properties

OUTPUT:

biomassValues: Biomass values for plotting targetUpperBounds: Target upper bounds for plotting (

biomassvalues x reactions)

targetLowerBounds: Target lower bounds for plotting (

biomassvalues x reactions)

plottedReactions: Reactions that led to relevant side product

multiProductionEnvelopeInorg(model, deletions, biomassRxn, geneDelFlag, nPts, plotAllFlag, plotTools)

Calculates the byproduct secretion envelopes for every product, including inorganic compounds

USAGE:

[biomassValues, targetValues] = multiProductionEnvelopeInorg(model, deletions, biomassRxn, geneDelFlag, nPts, plotAllFlag, plotTools)

INPUT:

model: COBRA model structure

OPTIONAL INPUT:
deletions: List of reaction or gene deletions (empty if wild type)

(Default = {})

biomassRxn: Biomass rxn name (Default = whatever is defined in model) geneDelFlag: Perform gene and not reaction deletions (Default = false) nPts: Number of points in the plot (Default = 20) plotAllFlag: Plots all envelopes, even ones that are not growth coupled

(Default = false)

plotTools: boolean (default = false) - add tools for editing the figure and its properties

OUTPUT:

biomassValues: Biomass values for plotting targetValues: Target upper and lower bounds for plotting

optGene(model, targetRxn, substrateRxn, generxnList, varargin)

Implements the optgene algorithm.

USAGE:

[x, population, scores, optGeneSol] = optGene(model, targetRxn, substrateRxn, generxnList, MaxKOs, population)

INPUTS:

model: Model of reconstruction targetRxn: (char) String name of reaction which is to

be maximized

substrateRxn: (char) Substrate reactions generxnList: (cell array) List of genes or rxns which

can be knocked out. The program will guess which of the two it is, based on the content in model.

OPTIONAL INPUTS:

MaxKOs: (double) Maximal KnockOuts population: (logical matrix) population matrix. Use this

parameter to interrupt simulation and resume afterwards.

mutationRate: (double) the rate of mutation. crossovermutationRate: (double) the rate of mutation after a

crossover. This value should probably be fairly low. It is only there to ensure that not every member of the population ends up with the same genotype.

CrossoverFraction: (double) Percentage of offspring created by

crossing over (as opposed to mutation). 0.7 - 0.8 were found to generate the highest mean, but this can be adjusted.

PopulationSize: (double) Number of individuals Generations: (double) Maximum number of generations TimeLimit: (double) global time limit in seconds StallTimeLimit: (double) Stall time limit (terminate after

this much time of not finding an improvement in fitness)

StallGenLimit: (double) terminate after this many

generations of not finding an improvement

MigrationFraction: (double). how many individuals migrate MigrationInterval: (double). how often individuals migrate from

one population to another.

saveFile: (double or boolean) saving a file with

inputs and outputs. Default = false;

outputFolder: (char) name of folder where

files will be generated

OUTPUTS:

x: best optimized value found population: Population of individuals. Pass this back

into optgene to continue simulating where you left off.

scores: An array of scores optGeneSol: optGene solution strcture

optGeneFitness(rxn_vector_matrix, model, targetRxn, rxnListInput, isGeneList)

The fitness function

USAGE:

[val] = optGeneFitness(rxn_vector_matrix, model, targetRxn, rxnListInput, isGeneList)

INPUTS:

rxn_vector_matrix: reactions vectors in a matrix model: model structure targetRxn: target reactions rxnListInput: list of reactions isGeneList: bolean checking if it is a gene list

OUTPUT:

val: fitness value

optGeneFitnessTilt(rxn_vector_matrix, model, targetRxn, rxnListInput, isGeneList)

The fitness function

USAGE:

[val] = optGeneFitnessTilt(rxn_vector_matrix, model, targetRxn, rxnListInput, isGeneList)

INPUTS:

rxn_vector_matrix: reactions vectors in a matrix model: model structure targetRxn: target reactions rxnListInput: list of reactions isGeneList: bolean checking if it is a gene list

OUTPUT:

val: fitness value

productionEnvelope(model, deletions, lineColor, targetRxn, biomassRxn, geneDelFlag, nPts)

Calculates the byproduct secretion envelope

USAGE:

[biomassValues, targetValues, lineHandle] = productionEnvelope(model, deletions, lineColor, targetRxn, biomassRxn, geneDelFlag, nPts)

INPUTS:

model: COBRA model structure

OPTIONAL INPUTS:

deletions: List of reaction or gene deletions (empty if wild type) lineColor: Line color for plotting (see help plot for color codes) targetRxn: Target metabolite production reaction name biomassRxn: Biomass rxn name geneDelFlag: Perform gene and not reaction deletions nPts: Number of points in the plot

OUTPUTS:

biomassValues: Biomass values for plotting targetValues: Target upper and lower bounds for plotting lineHandle: Handle to lineseries object

randomKO(modelRed, selectedRxns, N)

Knocks out N random genes and reports products from FBA

USAGE:

[products, productRates, KOrxns, BOF] = randomKO(modelRed, selectedRxns, N)

INPUTS:

modelRed: a reduced model (from the reduceModel.m function) selectedRxns: the reactions eligible for deletion (from the

getOptKnockTargets.m function)

N: the number of reactions to randomly knockout

OUTPUTS:

products: the exchange reactions that produce a siginifcant output productRates: the rates of those exhange reactions KOrxns: the N reactions randomly knocked out BOF: the value of the biomass objective function of the

knockout strain

simpleOptKnock(model, targetRxn, deletions, geneDelFlag, minGrowth, doubleDelFlag)

Simple OptKnock is used to check all one gene or reaction deletions for growth-coupled metabolite production

USAGE:

[wtRes, delRes] = simpleOptKnock(model, targetRxn, deletions, geneDelFlag, minGrowth, doubleDelFlag)

INPUTS:

model: COBRA model structure targetRxn: Target metabolite production reaction

OPTIONAL INPUTS:
deletions: Set of gene or reaction deletions to consider for KO

(Default = all reactions)

geneDelFlag: Gene deletion flag (Default = false) minGrowth: Minimum KO growth rate (Default = 0.05) doubleDelFlag: Double deletions (Default = false)

OUTPUTS:

wtRes: Wild type results delRes: Deletion strain results

The results structures have fields:

  • growth - Growth rate of strain

  • minProd - Minimum prod rate of target metabolite

  • maxProd - Maximum prod rate of target metabolite

singleProductionEnvelope(model, deletions, product, biomassRxn, varargin)

Plots maximum growth rate as a function of the output of one specified products

USAGE:

singleProductionEnvelope(model, deletions, product, biomassRxn, geneDelFlag, nPts)

INPUTS:
model: (struct) a metabolic model with at least

the following fields:

  • .rxns - Reaction IDs in the model

  • .mets - Metabolite IDs in the model

  • .S - Stoichiometric matrix (sparse)

  • .b - RHS of Sv = b (usually zeros)

  • .c - Objective coefficients

  • .lb - Lower bounds for fluxes

  • .ub - Upper bounds for fluxes

deletions: (cell array) the reactions or genes to knockout of

the model

product: (char) the product to investigate biomassRxn: (char) the biomass objective function

OPTIONAL INPUTS:
geneDelFlag: (double or boolean) perform gene and not reaction

deletions Default = false

nPts: (double) number of points to plot for each product

Default = 20

savePlot: (double or boolean) boolean for saving plot in a file

Default = false

showPlot: (double or boolean) boolean for showing the plot

Default = false

fileName: (char) name of the file where the plot is saved.

Default = product

outputFolder: (char) name of the folder where files are saved

Default = ‘Results’

OUTPUTS;

x: x axis for the curves y1: y results for the first curve - ‘Minimum Wild-type’ y2: y results for the second curve - ‘Maximum Wild-type’

testOptKnockSol(model, targetRxn, deletions)

Tests an OptKnock knockout strain

USAGE:

[growthRate, minProd, maxProd] = testOptKnockSol(model, targetRxn, deletions)

INPUTS:

model: COBRA model structure targetRxn: Target reaction (e.g. ‘EX_etoh(e)’)

OPTIONAL INPUT:
deletions: Set of reaction deletions (e.g. {‘PGI’, ‘TPI’})

(Default = [])

OUTPUTS:

growthRate: Maximim growth rate of the strain minProd: Minimum production rate at max growth rate maxProd: Maximum production rate at max growth rate

theoretMaxProd(model, inputrxn, criterion, normalize, rxns)

Determines the max theoretical output for each exchange reaction

USAGE:

[ExRxns, MaxTheoOut]= theoreticalMaxProduction(model, inputrxn, criterion, normalize, rxns)

INPUT:

model: model structure inputrxn: the input reaction (‘EX_glu(e)’, etc.)

OPTIONAL INPUT:

criterion: One of

  • ‘pr_mol’ (default, return the maximal molar fluxes for all exchangers)

  • ‘pr_mw’ (return the maximal molecular weight fluxes)

  • ‘pr_other_mol’ (return the molar fluxes through other carbon Exporters)

  • ‘pr_other_mw’ (return the molecular weights exported through other carbon Exporters)

normalize: normalize by input flux with respect to the input Flux. Either the flux rate in mol or

in molecular weight (Default = false)

rxns: Selection Vector (1 for selected, 0 otherwise)

OUTPUTS:

ExRxns: Vector of exchange reactions MaxThroOut: The max theoretical output for each exchange reaction