17.1.1.2. cobra.flux_analysis

17.1.1.2.1. Submodules

17.1.1.2.2. Package Contents

17.1.1.2.2.1. Functions

double_gene_deletion(model, gene_list1=None, gene_list2=None, method=’fba’, solution=None, processes=None, **kwargs)

Knock out each gene pair from the combination of two given lists.

double_reaction_deletion(model, reaction_list1=None, reaction_list2=None, method=’fba’, solution=None, processes=None, **kwargs)

Knock out each reaction pair from the combinations of two given lists.

single_gene_deletion(model, gene_list=None, method=’fba’, solution=None, processes=None, **kwargs)

Knock out each gene from a given list.

single_reaction_deletion(model, reaction_list=None, method=’fba’, solution=None, processes=None, **kwargs)

Knock out each reaction from a given list.

fastcc(model, flux_threshold=1.0, zero_cutoff=None)

Check consistency of a metabolic network using FASTCC [1]_.

gapfill(model, universal=None, lower_bound=0.05, penalties=None, demand_reactions=True, exchange_reactions=False, iterations=1)

Perform gapfilling on a model.

geometric_fba(model, epsilon=1e-06, max_tries=200, processes=None)

Perform geometric FBA to obtain a unique, centered flux distribution.

loopless_solution(model, fluxes=None)

Convert an existing solution to a loopless one.

add_loopless(model, zero_cutoff=None)

Modify a model so all feasible flux distributions are loopless.

add_moma(model, solution=None, linear=True)

Add constraints and objective representing for MOMA.

moma(model, solution=None, linear=True)

Compute a single solution based on (linear) MOMA.

pfba(model, fraction_of_optimum=1.0, objective=None, reactions=None)

Perform basic pFBA (parsimonious Enzyme Usage Flux Balance Analysis)

find_blocked_reactions(model, reaction_list=None, zero_cutoff=None, open_exchanges=False, processes=None)

Find reactions that cannot carry any flux.

find_essential_genes(model, threshold=None, processes=None)

Return a set of essential genes.

find_essential_reactions(model, threshold=None, processes=None)

Return a set of essential reactions.

flux_variability_analysis(model, reaction_list=None, loopless=False, fraction_of_optimum=1.0, pfba_factor=None, processes=None)

Determine the minimum and maximum possible flux value for each reaction.

production_envelope(model, reactions, objective=None, carbon_sources=None, points=20, threshold=None)

Calculate the objective value conditioned on all combinations of

add_room(model, solution=None, linear=False, delta=0.03, epsilon=0.001)

Add constraints and objective for ROOM.

room(model, solution=None, linear=False, delta=0.03, epsilon=0.001)

Compute a single solution based on regulatory on/off minimization (ROOM).

cobra.flux_analysis.double_gene_deletion(model, gene_list1=None, gene_list2=None, method='fba', solution=None, processes=None, **kwargs)[source]

Knock out each gene pair from the combination of two given lists.

We say ‘pair’ here but the order order does not matter.

Parameters
  • model (cobra.Model) – The metabolic model to perform deletions in.

  • gene_list1 (iterable, optional) – First iterable of ``cobra.Gene``s to be deleted. If not passed, all the genes from the model are used.

  • gene_list2 (iterable, optional) – Second iterable of ``cobra.Gene``s to be deleted. If not passed, all the genes from the model are used.

  • method ({"fba", "moma", "linear moma", "room", "linear room"}, optional) – Method used to predict the growth rate.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference for (linear) MOMA or ROOM.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not passed, will be set to the number of CPUs found.

  • kwargs – Keyword arguments are passed on to underlying simulation functions such as add_room.

Returns

A representation of all combinations of gene deletions. The columns are ‘growth’ and ‘status’, where

indextuple(str)

The gene identifiers that were knocked out.

growthfloat

The growth rate of the adjusted model.

statusstr

The solution’s status.

Return type

pandas.DataFrame

