Skip to content

Commit

Permalink
refactor: allow an external solution for MOMA (#671)
Browse files Browse the repository at this point in the history
  • Loading branch information
@Midnighter
Midnighter committed May 17, 2018
1 parent 1c37dc0 commit 164a0ab
Show file tree
Hide file tree
Showing 3 changed files with 126 additions and 63 deletions.
51 changes: 22 additions & 29 deletions cobra/flux_analysis/deletion.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,25 +72,23 @@ def _init_worker(model):


def _multi_deletion(model, entity, element_lists, method="fba",
processes=None):
solution=None, processes=None):
"""
Provide a common interface for single or multiple knockouts.
Parameters
----------
model : cobra.Model
The metabolic model to perform deletions in.
entity : 'gene' or 'reaction'
The entity to knockout (``cobra.Gene`` or ``cobra.Reaction``).
element_lists : list
List of iterables ``cobra.Reaction``s or ``cobra.Gene``s (or their IDs)
to be deleted.
method: {"fba", "moma", "linear moma"}, optional
Method used to predict the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -110,7 +108,7 @@ def _multi_deletion(model, entity, element_lists, method="fba",
The solution's status.
"""
solver = sutil.interface_to_str(model.problem.__name__)
if "moma" in method and solver not in sutil.qp_solvers:
if method == "moma" and solver not in sutil.qp_solvers:
raise RuntimeError(
"Cannot use MOMA since '{}' is not QP-capable."
"Please choose a different solver or use FBA only.".format(solver))
Expand All @@ -124,7 +122,7 @@ def _multi_deletion(model, entity, element_lists, method="fba",

with model:
if "moma" in method:
add_moma(model, linear="linear" in method)
add_moma(model, solution=solution, linear="linear" in method)

args = set([frozenset(comb) for comb in product(*element_lists)])
processes = min(processes, len(args))
Expand Down Expand Up @@ -180,22 +178,21 @@ def _element_lists(entities, *ids):


def single_reaction_deletion(model, reaction_list=None, method="fba",
processes=None):
solution=None, processes=None):
"""
Knock out each reaction from a given list.
Parameters
----------
model : cobra.Model
The metabolic model to perform deletions in.
reaction_list : iterable
``cobra.Reaction``s to be deleted. If not passed,
all the reactions from the model are used.
method: {"fba", "moma", "linear moma"}, optional
Method used to predict the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -218,25 +215,25 @@ def single_reaction_deletion(model, reaction_list=None, method="fba",
return _multi_deletion(
model, 'reaction',
element_lists=_element_lists(model.reactions, reaction_list),
method=method, processes=processes)
method=method, solution=solution, processes=processes)


def single_gene_deletion(model, gene_list=None, method="fba", processes=None):
def single_gene_deletion(model, gene_list=None, method="fba", solution=None,
processes=None):
"""
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"}, optional
Method used to predict the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -258,11 +255,11 @@ def single_gene_deletion(model, gene_list=None, method="fba", processes=None):
"""
return _multi_deletion(
model, 'gene', element_lists=_element_lists(model.genes, gene_list),
method=method, processes=processes)
method=method, solution=solution, processes=processes)


def double_reaction_deletion(model, reaction_list1=None, reaction_list2=None,
method="fba", processes=None):
method="fba", solution=None, processes=None):
"""
Knock out each reaction pair from the combinations of two given lists.
Expand All @@ -272,18 +269,16 @@ def double_reaction_deletion(model, reaction_list1=None, reaction_list2=None,
----------
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"}, optional
Method used to predict the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -309,11 +304,11 @@ def double_reaction_deletion(model, reaction_list1=None, reaction_list2=None,
reaction_list2)
return _multi_deletion(
model, 'reaction', element_lists=[reaction_list1, reaction_list2],
method=method, processes=processes)
method=method, solution=solution, processes=processes)


def double_gene_deletion(model, gene_list1=None, gene_list2=None,
method="fba", processes=None):
method="fba", solution=None, processes=None):
"""
Knock out each gene pair from the combination of two given lists.
Expand All @@ -323,18 +318,16 @@ def double_gene_deletion(model, gene_list1=None, gene_list2=None,
----------
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"}, optional
Method used to predict the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -359,4 +352,4 @@ def double_gene_deletion(model, gene_list1=None, gene_list2=None,
gene_list2)
return _multi_deletion(
model, 'gene', element_lists=[gene_list1, gene_list2],
method=method, processes=processes)
method=method, solution=solution, processes=processes)
18 changes: 12 additions & 6 deletions cobra/flux_analysis/variability.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -146,7 +146,7 @@ def flux_variability_analysis(model, reaction_list=None, loopless=False,
model.solver.objective.set_linear_coefficients(
{rxn.forward_variable: 0, rxn.reverse_variable: 0})

return fva_results
return fva_results[["minimum", "maximum"]]


def find_blocked_reactions(model, reaction_list=None,
Expand Down Expand Up @@ -193,7 +193,7 @@ def find_blocked_reactions(model, reaction_list=None,
max(abs(min_max)) < zero_cutoff]


def find_essential_genes(model, threshold=None, processes=None):
def find_essential_genes(model, threshold=None, solution=None, processes=None):
"""Return a set of essential genes.
A gene is considered essential if restricting the flux of all reactions
Expand All @@ -207,6 +207,8 @@ def find_essential_genes(model, threshold=None, processes=None):
threshold : float, optional
Minimal objective flux to be considered viable. By default this is
0.01 times the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -219,13 +221,15 @@ def find_essential_genes(model, threshold=None, processes=None):
"""
if threshold is None:
threshold = model.slim_optimize(error_value=None) * 1E-02
deletions = single_gene_deletion(model, method='fba', processes=processes)
deletions = single_gene_deletion(model, method='fba', solution=solution,
processes=processes)
essential = deletions.loc[deletions['growth'].isna() |
(deletions['growth'] < threshold), :].index
return set(model.genes.get_by_id(g) for ids in essential for g in ids)


def find_essential_reactions(model, threshold=None, processes=None):
def find_essential_reactions(model, threshold=None, solution=None,
processes=None):
"""Return a set of essential reactions.
A reaction is considered essential if restricting its flux to zero
Expand All @@ -238,6 +242,8 @@ def find_essential_reactions(model, threshold=None, processes=None):
threshold : float, optional
Minimal objective flux to be considered viable. By default this is
0.01 times the growth rate.
solution : cobra.Solution
A previous solution to use as a reference for (linear) MOMA.
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,
Expand All @@ -250,8 +256,8 @@ def find_essential_reactions(model, threshold=None, processes=None):
"""
if threshold is None:
threshold = model.slim_optimize(error_value=None) * 1E-02
deletions = single_reaction_deletion(model, method='fba',
processes=processes)
deletions = single_reaction_deletion(
model, method='fba', solution=solution, processes=processes)
essential = deletions.loc[deletions['growth'].isna() |
(deletions['growth'] < threshold), :].index
return set(model.reactions.get_by_id(r) for ids in essential for r in ids)

0 comments on commit 164a0ab

Please sign in to comment.