[med-svn] [python-cobra] 01/03: New upstream version 0.9.0
Afif Elghraoui
afif at moszumanska.debian.org
Wed Oct 18 02:23:12 UTC 2017
This is an automated email from the git hooks/post-receive script.
afif pushed a commit to branch master
in repository python-cobra.
commit 544431853b10e6b06b098523786ac3f681a9df9d
Author: Afif Elghraoui <afif at debian.org>
Date: Tue Oct 17 22:11:35 2017 -0400
New upstream version 0.9.0
---
.travis.yml | 2 +-
appveyor.yml | 2 +-
cobra/__init__.py | 2 +-
cobra/core/gene.py | 4 +-
cobra/core/metabolite.py | 7 +-
cobra/core/model.py | 65 ++++--
cobra/core/reaction.py | 43 ++--
cobra/core/solution.py | 11 +-
cobra/flux_analysis/phenotype_phase_plane.py | 305 +++++++++++++++------------
cobra/flux_analysis/summary.py | 10 +-
cobra/manipulation/validate.py | 4 -
cobra/test/data/iJO1366.pickle | Bin 1810497 -> 1810487 bytes
cobra/test/data/mini.pickle | Bin 32825 -> 33312 bytes
cobra/test/data/raven.pickle | Bin 13758 -> 13310 bytes
cobra/test/data/salmonella.pickle | Bin 2151196 -> 2151186 bytes
cobra/test/data/textbook_solution.pickle | Bin 143584 -> 7330 bytes
cobra/test/test_flux_analysis.py | 22 +-
cobra/test/test_manipulation.py | 3 +-
cobra/test/test_model.py | 27 +++
cobra/util/util.py | 7 +
release-notes/0.8.2.md | 13 +-
release-notes/0.9.0.md | 30 +++
release-notes/next-release.md | 2 +
scripts/compare-benchmark.py | 16 +-
setup.cfg | 2 +-
setup.py | 2 +-
26 files changed, 365 insertions(+), 214 deletions(-)
diff --git a/.travis.yml b/.travis.yml
index db8a64a..5bb203c 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -68,7 +68,7 @@ matrix:
before_install:
- if [[ -n "${MB_PYTHON_VERSION}" ]]; then
- (travis_retry git clone https://github.com/matthew-brett/multibuild.git && cd multibuild && git checkout edf5b691d0d565b4e65e655b983c11c883acbeca);
+ (travis_retry git clone https://github.com/matthew-brett/multibuild.git && cd multibuild && git checkout 37040e31b1044468027bd86991c805006a92bf09);
TEST_DEPENDS="swiglpk optlang sympy decorator cython codecov coverage numpy scipy jsonschema six pytest pytest-cov pytest-benchmark tabulate";
BUILD_DEPENDS="swiglpk optlang sympy cython numpy scipy";
source multibuild/common_utils.sh;
diff --git a/appveyor.yml b/appveyor.yml
index 2571dc9..017a3a4 100644
--- a/appveyor.yml
+++ b/appveyor.yml
@@ -84,7 +84,7 @@ install:
python -m pip install pip setuptools>=24.0 wheel --upgrade }
- if not defined CONDA %WITH_COMPILER% python -m pip install --upgrade pytest
- if defined CONDA conda install -q setuptools pytest numpy scipy
- - python -m pip install pytest-cov pytest-benchmark pandas swiglpk optlang python-libsbml decorator Cython jsonschema twine pypandoc==1.1.3
+ - python -m pip install pytest-cov pytest-benchmark pandas swiglpk optlang python-libsbml decorator Cython jsonschema twine pypandoc
- "%WITH_COMPILER% python appveyor/build_glpk.py"
build: off
diff --git a/cobra/__init__.py b/cobra/__init__.py
index 21666c1..ce83d7a 100644
--- a/cobra/__init__.py
+++ b/cobra/__init__.py
@@ -13,7 +13,7 @@ from cobra.core import (
DictList, Gene, Metabolite, Model, Object, Reaction, Species)
from cobra.util.version_info import show_versions
-__version__ = "0.8.2"
+__version__ = "0.9.0"
# set the warning format to be prettier and fit on one line
_cobra_path = _dirname(_abspath(__file__))
diff --git a/cobra/core/gene.py b/cobra/core/gene.py
index 4ea1f1b..b869b06 100644
--- a/cobra/core/gene.py
+++ b/cobra/core/gene.py
@@ -11,6 +11,7 @@ from warnings import warn
from cobra.core.species import Species
from cobra.util import resettable
+from cobra.util.util import format_long_string
keywords = list(kwlist)
keywords.remove("and")
@@ -289,4 +290,5 @@ class Gene(Species):
functional=self.functional,
address='0x0%x' % id(self),
n_reactions=len(self.reactions),
- reactions=', '.join(r.id for r in self.reactions))
+ reactions=format_long_string(
+ ', '.join(r.id for r in self.reactions), 200))
diff --git a/cobra/core/metabolite.py b/cobra/core/metabolite.py
index 8f78e39..c3cd0f7 100644
--- a/cobra/core/metabolite.py
+++ b/cobra/core/metabolite.py
@@ -12,6 +12,7 @@ from cobra.exceptions import OptimizationError
from cobra.core.formula import elements_and_molecular_weights
from cobra.core.species import Species
from cobra.util.solver import check_solver_status
+from cobra.util.util import format_long_string
# Numbers are not required because of the |(?=[A-Z])? block. See the
@@ -257,8 +258,10 @@ class Metabolite(Species):
<td><strong>In {n_reactions} reaction(s)</strong></td><td>
{reactions}</td>
</tr>
- </table>""".format(id=self.id, name=self.name, formula=self.formula,
+ </table>""".format(id=self.id, name=format_long_string(self.name),
+ formula=self.formula,
address='0x0%x' % id(self),
compartment=self.compartment,
n_reactions=len(self.reactions),
- reactions=', '.join(r.id for r in self.reactions))
+ reactions=format_long_string(
+ ', '.join(r.id for r in self.reactions), 200))
diff --git a/cobra/core/model.py b/cobra/core/model.py
index fb0be9a..7eabf78 100644
--- a/cobra/core/model.py
+++ b/cobra/core/model.py
@@ -24,7 +24,7 @@ from cobra.util.solver import (
SolverNotFound, get_solver_name, interface_to_str, set_objective, solvers,
add_cons_vars_to_problem, remove_cons_vars_from_problem, choose_solver,
check_solver_status, assert_optimal)
-from cobra.util.util import AutoVivification
+from cobra.util.util import AutoVivification, format_long_string
LOGGER = logging.getLogger(__name__)
@@ -63,8 +63,6 @@ class Model(Object):
for y in ['reactions', 'genes', 'metabolites']:
for x in getattr(self, y):
x._model = self
- if y == 'reactions':
- x._reset_var_cache()
if not hasattr(self, "name"):
self.name = None
@@ -94,7 +92,7 @@ class Model(Object):
self.reactions = DictList() # A list of cobra.Reactions
self.metabolites = DictList() # A list of cobra.Metabolites
# genes based on their ids {Gene.id: Gene}
- self.compartments = dict()
+ self._compartments = dict()
self._contexts = []
# from cameo ...
@@ -148,9 +146,6 @@ class Model(Object):
# Do nothing if the solver did not change
if self.problem == interface:
return
-
- for reaction in self.reactions:
- reaction._reset_var_cache()
self._solver = interface.Model.clone(self._solver)
@property
@@ -165,10 +160,37 @@ class Model(Object):
def get_metabolite_compartments(self):
"""Return all metabolites' compartments."""
+ warn('use Model.compartments instead', DeprecationWarning)
return {met.compartment for met in self.metabolites
if met.compartment is not None}
@property
+ def compartments(self):
+ return {met.compartment: self._compartments.get(met.compartment, '')
+ for met in self.metabolites if met.compartment is not None}
+
+ @compartments.setter
+ def compartments(self, value):
+ """Get or set the dictionary of current compartment descriptions.
+
+ Assigning a dictionary to this property updates the model's
+ dictionary of compartment descriptions with the new values.
+
+ Parameters
+ ----------
+ value : dict
+ Dictionary mapping compartments abbreviations to full names.
+
+ Examples
+ --------
+ >>> import cobra.test
+ >>> model = cobra.test.create_test_model("textbook")
+ >>> model.compartments = {'c': 'the cytosol'}
+ {'c': 'the cytosol', 'e': 'extracellular'}
+ """
+ self._compartments.update(value)
+
+ @property
def medium(self):
def is_active(reaction):
@@ -300,9 +322,6 @@ class Model(Object):
new_gene = new.genes.get_by_id(gene.id)
new_reaction._genes.add(new_gene)
new_gene._reaction.add(new_reaction)
-
- for reaction in new.reactions:
- reaction._reset_var_cache()
try:
new._solver = deepcopy(self.solver)
# Cplex has an issue with deep copies
@@ -521,7 +540,6 @@ class Model(Object):
# Add reactions. Also take care of genes and metabolites in the loop
for reaction in reaction_list:
- reaction._reset_var_cache()
reaction._model = self # the reaction now points to the model
# keys() is necessary because the dict will be modified during
# the loop
@@ -604,7 +622,6 @@ class Model(Object):
reaction._model = None
if context:
- context(reaction._reset_var_cache)
context(partial(setattr, reaction, '_model', self))
context(partial(self.reactions.add, reaction))
@@ -940,6 +957,28 @@ class Model(Object):
value = {rxn: 1 for rxn in reactions}
set_objective(self, value, additive=False)
+ @property
+ def objective_direction(self):
+ """
+ Get or set the objective direction.
+
+ When using a `HistoryManager` context, this attribute can be set
+ temporarily, reversed when exiting the context.
+
+ """
+ return self.solver.objective.direction
+
+ @objective_direction.setter
+ @resettable
+ def objective_direction(self, value):
+ value = value.lower()
+ if value.startswith("max"):
+ self.solver.objective.direction = "max"
+ elif value.startswith("min"):
+ self.solver.objective.direction = "min"
+ else:
+ raise ValueError("Unknown objective direction '{}'.".format(value))
+
def summary(self, solution=None, threshold=1E-8, fva=None, floatfmt='.3g'):
"""Print a summary of the input and output fluxes of the model. This
method requires the model to have been previously solved.
@@ -1064,7 +1103,7 @@ class Model(Object):
address='0x0%x' % id(self),
num_metabolites=len(self.metabolites),
num_reactions=len(self.reactions),
- objective=str(self.objective.expression),
+ objective=format_long_string(str(self.objective.expression), 100),
compartments=", ".join(
v if v else k for k, v in iteritems(self.compartments)
))
diff --git a/cobra/core/reaction.py b/cobra/core/reaction.py
index e2d424f..e5c3c42 100644
--- a/cobra/core/reaction.py
+++ b/cobra/core/reaction.py
@@ -20,6 +20,7 @@ from cobra.core.object import Object
from cobra.util.context import resettable, get_context
from cobra.util.solver import (
linear_reaction_coefficients, set_objective, check_solver_status)
+from cobra.util.util import format_long_string
# precompiled regular expressions
# Matches and/or in a gene reaction rule
@@ -137,11 +138,7 @@ class Reaction(Object):
not associated with a model.
"""
if self.model is not None:
- if self._forward_variable is None:
- self._forward_variable = self.model.variables[
- self.id]
- assert self._forward_variable.problem is self.model.solver
- return self._forward_variable
+ return self.model.variables[self.id]
else:
return None
@@ -155,20 +152,12 @@ class Reaction(Object):
An optlang variable for the reverse flux or None if reaction is
not associated with a model.
"""
- model = self.model
- if model is not None:
- if self._reverse_variable is None:
- self._reverse_variable = model.variables[
- self.reverse_id]
- assert self._reverse_variable.problem is self.model.solver
- return self._reverse_variable
+
+ if self.model is not None:
+ return self.model.variables[self.reverse_id]
else:
return None
- def _reset_var_cache(self):
- self._forward_variable = None
- self._reverse_variable = None
-
@property
def objective_coefficient(self):
""" Get the coefficient for this reaction in a linear
@@ -190,12 +179,10 @@ class Reaction(Object):
def __copy__(self):
cop = copy(super(Reaction, self))
- cop._reset_var_cache()
return cop
def __deepcopy__(self, memo):
cop = deepcopy(super(Reaction, self), memo)
- cop._reset_var_cache()
return cop
@property
@@ -621,9 +608,15 @@ class Reaction(Object):
"""
new_reaction = self.copy()
- new_reaction += other
+ if other == 0:
+ return new_reaction
+ else:
+ new_reaction += other
+
return new_reaction
+ __radd__ = __add__
+
def __iadd__(self, other):
self.add_metabolites(other._metabolites, combine=True)
gpr1 = self.gene_reaction_rule.strip()
@@ -911,6 +904,7 @@ class Reaction(Object):
def get_compartments(self):
"""lists compartments the metabolites are in"""
+ warn('use Reaction.compartments instead', DeprecationWarning)
return list(self.compartments)
def _associate_gene(self, cobra_gene):
@@ -1064,11 +1058,14 @@ class Reaction(Object):
<td><strong>Upper bound</strong></td><td>{ub}</td>
</tr>
</table>
- """.format(id=self.id, name=self.name,
+ """.format(id=format_long_string(self.id, 100),
+ name=format_long_string(self.name, 100),
address='0x0%x' % id(self),
- stoich_id=self.build_reaction_string(),
- stoich_name=self.build_reaction_string(True),
- gpr=self.gene_reaction_rule,
+ stoich_id=format_long_string(
+ self.build_reaction_string(), 200),
+ stoich_name=format_long_string(
+ self.build_reaction_string(True), 200),
+ gpr=format_long_string(self.gene_reaction_rule, 100),
lb=self.lower_bound, ub=self.upper_bound)
diff --git a/cobra/core/solution.py b/cobra/core/solution.py
index e730565..05bba87 100644
--- a/cobra/core/solution.py
+++ b/cobra/core/solution.py
@@ -10,7 +10,7 @@ from warnings import warn
from numpy import empty, nan
from optlang.interface import OPTIMAL
-from pandas import Series, DataFrame
+from pandas import Series, DataFrame, option_context
from cobra.util.solver import check_solver_status
@@ -89,10 +89,11 @@ class Solution(object):
def _repr_html_(self):
if self.status == OPTIMAL:
- html = ('<strong><em>Optimal</em> solution with objective value '
- '{:.3f}</strong><br>{}'
- .format(self.objective_value,
- self.to_frame()._repr_html_()))
+ with option_context('display.max_rows', 10):
+ html = ('<strong><em>Optimal</em> solution with objective '
+ 'value {:.3f}</strong><br>{}'
+ .format(self.objective_value,
+ self.to_frame()._repr_html_()))
else:
html = '<strong><em>{}</em> solution</strong>'.format(self.status)
return html
diff --git a/cobra/flux_analysis/phenotype_phase_plane.py b/cobra/flux_analysis/phenotype_phase_plane.py
index 72d36ac..5eacbd2 100644
--- a/cobra/flux_analysis/phenotype_phase_plane.py
+++ b/cobra/flux_analysis/phenotype_phase_plane.py
@@ -3,17 +3,15 @@
from __future__ import absolute_import, division
from warnings import warn
-from collections import defaultdict
-from operator import itemgetter
-from itertools import product
from six import iteritems
+from itertools import product
from multiprocessing import Pool
+
import pandas as pd
from optlang.interface import OPTIMAL
-
from numpy import (
nan, abs, arange, dtype, empty, int32, linspace, meshgrid, unravel_index,
- zeros, array)
+ zeros, full)
from cobra.solvers import get_solver_name, solver_dict
import cobra.util.solver as sutil
@@ -307,17 +305,18 @@ def calculate_phenotype_phase_plane(
return data
-def production_envelope(model, reactions, objective=None, c_source=None,
- points=20, solver=None):
+def production_envelope(model, reactions, objective=None, carbon_sources=None,
+ points=20, threshold=1e-7, solver=None):
"""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 models ability to
+ 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
- reaction, such as the uptake rates. The model is alternately optimize
- with respect to minimizing and maximizing the objective and record the
- obtained fluxes. Ranges to compute production is set to the effective
- bounds, i.e. the minimum / maximum fluxes that can be obtained given
+ 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
@@ -325,21 +324,24 @@ def production_envelope(model, reactions, objective=None, c_source=None,
model : cobra.Model
The model to compute the production envelope for.
reactions : list or string
- A list of reactions, reaction identifiers or single reaction
- objective : string, dict, model.solver.interface.Objective
+ 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 is left missing.
- c_source : cobra.Reaction or string
- A reaction or reaction identifier that is 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.
- points : int
+ 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.
- solver : string
+ threshold : float, optional
+ A cut-off under which flux values will be considered to be zero.
+ solver : string, optional
The solver to use - only here for consistency with older
implementations (this argument will be removed in the future). The
- solver should be set using `model.solver` directly. Only optlang
+ solver should be set using ``model.solver`` directly. Only optlang
based solvers are supported.
Returns
@@ -349,18 +351,16 @@ def production_envelope(model, reactions, objective=None, c_source=None,
- reaction id : one column per input reaction indicating the flux at
each given point,
- - flux: the objective flux
+ - 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)
- - direction: the direction of the optimization.
-
- Only points that give a valid solution are returned.
-
Examples
--------
>>> import cobra.test
@@ -374,168 +374,193 @@ def production_envelope(model, reactions, objective=None, c_source=None,
'based solver interfaces.')
reactions = model.reactions.get_by_any(reactions)
objective = model.solver.objective if objective is None else objective
- result = None
+ data = dict()
+ if carbon_sources is None:
+ c_input = find_carbon_sources(model)
+ else:
+ c_input = model.reactions.get_by_any(carbon_sources)
+ if c_input is None:
+ data['carbon_source'] = None
+ elif hasattr(c_input, 'id'):
+ data['carbon_source'] = c_input.id
+ else:
+ data['carbon_source'] = ', '.join(rxn.id for rxn in c_input)
+
+ size = points ** len(reactions)
+ for direction in ('minimum', 'maximum'):
+ data['flux_{}'.format(direction)] = full(size, nan, dtype=float)
+ data['carbon_yield_{}'.format(direction)] = full(
+ size, nan, dtype=float)
+ data['mass_yield_{}'.format(direction)] = full(
+ size, nan, dtype=float)
+ grid = pd.DataFrame(data)
with model:
model.objective = objective
- if c_source is None:
- c_input = get_c_input(model)
- else:
- c_input = model.reactions.get_by_any(c_source)[0]
objective_reactions = list(sutil.linear_reaction_coefficients(model))
if len(objective_reactions) != 1:
raise ValueError('cannot calculate yields for objectives with '
'multiple reactions')
- carbon_io = c_input, objective_reactions[0]
+ c_output = objective_reactions[0]
min_max = fva(model, reactions, fraction_of_optimum=0)
- grid = [linspace(min_max.minimum[rxn.id], min_max.maximum[rxn.id],
- points, endpoint=True) for rxn in reactions]
- grid_list = list(product(*grid))
- result = envelope_for_points(model, reactions, grid_list, carbon_io)
-
- return pd.DataFrame(result)
-
+ min_max[min_max.abs() < threshold] = 0.0
+ points = list(product(*[
+ linspace(min_max.at[rxn.id, "minimum"],
+ min_max.at[rxn.id, "maximum"],
+ points, endpoint=True) for rxn in reactions]))
+ tmp = pd.DataFrame(points, columns=[rxn.id for rxn in reactions])
+ grid = pd.concat([grid, tmp], axis=1, copy=False)
+ add_envelope(model, reactions, grid, c_input, c_output, threshold)
+ return grid
+
+
+def add_envelope(model, reactions, grid, c_input, c_output, threshold):
+ if c_input is not None:
+ input_components = [reaction_elements(rxn) for rxn in c_input]
+ output_components = reaction_elements(c_output)
+ try:
+ input_weights = [reaction_weight(rxn) for rxn in c_input]
+ output_weight = reaction_weight(c_output)
+ except ValueError:
+ input_weights = []
+ output_weight = []
+ else:
+ input_components = []
+ output_components = []
+ input_weights = []
+ output_weight = []
-def envelope_for_points(model, reactions, grid, carbon_io):
- results = defaultdict(list)
for direction in ('minimum', 'maximum'):
- sense = "min" if direction == "minimum" else "max"
- for point in grid:
- with model:
- model.solver.objective.direction = sense
- for reaction, coordinate in zip(reactions, point):
- reaction.bounds = coordinate, coordinate
- model.slim_optimize()
- if model.solver.status == OPTIMAL:
- for reaction, coordinate in zip(reactions, point):
- results[reaction.id].append(coordinate)
- results['direction'].append(direction)
- results['flux'].append(model.solver.objective.value)
- if carbon_io[0] is not None:
- results['carbon_yield'].append(carbon_yield(carbon_io))
- results['mass_yield'].append(mass_yield(carbon_io))
- for key, value in results.items():
- results[key] = array(value)
- if carbon_io[0] is not None:
- results['carbon_source'] = carbon_io[0].id
- return results
+ with model:
+ model.objective_direction = direction
+ for i in range(len(grid)):
+ with model:
+ for rxn in reactions:
+ point = grid.at[i, rxn.id]
+ rxn.bounds = point, point
+ obj_val = model.slim_optimize()
+ if model.solver.status != OPTIMAL:
+ continue
+ grid.at[i, 'flux_{}'.format(direction)] = \
+ 0.0 if abs(obj_val) < threshold else obj_val
+ if c_input is not None:
+ grid.at[i, 'carbon_yield_{}'.format(direction)] = \
+ total_yield([rxn.flux for rxn in c_input],
+ input_components,
+ obj_val,
+ output_components)
+ grid.at[i, 'mass_yield_{}'.format(direction)] = \
+ total_yield([rxn.flux for rxn in c_input],
+ input_weights,
+ obj_val,
+ output_weight)
+
+
+def total_yield(input_fluxes, input_elements, output_flux, output_elements):
+ """
+ Compute total output per input unit.
+ Units are typically mol carbon atoms or gram of source and product.
-def carbon_yield(c_input_output):
- """ mol product per mol carbon input
+ Parameters
+ ----------
+ input_fluxes : list
+ A list of input reaction fluxes in the same order as the
+ ``input_components``.
+ input_elements : list
+ A list of reaction components which are in turn list of numbers.
+ output_flux : float
+ The output flux value.
+ output_elements : list
+ A list of stoichiometrically weighted output reaction components.
Returns
-------
float
- the mol carbon atoms in the product (as defined by the model
- objective) divided by the mol carbon in the input reactions (as
- defined by the model medium) or zero in case of division by zero
- arises
+ The ratio between output (mol carbon atoms or grams of product) and
+ input (mol carbon atoms or grams of source compounds).
"""
- c_input, c_output = c_input_output
- if c_input is None:
- return nan
- carbon_input_flux = total_carbon_flux(c_input, consumption=True)
- carbon_output_flux = total_carbon_flux(c_output, consumption=False)
+ carbon_input_flux = sum(
+ total_components_flux(flux, components, consumption=True)
+ for flux, components in zip(input_fluxes, input_elements))
+ carbon_output_flux = total_components_flux(
+ output_flux, output_elements, consumption=False)
try:
return carbon_output_flux / carbon_input_flux
except ZeroDivisionError:
return nan
-def mass_yield(c_input_output):
- """Gram product divided by gram of carbon input source
+def reaction_elements(reaction):
+ """
+ Split metabolites into the atoms times their stoichiometric coefficients.
Parameters
----------
- c_input_output : tuple
- Two reactions, the one that feeds carbon to the system and the one
- that produces carbon containing compound.
+ reaction : Reaction
+ The metabolic reaction whose components are desired.
Returns
-------
- float
- gram product per 1 g of feeding source
+ list
+ Each of the reaction's metabolites' desired carbon elements (if any)
+ times that metabolite's stoichiometric coefficient.
"""
- c_input, c_output = c_input_output
- if input is None:
- return nan
- try:
- c_source, source_flux = single_flux(c_input, consumption=True)
- c_product, product_flux = single_flux(c_output, consumption=False)
- except ValueError:
- return nan
- mol_prod_mol_src = product_flux / source_flux
- x = mol_prod_mol_src * c_product.formula_weight
- return x / c_source.formula_weight
+ c_elements = [coeff * met.elements.get('C', 0)
+ for met, coeff in iteritems(reaction.metabolites)]
+ return [elem for elem in c_elements if elem != 0]
-def total_carbon_flux(reaction, consumption=True):
- """summed product carbon flux for a reaction
+def reaction_weight(reaction):
+ """Return the metabolite weight times its stoichiometric coefficient."""
+ if len(reaction.metabolites) != 1:
+ raise ValueError('Reaction weight is only defined for single '
+ 'metabolite products or educts.')
+ met, coeff = next(iteritems(reaction.metabolites))
+ return [coeff * met.formula_weight]
+
+
+def total_components_flux(flux, components, consumption=True):
+ """
+ Compute the total components consumption or production flux.
Parameters
----------
- reaction : Reaction
- the reaction to carbon return flux for
- consumption : bool
- flux for consumed metabolite, else produced
+ flux : float
+ The reaction flux for the components.
+ components : list
+ List of stoichiometrically weighted components.
+ consumption : bool, optional
+ Whether to sum up consumption or production fluxes.
- Returns
- -------
- float
- reaction flux multiplied by number of carbon for the products of the
- reaction
"""
direction = 1 if consumption else -1
- c_flux = [reaction.flux * coeff * met.elements.get('C', 0) * direction
- for met, coeff in reaction.metabolites.items()]
+ c_flux = [elem * flux * direction for elem in components]
return sum([flux for flux in c_flux if flux > 0])
-def single_flux(reaction, consumption=True):
- """flux into single product for a reaction
-
- only defined for reactions with single products
+def find_carbon_sources(model):
+ """
+ Find all active carbon source reactions.
Parameters
----------
- reaction : Reaction
- the reaction to product flux for
- consumption : bool
- flux for consumed metabolite, else produced
+ model : Model
+ A genome-scale metabolic model.
Returns
-------
- tuple
- metabolite, flux for the metabolite
- """
- if len(list(reaction.metabolites)) != 1:
- raise ValueError('product flux only defined for single metabolite '
- 'reactions')
- met, coeff = next(iteritems(reaction.metabolites))
- direction = 1 if consumption else -1
- return met, reaction.flux * coeff * direction
-
-
-def get_c_input(model):
- """ carbon source reactions
+ list
+ The medium reactions with carbon input flux.
- Returns
- -------
- Reaction
- The medium reaction with highest input carbon flux
"""
try:
model.slim_optimize(error_value=None)
except OptimizationError:
- return None
+ return []
reactions = model.reactions.get_by_any(list(model.medium))
- reactions_fluxes = [(rxn, total_carbon_flux(rxn, consumption=True))
- for rxn in reactions]
- source_reactions = [(rxn, c_flux) for rxn, c_flux
- in reactions_fluxes if c_flux > 0]
- try:
- return max(source_reactions, key=itemgetter(1))[0]
- except ValueError:
- return None
+ reactions_fluxes = [
+ (rxn, total_components_flux(rxn.flux, reaction_elements(rxn),
+ consumption=True)) for rxn in reactions]
+ return [rxn for rxn, c_flux in reactions_fluxes if c_flux > 0]
diff --git a/cobra/flux_analysis/summary.py b/cobra/flux_analysis/summary.py
index c75f4b6..a53b205 100644
--- a/cobra/flux_analysis/summary.py
+++ b/cobra/flux_analysis/summary.py
@@ -10,16 +10,10 @@ from tabulate import tabulate
from cobra.flux_analysis.variability import flux_variability_analysis
from cobra.util.solver import linear_reaction_coefficients
+from cobra.util.util import format_long_string
from cobra.core import get_solution
-def format_long_string(string, max_length):
- if len(string) > max_length:
- string = string[:max_length - 3]
- string += '...'
- return string
-
-
def metabolite_summary(met, solution=None, threshold=0.01, fva=False,
floatfmt='.3g'):
"""Print a summary of the reactions which produce and consume this
@@ -145,7 +139,7 @@ def model_summary(model, solution=None, threshold=1E-8, fva=None,
"""
objective_reactions = linear_reaction_coefficients(model)
boundary_reactions = model.exchanges
- summary_rxns = list(objective_reactions.keys()) + boundary_reactions
+ summary_rxns = set(objective_reactions.keys()).union(boundary_reactions)
if solution is None:
model.slim_optimize(error_value=None)
diff --git a/cobra/manipulation/validate.py b/cobra/manipulation/validate.py
index cf39969..e9f217f 100644
--- a/cobra/manipulation/validate.py
+++ b/cobra/manipulation/validate.py
@@ -50,10 +50,6 @@ def check_reaction_bounds(model):
def check_metabolite_compartment_formula(model):
errors = []
for met in model.metabolites:
- if met.compartment is not None and \
- met.compartment not in model.compartments:
- errors.append("Metabolite '%s' compartment '%s' not found" %
- (met.id, met.compartment))
if met.formula is not None and len(met.formula) > 0:
if not met.formula.isalnum():
errors.append("Metabolite '%s' formula '%s' not alphanumeric" %
diff --git a/cobra/test/data/iJO1366.pickle b/cobra/test/data/iJO1366.pickle
index 6346508..2e89067 100644
Binary files a/cobra/test/data/iJO1366.pickle and b/cobra/test/data/iJO1366.pickle differ
diff --git a/cobra/test/data/mini.pickle b/cobra/test/data/mini.pickle
index eb9d561..af9ef9c 100644
Binary files a/cobra/test/data/mini.pickle and b/cobra/test/data/mini.pickle differ
diff --git a/cobra/test/data/raven.pickle b/cobra/test/data/raven.pickle
index 9750952..adfaf9b 100644
Binary files a/cobra/test/data/raven.pickle and b/cobra/test/data/raven.pickle differ
diff --git a/cobra/test/data/salmonella.pickle b/cobra/test/data/salmonella.pickle
index 3848bdb..ba0568c 100644
Binary files a/cobra/test/data/salmonella.pickle and b/cobra/test/data/salmonella.pickle differ
diff --git a/cobra/test/data/textbook_solution.pickle b/cobra/test/data/textbook_solution.pickle
index f5e453d..a9633a3 100644
Binary files a/cobra/test/data/textbook_solution.pickle and b/cobra/test/data/textbook_solution.pickle differ
diff --git a/cobra/test/test_flux_analysis.py b/cobra/test/test_flux_analysis.py
index 773d10b..aef830e 100644
--- a/cobra/test/test_flux_analysis.py
+++ b/cobra/test/test_flux_analysis.py
@@ -620,6 +620,10 @@ class TestCobraFluxAnalysis:
model.summary()
self.check_in_line(out.getvalue(), expected_entries)
+ with model:
+ model.objective = model.exchanges[0]
+ model.summary()
+
@pytest.mark.parametrize("fraction", [0.95])
def test_model_summary_with_fva(self, model, opt_solver, fraction):
if opt_solver == "optlang-gurobi":
@@ -834,7 +838,7 @@ class TestProductionEnvelope:
def test_envelope_one(self, model):
df = production_envelope(model, ["EX_o2_e"])
- assert abs(sum(df.flux) - 9.342) < 0.001
+ assert numpy.isclose(df["flux_maximum"].sum(), 9.342, atol=1e-3)
def test_envelope_multi_reaction_objective(self, model):
obj = {model.reactions.EX_ac_e: 1,
@@ -842,12 +846,22 @@ class TestProductionEnvelope:
with pytest.raises(ValueError):
production_envelope(model, "EX_o2_e", obj)
+ @pytest.mark.parametrize("variables, num", [
+ (["EX_glc__D_e"], 30),
+ (["EX_glc__D_e", "EX_o2_e"], 20),
+ (["EX_glc__D_e", "EX_o2_e", "EX_ac_e"], 10)
+ ])
+ def test_multi_variable_envelope(self, model, variables, num):
+ df = production_envelope(model, variables, points=num)
+ assert len(df) == num ** len(variables)
+
def test_envelope_two(self, model):
df = production_envelope(model, ["EX_glc__D_e", "EX_o2_e"],
objective="EX_ac_e")
- assert abs(numpy.sum(df.carbon_yield) - 83.579) < 0.001
- assert abs(numpy.sum(df.flux) - 1737.466) < 0.001
- assert abs(numpy.sum(df.mass_yield) - 82.176) < 0.001
+ assert numpy.isclose(df["flux_maximum"].sum(), 1737.466, atol=1e-3)
+ assert numpy.isclose(df["carbon_yield_maximum"].sum(), 83.579,
+ atol=1e-3)
+ assert numpy.isclose(df["mass_yield_maximum"].sum(), 82.176, atol=1e-3)
class TestReactionUtils:
diff --git a/cobra/test/test_manipulation.py b/cobra/test/test_manipulation.py
index 40b98d5..a7d37f3 100644
--- a/cobra/test/test_manipulation.py
+++ b/cobra/test/test_manipulation.py
@@ -225,10 +225,9 @@ class TestManipulation:
assert rxns.EX_h_e.annotation["SBO"] == "SBO:0000628"
def test_validate_formula_compartment(self, model):
- model.metabolites[1].compartment = "fake"
model.metabolites[1].formula = "(a*.bcde)"
errors = check_metabolite_compartment_formula(model)
- assert len(errors) == 2
+ assert len(errors) == 1
def test_validate_mass_balance(self, model):
assert len(check_mass_balance(model)) == 0
diff --git a/cobra/test/test_model.py b/cobra/test/test_model.py
index 745d5d7..b272264 100644
--- a/cobra/test/test_model.py
+++ b/cobra/test/test_model.py
@@ -287,6 +287,11 @@ class TestReactions:
assert gene is not model.genes.get_by_id(gene.id)
assert gene.model is not model
+ def test_radd(self, model):
+ new = sum([model.reactions.PGI, model.reactions.EX_h2o_e])
+ assert new._model is not model
+ assert len(new.metabolites) == 3
+
def test_mul(self, model):
new = model.reactions.PGI * 2
assert set(new.metabolites.values()) == {-2, 2}
@@ -395,6 +400,15 @@ class TestCobraModel:
def test_compartments(self, model):
assert set(model.compartments) == {"c", "e"}
+ model = Model("test", "test")
+ met_c = Metabolite("a_c", compartment="c")
+ met_e = Metabolite("a_e", compartment="e")
+ rxn = Reaction("foo")
+ rxn.add_metabolites({met_e: -1, met_c: 1})
+ model.add_reactions([rxn])
+ assert model.compartments == {'c': '', 'e': ''}
+ model.compartments = {'c': 'cytosol'}
+ assert model.compartments == {'c': 'cytosol', 'e': ''}
def test_add_reaction(self, model):
old_reaction_count = len(model.reactions)
@@ -766,6 +780,19 @@ class TestCobraModel:
benchmark(benchmark_change_objective)
+ def test_get_objective_direction(self, model):
+ assert model.objective_direction == "max"
+ value = model.slim_optimize()
+ assert numpy.isclose(value, 0.874, 1e-3)
+
+ def test_set_objective_direction(self, model):
+ with model:
+ model.objective_direction = "min"
+ assert model.objective_direction == "min"
+ value = model.slim_optimize()
+ assert value == 0.0
+ assert model.objective_direction == "max"
+
def test_slim_optimize(self, model):
with model:
assert model.slim_optimize() > 0.872
diff --git a/cobra/util/util.py b/cobra/util/util.py
index 3aeb377..6299a50 100644
--- a/cobra/util/util.py
+++ b/cobra/util/util.py
@@ -3,6 +3,13 @@
from __future__ import absolute_import
+def format_long_string(string, max_length=50):
+ if len(string) > max_length:
+ string = string[:max_length - 3]
+ string += '...'
+ return string
+
+
class AutoVivification(dict):
"""Implementation of perl's autovivification feature. Checkout
http://stackoverflow.com/a/652284/280182 """
diff --git a/release-notes/0.8.2.md b/release-notes/0.8.2.md
index d2ecece..a5251a2 100644
--- a/release-notes/0.8.2.md
+++ b/release-notes/0.8.2.md
@@ -2,8 +2,6 @@
## Fixes
-- the Solution class no longer contains links progenitor model's
- reactions and metabolites
- Guarantee that sampler._reproject always returns a feasible point
and will not attempt to reproject already feasible
points. [#564](https://github.com/opencobra/cobrapy/pull/564)
@@ -15,3 +13,14 @@
model or `ValueError` is raised.
- Fix use of underscores in key/value pairs in legacy sbml
notes. [#547](https://github.com/opencobra/cobrapy/issues/547)
+
+## Backwards incompatible changes
+
+- the Solution class no longer contains links progenitor model's
+ reactions and metabolites. Removed since it those can change after
+ the solution has been computed making them erroneous. This of course
+ implies that `Solution` constructor changes to:
+ ```
+ def __init__(self, objective_value, status, fluxes, reduced_costs=None,
+ shadow_prices=None, **kwargs):
+ ```
diff --git a/release-notes/0.9.0.md b/release-notes/0.9.0.md
new file mode 100644
index 0000000..178a215
--- /dev/null
+++ b/release-notes/0.9.0.md
@@ -0,0 +1,30 @@
+# Release notes for cobrapy 0.9.0
+
+## Fixes
+
+- `Model.compartment` is now a dynamic property fetching the
+ compartments from all metabolites therefore always
+ up-to-date. Assigning a dictionary to the same property updates the
+ internal dictionary of compartment descriptions. This change removes
+ the need for the check for missing compartments from
+ `validation.check_metabolite_compartment_formula`.
+- Excessively long output of html representations in jupyter notebooks are now abbreviated [#577](https://github.com/opencobra/cobrapy/pull/577).
+- Reaction forward and reverse variables are no longer cached with those object. No visible effect but simplifies the code.
+- Fix bug in summary methods when used with exchange reaction sas objective. [#595](https://github.com/opencobra/cobrapy/pull/595).
+
+## New features
+
+- `Model.objective_direction` is a new revertible property to set maximization / minimization using the context manager.
+- Change output of `production_envelope` to wide data frame format [#587](https://github.com/opencobra/cobrapy/pull/587). Also allow multiple carbon source reactions and better handling of zero-division exceptions.
+- Enable summing lists of reactions, see [#596](https://github.com/opencobra/cobrapy/pull/596)
+
+## Deprecated features
+
+- `Model.get_metabolite_compartments` is deprecated (use
+ `Model.compartments` instead).
+- `Reaction.get_compartments` is deprecated (use
+ `Reaction.compartments` instead).
+
+## Backwards incompatible changes
+
+- The format of the dataframe `production_envelope` changed listing max and min on different columns instead of the same.
diff --git a/release-notes/next-release.md b/release-notes/next-release.md
index c622d35..7d5acdd 100644
--- a/release-notes/next-release.md
+++ b/release-notes/next-release.md
@@ -5,3 +5,5 @@
## New features
## Deprecated features
+
+## Backwards incompatible changes
diff --git a/scripts/compare-benchmark.py b/scripts/compare-benchmark.py
index c203dc4..686c659 100644
--- a/scripts/compare-benchmark.py
+++ b/scripts/compare-benchmark.py
@@ -4,6 +4,8 @@ import argparse
import re
from os.path import basename
+pd.set_option('display.width', 200)
+
def benchmark_to_df(json_file):
with open(json_file) as jf:
@@ -18,9 +20,13 @@ def benchmark_to_df(json_file):
if __name__ == "__main__":
parser = argparse.ArgumentParser(description="""
- compare cobrapy benchmarks.""")
- parser.add_argument('first', help='first branch to compare')
- parser.add_argument('second', help='second branch to compare')
+ compare cobrapy benchmarks.
+ Run pytest with
+ pytest --benchmark-save=without-cache --benchmark-min-rounds=20
+ then compare saved json files with this script.
+ """)
+ parser.add_argument('first', help='first json file')
+ parser.add_argument('second', help='second json file')
args = parser.parse_args()
first = benchmark_to_df(args.first)
@@ -31,5 +37,5 @@ if __name__ == "__main__":
both = pd.merge(first, second, how="inner", on="test",
suffixes=(first_name, second_name))
both["fraction"] = both.iloc[:, 2] / both.iloc[:, 1]
- both.sort_values(by="fraction")
- print(both)
+ print(both.sort_values(by="fraction"))
+
diff --git a/setup.cfg b/setup.cfg
index 6c7b437..73d72f1 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -1,5 +1,5 @@
[bumpversion]
-current_version = 0.8.2
+current_version = 0.9.0
commit = True
tag = True
parse = (?P<major>\d+)
diff --git a/setup.py b/setup.py
index be73b98..bc1dcb6 100644
--- a/setup.py
+++ b/setup.py
@@ -144,7 +144,7 @@ except:
setup(
name="cobra",
- version="0.8.2",
+ version="0.9.0",
packages=find_packages(),
setup_requires=setup_requirements,
install_requires=["future", "swiglpk", "optlang>=1.2.1",
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-cobra.git
More information about the debian-med-commit
mailing list