Merge branch 'release'

NEWS.md

@@ -1,3 +1,18 @@
+v0.22 (2016-07-01)
+------------------
+
+- Better unicode handling in commands.
+- When running the `gapfill` command the epsilon parameter can now be
+  specified on the command line.
+- When parsing reaction and compound entities from the YAML files, produce
+  better error messages when IDs are invalid.
+- Work around a bug in Cplex that in rare causes a segmentation fault when a
+  linear programming problem is solved repeatedly.
+- API: Add `fastgapfill` module which allows access to run the fastGapFill
+  algorithm.
+- API: Add `randomsparse` module which allows access to generate a random
+  minimal model which satisfies the flux threshold of the objective reaction.
+
 v0.21 (2016-06-09)
 ------------------
 

docs/api/fastgapfill.rst

@@ -0,0 +1,6 @@
+
+``psamm.fastgapfill`` -- FastGapFill algorithm
+==============================================
+
+.. automodule:: psamm.fastgapfill
+   :members:

docs/api/randomsparse.rst

@@ -0,0 +1,6 @@
+
+``psamm.randomsparse`` -- Find a random minimal network of model reactions
+===========================================================================
+
+.. automodule:: psamm.randomsparse
+   :members:

psamm/balancecheck.py

@@ -16,6 +16,8 @@
 # Copyright 2016  Jon Lund Steffensen <[email protected]>
 # Copyright 2016  Chao Liu <[email protected]>
 
+from __future__ import unicode_literals
+
 import operator
 import logging
 

psamm/command.py

