ENH: allows the user to choose which constraint method to use to calc…

…ulate lumps
EPFL-LCSB · May 13, 2019 · 0c26b16 · 0c26b16
1 parent b25a274
commit 0c26b16
Show file tree

Hide file tree

Showing 5 changed files with 51 additions and 20 deletions.
diff --git a/pytfa/redgem/lumpgem.py b/pytfa/redgem/lumpgem.py
@@ -48,8 +48,11 @@ def __init__(self, reaction, **kwargs):
                                   **kwargs)
 
 # Define a new constraint type:
-class UseOrKO(ReactionConstraint):
-    prefix = 'UK_'
+class UseOrKOInt(ReactionConstraint):
+    prefix = 'UKI_'
+# Define a new constraint type:
+class UseOrKOFlux(ReactionConstraint):
+    prefix = 'UKF_'
 
 
 class LumpGEM:
@@ -149,27 +152,48 @@ def init_params(self):
 
         self.timeout_limit = self._param_dict["timeout"]
 
+        self.constraint_method = self._param_dict["constraint_method"]
+
     def _generate_usage_constraints(self):
         """
         Generate carbon intake related constraints for each non-core reaction
         For each reaction rxn : rxn.forward_variable + rxn.reverse_variable + activation_var * C_uptake < C_uptake
         """
+        flux_methods = ['flux', 'fluxes', 'both']
+        int_methods = ['int', 'integer', 'both']
+
+        if self.constraint_method.lower() not in flux_methods + int_methods:
+            raise ArgumentError('{} is not a correct constraint method. '
+                                'Choose among [Flux, Integer, Both]. '
+                                'If you do not know what to choose, go for Flux.'
+                                'If it is too slow, go for integer.'
+                                'If you get strange lumps, go for both'
+                                .format(self.constraint_method))
 
         for rxn in self._rncore:
             activation_var = self._activation_vars[rxn]
-            bigM = 100
-            reac_var = rxn.forward_variable + rxn.reverse_variable + activation_var * bigM
-            # fu = self._tfa_model.forward_use_variable .get_by_id(rxn.id)
-            # bu = self._tfa_model.backward_use_variable.get_by_id(rxn.id)
-            # reac_var = fu + bu + activation_var
-            # adding the constraint to the model
-            self._tfa_model.add_constraint(kind=UseOrKO,
-                                           hook=rxn,
-                                           expr=reac_var,
-                                           # ub=1,
-                                           ub=bigM,
-                                           lb=0,
-                                           queue=True)
+            if self.constraint_method.lower() in flux_methods:
+                bigM = 100
+                reac_var = rxn.forward_variable + rxn.reverse_variable + activation_var * bigM
+                # adding the constraint to the model
+                self._tfa_model.add_constraint(kind=UseOrKOFlux,
+                                               hook=rxn,
+                                               expr=reac_var,
+                                               ub=bigM,
+                                               lb=0,
+                                               queue=True)
+            if self.constraint_method.lower() in int_methods:
+                fu = self._tfa_model.forward_use_variable .get_by_id(rxn.id)
+                bu = self._tfa_model.backward_use_variable.get_by_id(rxn.id)
+                reac_var = fu + bu + activation_var
+                # adding the constraint to the model
+                self._tfa_model.add_constraint(kind=UseOrKOInt,
+                                               hook=rxn,
+                                               expr=reac_var,
+                                               ub=1,
+                                               lb=0,
+                                               queue=True)
+
 
         # push constraints in one bulk (faster)
         self._tfa_model._push_queue()

diff --git a/tests/redgem_params.yml b/tests/redgem_params.yml
@@ -61,3 +61,4 @@ solver: auto
 force_solve: False
 timeout: 300
 feasibility: 1e-9
+constraint_method: integer # flux, integer, or both
diff --git a/tests/redgem_params_textbook.yaml b/tests/redgem_params_textbook.yaml
@@ -53,3 +53,4 @@ solver: auto
 force_solve: False
 timeout: 300
 feasibility: 1e-9
+constraint_method: flux # flux, integer, or both
diff --git a/tests/test_redgem.py b/tests/test_redgem.py
@@ -39,8 +39,6 @@
 
 else:
     path_to_model = join(this_directory, '..', 'models/small_ecoli.mat')
-    # path_to_model = join(this_directory,'..','models/GSmodel_Ecoli.mat')
-
     model = import_matlab_model(path_to_model)
     path_to_params = join(this_directory, '..', 'tests/redgem_params.yml')
 
@@ -50,9 +48,9 @@
 
 
 # Scaling to avoid numerical errors with bad lumps
-# for rxn in model.reactions:
-#     if rxn.id.startswith('LMPD_'):
-#         rxn.add_metabolites({x:v*(0.1 - 1) for x,v in rxn.metabolites.items()})
+for rxn in model.reactions:
+    if rxn.id.startswith('LMPD_'):
+        rxn.add_metabolites({x:v*(0.1 - 1) for x,v in rxn.metabolites.items()})
 
 thermo_data = load_thermoDB(thermoDB)
 lexicon = read_lexicon(path_to_lexicon)

diff --git a/tutorials/tuto_redgem_params.yaml b/tutorials/tuto_redgem_params.yaml
@@ -62,3 +62,10 @@ solver: auto
 force_solve: False
 timeout: 300
 feasibility: 1e-9
+
+# Optimisation shenanigans
+# Integer is faster with a good parallel solver, but might yield some numerical errors
+# Flux is slower but more reliable
+# Both is good if you can afford it
+constraint_method: both # flux, integer, or both.
+