cobra.flux_analysis.double_reaction_deletion(model, reaction_list1=None, reaction_list2=None, method='fba', solution=None, processes=None, **kwargs)[source]

Knock out each reaction pair from the combinations of two given lists.

We say ‘pair’ here but the order order does not matter.

Parameters
  • model (cobra.Model) – The metabolic model to perform deletions in.

  • reaction_list1 (iterable, optional) – First iterable of ``cobra.Reaction``s to be deleted. If not passed, all the reactions from the model are used.

  • reaction_list2 (iterable, optional) – Second iterable of ``cobra.Reaction``s to be deleted. If not passed, all the reactions from the model are used.

  • method ({"fba", "moma", "linear moma", "room", "linear room"}, optional) – Method used to predict the growth rate.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference for (linear) MOMA or ROOM.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not passed, will be set to the number of CPUs found.

  • kwargs – Keyword arguments are passed on to underlying simulation functions such as add_room.

Returns

A representation of all combinations of reaction deletions. The columns are ‘growth’ and ‘status’, where

indextuple(str)

The reaction identifiers that were knocked out.

growthfloat

The growth rate of the adjusted model.

statusstr

The solution’s status.

Return type

pandas.DataFrame

cobra.flux_analysis.single_gene_deletion(model, gene_list=None, method='fba', solution=None, processes=None, **kwargs)[source]

Knock out each gene from a given list.

Parameters
  • model (cobra.Model) – The metabolic model to perform deletions in.

  • gene_list (iterable) – ``cobra.Gene``s to be deleted. If not passed, all the genes from the model are used.

  • method ({"fba", "moma", "linear moma", "room", "linear room"}, optional) – Method used to predict the growth rate.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference for (linear) MOMA or ROOM.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not passed, will be set to the number of CPUs found.

  • kwargs – Keyword arguments are passed on to underlying simulation functions such as add_room.

Returns

A representation of all single gene deletions. The columns are ‘growth’ and ‘status’, where

indextuple(str)

The gene identifier that was knocked out.

growthfloat

The growth rate of the adjusted model.

statusstr

The solution’s status.

Return type

pandas.DataFrame

cobra.flux_analysis.single_reaction_deletion(model, reaction_list=None, method='fba', solution=None, processes=None, **kwargs)[source]

Knock out each reaction from a given list.

Parameters
  • model (cobra.Model) – The metabolic model to perform deletions in.

  • reaction_list (iterable, optional) – ``cobra.Reaction``s to be deleted. If not passed, all the reactions from the model are used.

  • method ({"fba", "moma", "linear moma", "room", "linear room"}, optional) – Method used to predict the growth rate.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference for (linear) MOMA or ROOM.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not passed, will be set to the number of CPUs found.

  • kwargs – Keyword arguments are passed on to underlying simulation functions such as add_room.

Returns

A representation of all single reaction deletions. The columns are ‘growth’ and ‘status’, where

indextuple(str)

The reaction identifier that was knocked out.

growthfloat

The growth rate of the adjusted model.

statusstr

The solution’s status.

Return type

pandas.DataFrame

cobra.flux_analysis.fastcc(model, flux_threshold=1.0, zero_cutoff=None)[source]

Check consistency of a metabolic network using FASTCC [1]_.

FASTCC (Fast Consistency Check) is an algorithm for rapid and efficient consistency check in metabolic networks. FASTCC is a pure LP implementation and is low on computation resource demand. FASTCC also circumvents the problem associated with reversible reactions for the purpose. Given a global model, it will generate a consistent global model i.e., remove blocked reactions. For more details on FASTCC, please check [1]_.

Parameters
  • model (cobra.Model) – The constraint-based model to operate on.

  • flux_threshold (float, optional (default 1.0)) – The flux threshold to consider.

  • zero_cutoff (float, optional) – The cutoff to consider for zero flux (default model.tolerance).

Returns

The consistent constraint-based model.

Return type

cobra.Model

Notes

The LP used for FASTCC is like so: maximize: sum_{i in J} z_i s.t. : z_i in [0, varepsilon] forall i in J, z_i in mathbb{R}_+

v_i ge z_i forall i in J Sv = 0 v in B

References

1

Vlassis N, Pacheco MP, Sauter T (2014) Fast Reconstruction of Compact Context-Specific Metabolic Network Models. PLoS Comput Biol 10(1): e1003424. doi:10.1371/journal.pcbi.1003424

cobra.flux_analysis.gapfill(model, universal=None, lower_bound=0.05, penalties=None, demand_reactions=True, exchange_reactions=False, iterations=1)[source]

Perform gapfilling on a model.

See documentation for the class GapFiller.

modelcobra.Model

The model to perform gap filling on.

universalcobra.Model, None

A universal model with reactions that can be used to complete the model. Only gapfill considering demand and exchange reactions if left missing.

lower_boundfloat

The minimally accepted flux for the objective in the filled model.

penaltiesdict, None

A dictionary with keys being ‘universal’ (all reactions included in the universal model), ‘exchange’ and ‘demand’ (all additionally added exchange and demand reactions) for the three reaction types. Can also have reaction identifiers for reaction specific costs. Defaults are 1, 100 and 1 respectively.

iterationsint

The number of rounds of gapfilling to perform. For every iteration, the penalty for every used reaction increases linearly. This way, the algorithm is encouraged to search for alternative solutions which may include previously used reactions. I.e., with enough iterations pathways including 10 steps will eventually be reported even if the shortest pathway is a single reaction.

exchange_reactionsbool

Consider adding exchange (uptake) reactions for all metabolites in the model.

demand_reactionsbool

Consider adding demand reactions for all metabolites.

iterable

list of lists with on set of reactions that completes the model per requested iteration.

import cobra as ct
>>> from cobra import Model
>>> from cobra.flux_analysis import gapfill
>>> model = ct.create_test_model("salmonella")
>>> universal = Model('universal')
>>> universal.add_reactions(model.reactions.GF6PTA.copy())
>>> model.remove_reactions([model.reactions.GF6PTA])
>>> gapfill(model, universal)
cobra.flux_analysis.geometric_fba(model, epsilon=1e-06, max_tries=200, processes=None)[source]

Perform geometric FBA to obtain a unique, centered flux distribution.

Geometric FBA [1]_ formulates the problem as a polyhedron and then solves it by bounding the convex hull of the polyhedron. The bounding forms a box around the convex hull which reduces with every iteration and extracts a unique solution in this way.

Parameters
  • model (cobra.Model) – The model to perform geometric FBA on.

  • epsilon (float, optional) – The convergence tolerance of the model (default 1E-06).

  • max_tries (int, optional) – Maximum number of iterations (default 200).

  • processes (int, optional) – The number of parallel processes to run. If not explicitly passed, will be set from the global configuration singleton.

Returns

The solution object containing all the constraints required for geometric FBA.

Return type

cobra.Solution

References

1

Smallbone, Kieran & Simeonidis, Vangelis. (2009). Flux balance analysis: A geometric perspective. Journal of theoretical biology.258. 311-5. 10.1016/j.jtbi.2009.01.027.

cobra.flux_analysis.loopless_solution(model, fluxes=None)[source]

Convert an existing solution to a loopless one.

Removes as many loops as possible (see Notes). Uses the method from CycleFreeFlux [1]_ and is much faster than add_loopless and should therefore be the preferred option to get loopless flux distributions.

Parameters
  • model (cobra.Model) – The model to which to add the constraints.

  • fluxes (dict) – A dictionary {rxn_id: flux} that assigns a flux to each reaction. If not None will use the provided flux values to obtain a close loopless solution.

Returns

A solution object containing the fluxes with the least amount of loops possible or None if the optimization failed (usually happening if the flux distribution in fluxes is infeasible).

Return type

cobra.Solution

Notes

The returned flux solution has the following properties:

  • it contains the minimal number of loops possible and no loops at all if all flux bounds include zero

  • it has an objective value close to the original one and the same objective value id the objective expression can not form a cycle (which is usually true since it consumes metabolites)

  • it has the same exact exchange fluxes as the previous solution

  • all fluxes have the same sign (flow in the same direction) as the previous solution

References

1

CycleFreeFlux: efficient removal of thermodynamically infeasible loops from flux distributions. Desouki AA, Jarre F, Gelius-Dietrich G, Lercher MJ. Bioinformatics. 2015 Jul 1;31(13):2159-65. doi: 10.1093/bioinformatics/btv096.

cobra.flux_analysis.add_loopless(model, zero_cutoff=None)[source]

Modify a model so all feasible flux distributions are loopless.

In most cases you probably want to use the much faster loopless_solution. May be used in cases where you want to add complex constraints and objecives (for instance quadratic objectives) to the model afterwards or use an approximation of Gibbs free energy directions in you model. Adds variables and constraints to a model which will disallow flux distributions with loops. The used formulation is described in [1]_. This function will modify your model.

Parameters
  • model (cobra.Model) – The model to which to add the constraints.

  • zero_cutoff (positive float, optional) – Cutoff used for null space. Coefficients with an absolute value smaller than zero_cutoff are considered to be zero (default model.tolerance).

Returns

Return type

Nothing

References

1

Elimination of thermodynamically infeasible loops in steady-state metabolic models. Schellenberger J, Lewis NE, Palsson BO. Biophys J. 2011 Feb 2;100(3):544-53. doi: 10.1016/j.bpj.2010.12.3707. Erratum in: Biophys J. 2011 Mar 2;100(5):1381.

cobra.flux_analysis.add_moma(model, solution=None, linear=True)[source]

Add constraints and objective representing for MOMA.

This adds variables and constraints for the minimization of metabolic adjustment (MOMA) to the model.

Parameters
  • model (cobra.Model) – The model to add MOMA constraints and objective to.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference. If no solution is given, one will be computed using pFBA.

  • linear (bool, optional) – Whether to use the linear MOMA formulation or not (default True).

Notes

In the original MOMA [1]_ specification one looks for the flux distribution of the deletion (v^d) closest to the fluxes without the deletion (v). In math this means:

minimize sum_i (v^d_i - v_i)^2 s.t. Sv^d = 0

lb_i <= v^d_i <= ub_i

Here, we use a variable transformation v^t := v^d_i - v_i. Substituting and using the fact that Sv = 0 gives:

minimize sum_i (v^t_i)^2 s.t. Sv^d = 0

v^t = v^d_i - v_i lb_i <= v^d_i <= ub_i

So basically we just re-center the flux space at the old solution and then find the flux distribution closest to the new zero (center). This is the same strategy as used in cameo.

In the case of linear MOMA [2]_, we instead minimize sum_i abs(v^t_i). The linear MOMA is typically significantly faster. Also quadratic MOMA tends to give flux distributions in which all fluxes deviate from the reference fluxes a little bit whereas linear MOMA tends to give flux distributions where the majority of fluxes are the same reference with few fluxes deviating a lot (typical effect of L2 norm vs L1 norm).

The former objective function is saved in the optlang solver interface as "moma_old_objective" and this can be used to immediately extract the value of the former objective after MOMA optimization.

See also

pfba()

parsimonious FBA

References

1

Segrè, Daniel, Dennis Vitkup, and George M. Church. “Analysis of Optimality in Natural and Perturbed Metabolic Networks.” Proceedings of the National Academy of Sciences 99, no. 23 (November 12, 2002): 15112. https://doi.org/10.1073/pnas.232349399.

2

Becker, Scott A, Adam M Feist, Monica L Mo, Gregory Hannum, Bernhard Ø Palsson, and Markus J Herrgard. “Quantitative Prediction of Cellular Metabolism with Constraint-Based Models: The COBRA Toolbox.” Nature Protocols 2 (March 29, 2007): 727.