@@ -25,17 +25,18 @@
 The :func:`.main` function is the entry point of command line interface.
 """
 
-from __future__ import division
+from __future__ import division, unicode_literals
 
 import os
+import sys
 import argparse
 import logging
 import abc
 from itertools import islice
 import multiprocessing as mp
 
 import pkg_resources
-from six import add_metaclass, iteritems, itervalues
+from six import add_metaclass, iteritems, itervalues, text_type
 
 from . import __version__ as package_version
 from .datasource.native import NativeModel
@@ -76,7 +77,7 @@ def __init__(self, model, args):
 
         name = self._model.get_name()
         if name is None:
-            name = str(self._model.context)
+            name = text_type(self._model.context)
         logger.info('Model: {}'.format(name))
 
         version = util.git_try_describe(self._model.context.basepath)
@@ -91,6 +92,17 @@ def init_parser(cls, parser):
     def run(self):
         """Execute command"""
 
+    def argument_error(self, msg):
+        """Raise error indicating error parsing an argument."""
+        raise CommandError(msg)
+
+    def fail(self, msg, exc=None):
+        """Exit command as a result of a failure."""
+        logger.error(msg)
+        if exc is not None:
+            logger.debug('Command failure caused by exception!', exc_info=exc)
+        sys.exit(1)
+
 
 class MetabolicMixin(object):
     """Mixin for commands that use a metabolic model representation."""
@@ -396,4 +408,4 @@ def main(command_class=None, args=None):
     try:
         command.run()
     except CommandError as e:
-        parser.error(str(e))
+        parser.error(text_type(e))

psamm/commands/excelexport.py

@@ -19,7 +19,9 @@
 
 import logging
 
-from ..command import Command, CommandError
+from six import text_type
+
+from ..command import Command
 
 try:
     import xlsxwriter
@@ -40,8 +42,7 @@ def init_parser(cls, parser):
     def run(self):
         model = self._model
         if xlsxwriter is None:
-            raise CommandError(
-                'Excel export requires the XlsxWriter python module')
+            self.fail('Excel export requires the XlsxWriter python module')
         workbook = xlsxwriter.Workbook(self._args.file)
         reaction_sheet = workbook.add_worksheet(name='Reactions')
 
@@ -53,10 +54,11 @@ def run(self):
                                       key=lambda x: (x != 'id',
                                                      x != 'equation', x))
         for z, i in enumerate(property_list_sorted):
-            reaction_sheet.write_string(0, z, str(i))
+            reaction_sheet.write_string(0, z, text_type(i))
         for x, i in enumerate(model.parse_reactions()):
             for y, j in enumerate(property_list_sorted):
-                reaction_sheet.write_string(x+1, y, str(i.properties.get(j)))
+                reaction_sheet.write_string(
+                    x+1, y, text_type(i.properties.get(j)))
 
         compound_sheet = workbook.add_worksheet(name='Compounds')
 
@@ -69,10 +71,11 @@ def run(self):
                                                      x != 'name', x))
 
         for z, i in enumerate(compound_list_sorted):
-            compound_sheet.write_string(0, z, str(i))
+            compound_sheet.write_string(0, z, text_type(i))
         for x, i in enumerate(model.parse_compounds()):
             for y, j in enumerate(compound_list_sorted):
-                compound_sheet.write_string(x+1, y, str(i.properties.get(j)))
+                compound_sheet.write_string(
+                    x+1, y, text_type(i.properties.get(j)))
 
         media_sheet = workbook.add_worksheet(name='Medium')
 
@@ -82,8 +85,6 @@ def run(self):
         media_sheet.write_string(0, 3, 'Upper Limit')
 
         default_flux = model.get_default_flux_limit()
-        if default_flux is None:
-            default_flux = 1000
 
         for x, (compound, reaction, lower, upper) in enumerate(
                 model.parse_medium()):
@@ -93,16 +94,16 @@ def run(self):
             if upper is None:
                 upper = default_flux
 
-            media_sheet.write(x+1, 0, str(compound))
-            media_sheet.write(x+1, 1, str(reaction))
-            media_sheet.write(x+1, 2, str(lower))
-            media_sheet.write(x+1, 3, str(upper))
+            media_sheet.write(x+1, 0, text_type(compound))
+            media_sheet.write(x+1, 1, text_type(reaction))
+            media_sheet.write(x+1, 2, text_type(lower))
+            media_sheet.write(x+1, 3, text_type(upper))
 
         limits_sheet = workbook.add_worksheet(name='Limits')
 
         limits_sheet.write_string(0, 0, 'Reaction ID')
         limits_sheet.write_string(0, 1, 'Lower Limit')
-        limits_sheet.write_string(0, 2, str('Upper Limit'))
+        limits_sheet.write_string(0, 2, 'Upper Limit')
 
         for x, i in enumerate(model.parse_limits()):
             limits_sheet.write(x, 0, (i[0]))

psamm/commands/fastgapfill.py

@@ -20,8 +20,9 @@
 import argparse
 import logging
 
-from ..command import Command, MetabolicMixin, SolverCommandMixin, CommandError
+from ..command import Command, MetabolicMixin, SolverCommandMixin
 from .. import fastcore, fluxanalysis
+from ..fastgapfill import create_extended_model, fastgapfill
 
 logger = logging.getLogger(__name__)
 
@@ -69,71 +70,53 @@ def run(self):
             compound_name[compound.id] = (
                 compound.name if compound.name is not None else compound.id)
 
-        epsilon = self._args.epsilon
-        model_compartments = set(self._mm.compartments)
-        extra_comp = self._model.get_extracellular_compartment()
-
-        # Add exchange and transport reactions to database
-        model_complete = self._mm.copy()
-        logger.info('Adding database, exchange and transport reactions')
-        db_added = model_complete.add_all_database_reactions(
-            model_compartments)
-        ex_added = model_complete.add_all_exchange_reactions(extra_comp)
-        tp_added = model_complete.add_all_transport_reactions(extra_comp)
-
         # TODO: The exchange and transport reactions have tuple names. This
         # means that in Python 3 the reactions can no longer be directly
         # compared (e.g. while sorting) so define this helper function as a
         # workaround.
         def reaction_key(r):
             return r if isinstance(r, tuple) else (r,)
 
-        # Add penalty weights on reactions
-        weights = {}
-        if self._args.db_weight is not None:
-            weights.update((rxnid, self._args.db_weight) for rxnid in db_added)
-        if self._args.tp_weight is not None:
-            weights.update((rxnid, self._args.tp_weight) for rxnid in tp_added)
-        if self._args.ex_weight is not None:
-            weights.update((rxnid, self._args.ex_weight) for rxnid in ex_added)
-
+        # Calculate penalty if penalty file exists
+        penalties = {}
         if self._args.penalty is not None:
             for line in self._args.penalty:
                 line, _, comment = line.partition('#')
                 line = line.strip()
                 if line == '':
                     continue
                 rxnid, weight = line.split(None, 1)
-                weights[rxnid] = float(weight)
+                penalties[rxnid] = float(weight)
 
-        # Run Fastcore and print the induced reaction set
-        logger.info('Calculating Fastcore induced set on model')
-        core = set(self._mm.reactions)
+        model_extended, weights = create_extended_model(
+            self._model,
+            db_weight=self._args.db_weight,
+            ex_weight=self._args.ex_weight,
+            tp_weight=self._args.tp_weight,
+            penalties=penalties)
 
-        induced = fastcore.fastcore(model_complete, core, epsilon,
-                                    weights=weights, solver=solver)
-        logger.info('Result: |A| = {}, A = {}'.format(len(induced), induced))
-        added_reactions = induced - core
-        logger.info('Extended: |E| = {}, E = {}'.format(
-            len(added_reactions), added_reactions))
+        epsilon = self._args.epsilon
+        core = set(self._mm.reactions)
+        induced = fastgapfill(model_extended, core, weights=weights,
+                              epsilon=epsilon, solver=solver)
 
         if self._args.reaction is not None:
             maximized_reaction = self._args.reaction
         else:
             maximized_reaction = self._model.get_biomass_reaction()
             if maximized_reaction is None:
-                raise CommandError('The maximized reaction was not specified')
+                self.argument_error(
+                    'The maximized reaction was not specified')
 
         if not self._mm.has_reaction(maximized_reaction):
-            raise CommandError(
-                'The biomass reaction is not a valid model'
-                ' reaction: {}'.format(maximized_reaction))
+            self.fail('The biomass reaction is not a valid model'
+                      ' reaction: {}'.format(maximized_reaction))
 
         logger.info('Flux balance on induced model maximizing {}'.format(
             maximized_reaction))
-        model_induced = self._mm.copy()
-        for rxnid in induced:
-            model_induced.add_reaction(rxnid)
+        model_induced = model_extended.copy()
+        for rxnid in set(model_extended.reactions) - induced:
+            model_induced.remove_reaction(rxnid)
         for rxnid, flux in sorted(fluxanalysis.flux_balance(
                 model_induced, maximized_reaction, tfba=enable_tfba,
                 solver=solver), key=lambda x: (reaction_key(x[0]), x[1])):
@@ -142,7 +125,7 @@ def reaction_key(r):
             if self._mm.has_reaction(rxnid):
                 reaction_class = 'Model'
                 weight = 0
-            rx = model_complete.get_reaction(rxnid)
+            rx = model_induced.get_reaction(rxnid)
             rxt = rx.translated_compounds(lambda x: compound_name.get(x, x))
             print('{}\t{}\t{}\t{}\t{}'.format(
                 rxnid, reaction_class, weight, flux, rxt))

psamm/commands/fba.py

@@ -20,7 +20,7 @@
 import time
 import logging
 
-from ..command import SolverCommandMixin, MetabolicMixin, Command, CommandError
+from ..command import SolverCommandMixin, MetabolicMixin, Command
 from .. import fluxanalysis
 
 logger = logging.getLogger(__name__)
@@ -65,11 +65,11 @@ def run(self):
         else:
             reaction = self._model.get_biomass_reaction()
             if reaction is None:
-                raise CommandError('The biomass reaction was not specified')
+                self.argument_error('The biomass reaction was not specified')
 
         if not self._mm.has_reaction(reaction):
-            raise CommandError('Specified reaction is not in model: {}'.format(
-                reaction))
+            self.fail(
+                'Specified reaction is not in model: {}'.format(reaction))
 
         logger.info('Using {} as objective'.format(reaction))
 
@@ -84,8 +84,8 @@ def run(self):
             logger.info('Loop removal using thermodynamic constraints')
             result = self.run_tfba(reaction)
         else:
-            raise CommandError('Invalid loop constraint mode: {}'.format(
-                loop_removal))
+            self.argument_error(
+                'Invalid loop constraint mode: {}'.format(loop_removal))
 
         optimum = None
         total_reactions = 0

psamm/commands/fluxcheck.py

@@ -22,7 +22,7 @@
 from itertools import product
 
 from ..command import (Command, MetabolicMixin, SolverCommandMixin,
-                       ParallelTaskMixin, CommandError)
+                       ParallelTaskMixin)
 from .. import fluxanalysis, fastcore
 
 logger = logging.getLogger(__name__)
@@ -82,10 +82,9 @@ def run(self):
         enable_fastcore = self._args.fastcore
 
         if enable_tfba and enable_fastcore:
-            raise CommandError(
+            self.argument_error(
                 'Using Fastcore with thermodynamic constraints'
                 ' is not supported!')
-
         start_time = time.time()
 
         if enable_fastcore:

psamm/commands/fluxcoupling.py

@@ -20,7 +20,7 @@
 import logging
 
 from ..command import (Command, SolverCommandMixin, MetabolicMixin,
-                       ParallelTaskMixin, CommandError)
+                       ParallelTaskMixin)
 from .. import fluxanalysis, fluxcoupling
 
 logger = logging.getLogger(__name__)
@@ -43,7 +43,7 @@ def run(self):
 
         max_reaction = self._model.get_biomass_reaction()
         if max_reaction is None:
-            raise CommandError('The biomass reaction was not specified')
+            self.fail('The biomass reaction was not specified')
 
         fba_fluxes = dict(fluxanalysis.flux_balance(
             self._mm, max_reaction, tfba=False, solver=solver))

psamm/commands/fva.py

@@ -22,7 +22,7 @@
 from itertools import product
 
 from ..command import (Command, SolverCommandMixin, MetabolicMixin,
-                       ParallelTaskMixin, CommandError)
+                       ParallelTaskMixin)
 from ..util import MaybeRelative
 from .. import fluxanalysis
 
@@ -60,11 +60,11 @@ def run(self):
         else:
             reaction = self._model.get_biomass_reaction()
             if reaction is None:
-                raise CommandError('The biomass reaction was not specified')
+                self.argument_error('The biomass reaction was not specified')
 
         if not self._mm.has_reaction(reaction):
-            raise CommandError('Specified reaction is not in model: {}'.format(
-                reaction))
+            self.fail(
+                'Specified reaction is not in model: {}'.format(reaction))
 
         enable_tfba = self._args.tfba
         if enable_tfba: