[WIP] Add Turbo Compass.

compass/compass/algorithm.py

@@ -2,6 +2,7 @@
 Run the procedure for COMPASS
 """
 from __future__ import print_function, division, absolute_import
+import numpy as np
 import pandas as pd
 from tqdm import tqdm
 from random import shuffle
@@ -13,6 +14,7 @@
 from .. import models
 from . import cache
 from ..globals import BETA, EXCHANGE_LIMIT
+import compass.global_state as global_state
 
 import cplex
 
@@ -75,7 +77,7 @@ def singleSampleCompass(data, model, media, directory, sample_index, args):
 
     logger.info("Processing Sample %i/%i: %s", sample_index,
                 len(expression.columns), sample_name)
-
+    global_state.set_current_cell_name(sample_name)
     # Run core compass algorithm
 
     # Evaluate reaction penalties
@@ -198,7 +200,7 @@ def compass_exchange(model, problem, reaction_penalties, select_reactions=None,
 
             uptake_scores[met_id] = 0.0
             secretion_scores[met_id] = 0.0
-            continue
+            continue  # In test mode this always continues!
 
         # Rectify exchange reactions
         # Either find existing pos and neg exchange reactions
@@ -226,7 +228,7 @@ def compass_exchange(model, problem, reaction_penalties, select_reactions=None,
             reactions = [r for r in reactions if ((r.id)[:-4] in selected_reaction_ids)]
 
 
-        for reaction in reactions:
+        for reaction in reactions:  # This only effectively loops over exchange reactions in fact.
             if reaction.is_exchange and met_id in reaction.products:
                 if uptake_rxn is None:
                     uptake_rxn = reaction.id
@@ -319,8 +321,8 @@ def compass_exchange(model, problem, reaction_penalties, select_reactions=None,
             problem.objective.set_name('reaction_penalties')
             problem.objective.set_sense(problem.objective.sense.minimize)
 
-        problem.solve()
-        value = problem.solution.get_objective_value()
+        global_state.set_current_reaction_id(secretion_rxn)
+        value = solve_problem_wrapper(problem)
         secretion_scores[met_id] = value
 
         # Clear Secretion constraint
@@ -371,8 +373,8 @@ def compass_exchange(model, problem, reaction_penalties, select_reactions=None,
             problem.objective.set_name('reaction_penalties')
             problem.objective.set_sense(problem.objective.sense.minimize)
 
-        problem.solve()
-        value = problem.solution.get_objective_value()
+        global_state.set_current_reaction_id(uptake_rxn)
+        value = solve_problem_wrapper(problem)
         uptake_scores[met_id] = value
 
         # Clear Secretion constraint
@@ -436,7 +438,6 @@ def compass_reactions(model, problem, reaction_penalties, select_reactions=None,
         reactions = [r for r in reactions if ((r.id)[:-4] in selected_reaction_ids)]
 
 
-
     for reaction in tqdm(reactions, file=sys.stderr):
 
         if reaction.is_exchange:
@@ -474,8 +475,8 @@ def compass_reactions(model, problem, reaction_penalties, select_reactions=None,
                 problem.objective.set_name('reaction_penalties')
                 problem.objective.set_sense(problem.objective.sense.minimize)
 
-            problem.solve()
-            value = problem.solution.get_objective_value()
+            global_state.set_current_reaction_id(reaction.id)
+            value = solve_problem_wrapper(problem)
             reaction_scores[reaction.id] = value
 
             # Remove Constraint
@@ -552,11 +553,23 @@ def maximize_reaction(model, problem, rxn, use_cache=True):
     problem.objective.set_name(rxn)
     problem.objective.set_sense(problem.objective.sense.maximize)
 
-    problem.solve()
-    rxn_max = problem.solution.get_objective_value()
+    global_state.set_current_reaction_id(rxn)
+    rxn_max = solve_problem_wrapper(problem)
 
     # Save the result
     model_cache = cache.load(model)
     model_cache[rxn] = rxn_max
 
     return rxn_max
+
+
+def solve_problem_wrapper(problem) -> float:
+    r"""
+    Only solve the problem if the reaction is selected for the cell. Else,
+    skip the computation and return np.nan.
+    """
+    if global_state.current_reaction_is_selected_for_current_cell():
+        problem.solve()
+        return problem.solution.get_objective_value()
+    else:
+        return np.nan

compass/global_state.py

@@ -0,0 +1,76 @@
+r"""
+This module contains global state. It is used for e.g. determining whether
+to skip the computation of a reaction for a cell when --selected-reactions-for-each-cell
+is provided.
+"""
+import os
+from typing import Optional
+
+_current_cell_name = -1
+_current_reaction_id = ""
+# This will hold an object of the class SelectedReactionsForEachCell that
+# will allow us to know what reactions are selected for each cell.
+_selected_reactions_for_each_cell = None
+
+
+class _SelectedReactionsForEachCell():
+    def __init__(self, path: Optional[str]):
+        r"""
+        :param path: Path to file containing the reactions for each cell.
+            The file contains one line per cell of interest. Each line is
+            comma-separated. It starts with the name of the cell, and is
+            followed by the reaction ids of interest.
+            E.g. 'cell42,BUP2_pos,BUP2_neg,DHPM2_pos' is a valid line.
+            If no path is provided, then no (cells, reaction) pairs will
+            be excluded by this method.
+        """
+        if path is None:
+            reactions_for_cells = None
+        else:
+            reactions_for_cells = {}
+            if not os.path.exists(path):
+                raise Exception(f"cannot find selected reactions for each cell: file {path}")
+
+            with open(path) as f:
+                for line in f:
+                    # Must be careful to remove newlines from end of file.
+                    tokens = [token.rstrip("\n") for token in line.split(',')]
+                    cell_name, reactions_for_this_cell = tokens[0], tokens[1:]
+                    reactions_for_cells[cell_name] = reactions_for_this_cell
+        self.reactions_for_cells = reactions_for_cells
+
+    def is_selected(self, cell_name: str, reaction_id: str) -> bool:
+        r"""
+        Whether the reaction_id is selected for cell_name.
+        """
+        if self.reactions_for_cells is None:
+            # No file was supplied: everything is selected by default.
+            return True
+        return reaction_id in self.reactions_for_cells.get(cell_name, [])
+
+
+def init_selected_reactions_for_each_cell(path: Optional[str]) -> None:
+    global _selected_reactions_for_each_cell
+    _selected_reactions_for_each_cell = _SelectedReactionsForEachCell(path)
+
+
+def set_current_cell_name(cell_name: str) -> None:
+    global _current_cell_name
+    _current_cell_name = cell_name
+
+
+def set_current_reaction_id(reaction_id: str) -> None:
+    global _current_reaction_id
+    _current_reaction_id = reaction_id
+
+
+def get_current_cell_name() -> str:
+    return _current_cell_name
+
+
+def get_current_reaction_id() -> str:
+    return _current_reaction_id
+
+
+def current_reaction_is_selected_for_current_cell() -> bool:
+    return _selected_reactions_for_each_cell.is_selected(_current_cell_name, _current_reaction_id)

compass/globals.py

@@ -33,16 +33,15 @@ def init_logger(directory="."):
     fh = logging.FileHandler(log_file, mode='w')
     fh.name = 'logfile'
 
-    formatter = logging.Formatter(
-        '%(message)s')
+    formatter = logging.Formatter("%(levelname)s %(lineno)s %(name)s: %(message)s")
     fh.setLevel(logging.DEBUG)
     fh.setFormatter(formatter)
 
     logger.addHandler(fh)
 
     # Add stream to stdout
     sh = logging.StreamHandler(sys.stdout)
-    formatter = logging.Formatter('%(message)s')
+    formatter = logging.Formatter("%(levelname)s %(lineno)s %(name)s: %(message)s")
     sh.setFormatter(formatter)
     sh.setLevel(logging.INFO)
     logger.addHandler(sh)