cobra.flux_analysis.moma(model, solution=None, linear=True)[source]

Compute a single solution based on (linear) MOMA.

Compute a new flux distribution that is at a minimal distance to a previous reference solution. Minimization of metabolic adjustment (MOMA) is generally used to assess the impact of knock-outs. Thus the typical usage is to provide a wildtype flux distribution as reference and a model in knock-out state.

Parameters
  • model (cobra.Model) – The model state to compute a MOMA-based solution for.

  • solution (cobra.Solution, optional) – A (wildtype) reference solution.

  • linear (bool, optional) – Whether to use the linear MOMA formulation or not (default True).

Returns

A flux distribution that is at a minimal distance compared to the reference solution.

Return type

cobra.Solution

See also

add_moma()

add MOMA constraints and objective

cobra.flux_analysis.pfba(model, fraction_of_optimum=1.0, objective=None, reactions=None)[source]

Perform basic pFBA (parsimonious Enzyme Usage Flux Balance Analysis) to minimize total flux.

pFBA [1] adds the minimization of all fluxes the the objective of the model. This approach is motivated by the idea that high fluxes have a higher enzyme turn-over and that since producing enzymes is costly, the cell will try to minimize overall flux while still maximizing the original objective function, e.g. the growth rate.

Parameters
  • model (cobra.Model) – The model

  • fraction_of_optimum (float, optional) – Fraction of optimum which must be maintained. The original objective reaction is constrained to be greater than maximal_value * fraction_of_optimum.

  • objective (dict or model.problem.Objective) – A desired objective to use during optimization in addition to the pFBA objective. Dictionaries (reaction as key, coefficient as value) can be used for linear objectives.

  • reactions (iterable) – List of reactions or reaction identifiers. Implies return_frame to be true. Only return fluxes for the given reactions. Faster than fetching all fluxes if only a few are needed.

Returns

The solution object to the optimized model with pFBA constraints added.

Return type

cobra.Solution

References

1

Lewis, N. E., Hixson, K. K., Conrad, T. M., Lerman, J. A., Charusanti, P., Polpitiya, A. D., Palsson, B. O. (2010). Omic data from evolved E. coli are consistent with computed optimal growth from genome-scale models. Molecular Systems Biology, 6, 390. doi:10.1038/msb.2010.47

cobra.flux_analysis.find_blocked_reactions(model, reaction_list=None, zero_cutoff=None, open_exchanges=False, processes=None)[source]

Find reactions that cannot carry any flux.

The question whether or not a reaction is blocked is highly dependent on the current exchange reaction settings for a COBRA model. Hence an argument is provided to open all exchange reactions.

Notes

Sink and demand reactions are left untouched. Please modify them manually.

Parameters
  • model (cobra.Model) – The model to analyze.

  • reaction_list (list, optional) – List of reactions to consider, the default includes all model reactions.

  • zero_cutoff (float, optional) – Flux value which is considered to effectively be zero (default model.tolerance).

  • open_exchanges (bool, optional) – Whether or not to open all exchange reactions to very high flux ranges.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of reactions is large. If not explicitly passed, it will be set from the global configuration singleton.

Returns

List with the identifiers of blocked reactions.

Return type

list

cobra.flux_analysis.find_essential_genes(model, threshold=None, processes=None)[source]

Return a set of essential genes.

A gene is considered essential if restricting the flux of all reactions that depend on it to zero causes the objective, e.g., the growth rate, to also be zero, below the threshold, or infeasible.

Parameters
  • model (cobra.Model) – The model to find the essential genes for.

  • threshold (float, optional) – Minimal objective flux to be considered viable. By default this is 1% of the maximal objective.

  • processes (int, optional) – The number of parallel processes to run. If not passed, will be set to the number of CPUs found.

  • processes – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not explicitly passed, it will be set from the global configuration singleton.

Returns

Set of essential genes

Return type

set

cobra.flux_analysis.find_essential_reactions(model, threshold=None, processes=None)[source]

Return a set of essential reactions.

A reaction is considered essential if restricting its flux to zero causes the objective, e.g., the growth rate, to also be zero, below the threshold, or infeasible.

Parameters
  • model (cobra.Model) – The model to find the essential reactions for.

  • threshold (float, optional) – Minimal objective flux to be considered viable. By default this is 1% of the maximal objective.

  • processes (int, optional) – The number of parallel processes to run. Can speed up the computations if the number of knockouts to perform is large. If not explicitly passed, it will be set from the global configuration singleton.

Returns

Set of essential reactions

Return type

set

cobra.flux_analysis.flux_variability_analysis(model, reaction_list=None, loopless=False, fraction_of_optimum=1.0, pfba_factor=None, processes=None)[source]

Determine the minimum and maximum possible flux value for each reaction.

Parameters
  • model (cobra.Model) – The model for which to run the analysis. It will not be modified.

  • reaction_list (list of cobra.Reaction or str, optional) – The reactions for which to obtain min/max fluxes. If None will use all reactions in the model (default).

  • loopless (boolean, optional) – Whether to return only loopless solutions. This is significantly slower. Please also refer to the notes.

  • fraction_of_optimum (float, optional) – Must be <= 1.0. Requires that the objective value is at least the fraction times maximum objective value. A value of 0.85 for instance means that the objective has to be at least at 85% percent of its maximum.

  • pfba_factor (float, optional) – Add an additional constraint to the model that requires the total sum of absolute fluxes must not be larger than this value times the smallest possible sum of absolute fluxes, i.e., by setting the value to 1.1 the total sum of absolute fluxes must not be more than 10% larger than the pFBA solution. Since the pFBA solution is the one that optimally minimizes the total flux sum, the pfba_factor should, if set, be larger than one. Setting this value may lead to more realistic predictions of the effective flux bounds.

  • processes (int, optional) – The number of parallel processes to run. If not explicitly passed, will be set from the global configuration singleton.

Returns

A data frame with reaction identifiers as the index and two columns: - maximum: indicating the highest possible flux - minimum: indicating the lowest possible flux

Return type

pandas.DataFrame

Notes

This implements the fast version as described in [1]_. Please note that the flux distribution containing all minimal/maximal fluxes does not have to be a feasible solution for the model. Fluxes are minimized/maximized individually and a single minimal flux might require all others to be suboptimal.

Using the loopless option will lead to a significant increase in computation time (about a factor of 100 for large models). However, the algorithm used here (see [2]_) is still more than 1000x faster than the “naive” version using add_loopless(model). Also note that if you have included constraints that force a loop (for instance by setting all fluxes in a loop to be non-zero) this loop will be included in the solution.

References

1

Computationally efficient flux variability analysis. Gudmundsson S, Thiele I. BMC Bioinformatics. 2010 Sep 29;11:489. doi: 10.1186/1471-2105-11-489, PMID: 20920235

2

CycleFreeFlux: efficient removal of thermodynamically infeasible loops from flux distributions. Desouki AA, Jarre F, Gelius-Dietrich G, Lercher MJ. Bioinformatics. 2015 Jul 1;31(13):2159-65. doi: 10.1093/bioinformatics/btv096.

cobra.flux_analysis.production_envelope(model, reactions, objective=None, carbon_sources=None, points=20, threshold=None)[source]

Calculate the objective value conditioned on all combinations of fluxes for a set of chosen reactions

The production envelope can be used to analyze a model’s ability to produce a given compound conditional on the fluxes for another set of reactions, such as the uptake rates. The model is alternately optimized with respect to minimizing and maximizing the objective and the obtained fluxes are recorded. Ranges to compute production is set to the effective bounds, i.e., the minimum / maximum fluxes that can be obtained given current reaction bounds.

Parameters
  • model (cobra.Model) – The model to compute the production envelope for.

  • reactions (list or string) – A list of reactions, reaction identifiers or a single reaction.

  • objective (string, dict, model.solver.interface.Objective, optional) – The objective (reaction) to use for the production envelope. Use the model’s current objective if left missing.

  • carbon_sources (list or string, optional) – One or more reactions or reaction identifiers that are the source of carbon for computing carbon (mol carbon in output over mol carbon in input) and mass yield (gram product over gram output). Only objectives with a carbon containing input and output metabolite is supported. Will identify active carbon sources in the medium if none are specified.

  • points (int, optional) – The number of points to calculate production for.

  • threshold (float, optional) – A cut-off under which flux values will be considered to be zero (default model.tolerance).

Returns

A data frame with one row per evaluated point and

  • reaction id : one column per input reaction indicating the flux at each given point,

  • carbon_source: identifiers of carbon exchange reactions

A column for the maximum and minimum each for the following types:

  • flux: the objective flux

  • carbon_yield: if carbon source is defined and the product is a single metabolite (mol carbon product per mol carbon feeding source)

  • mass_yield: if carbon source is defined and the product is a single metabolite (gram product per 1 g of feeding source)

Return type

pandas.DataFrame

Examples

>>> import cobra.test
>>> from cobra.flux_analysis import production_envelope
>>> model = cobra.test.create_test_model("textbook")
>>> production_envelope(model, ["EX_glc__D_e", "EX_o2_e"])
cobra.flux_analysis.add_room(model, solution=None, linear=False, delta=0.03, epsilon=0.001)[source]

Add constraints and objective for ROOM.

This function adds variables and constraints for applying regulatory on/off minimization (ROOM) to the model.

Parameters
  • model (cobra.Model) – The model to add ROOM constraints and objective to.

  • solution (cobra.Solution, optional) – A previous solution to use as a reference. If no solution is given, one will be computed using pFBA.

  • linear (bool, optional) – Whether to use the linear ROOM formulation or not (default False).

  • delta (float, optional) – The relative tolerance range which is additive in nature (default 0.03).

  • epsilon (float, optional) – The absolute range of tolerance which is multiplicative (default 0.001).

Notes

The formulation used here is the same as stated in the original paper [1]_. The mathematical expression is given below:

minimize sum_{i=1}^m y^i s.t. Sv = 0

v_min <= v <= v_max v_j = 0 j ∈ A for 1 <= i <= m v_i - y_i(v_{max,i} - w_i^u) <= w_i^u (1) v_i - y_i(v_{min,i} - w_i^l) <= w_i^l (2) y_i ∈ {0,1} (3) w_i^u = w_i + delta|w_i| + epsilon w_i^l = w_i - delta|w_i| - epsilon

So, for the linear version of the ROOM , constraint (3) is relaxed to 0 <= y_i <= 1.

See also

pfba()

parsimonious FBA

References

1

Tomer Shlomi, Omer Berkman and Eytan Ruppin, “Regulatory on/off minimization of metabolic flux changes after genetic perturbations”, PNAS 2005 102 (21) 7695-7700; doi:10.1073/pnas.0406346102

cobra.flux_analysis.room(model, solution=None, linear=False, delta=0.03, epsilon=0.001)[source]

Compute a single solution based on regulatory on/off minimization (ROOM).

Compute a new flux distribution that minimizes the number of active reactions needed to accommodate a previous reference solution. Regulatory on/off minimization (ROOM) is generally used to assess the impact of knock-outs. Thus the typical usage is to provide a wildtype flux distribution as reference and a model in knock-out state.

Parameters
  • model (cobra.Model) – The model state to compute a ROOM-based solution for.

  • solution (cobra.Solution, optional) – A (wildtype) reference solution.

  • linear (bool, optional) – Whether to use the linear ROOM formulation or not (default False).

  • delta (float, optional) – The relative tolerance range (additive) (default 0.03).

  • epsilon (float, optional) – The absolute tolerance range (multiplicative) (default 0.001).

Returns

A flux distribution with minimal active reaction changes compared to the reference.

Return type

cobra.Solution

See also

add_room()

add ROOM constraints and objective