Merge pull request #703 from carrascomj/refactor-stoichiometry

refactor: split stoichiometric consistency into three tests

HISTORY.rst

@@ -4,6 +4,7 @@ History
 Next Release
 ------------
 * Fix bug caused from upstream changes -> cobrapy -> pandas.
+* Refactor stoichiometric consistency test into 3 separate ones.
 
 0.11.1 (2020-08-10)
 -------------------

src/memote/suite/templates/test_config.yml

@@ -63,6 +63,8 @@ cards:
     - test_genes_presence
     - test_compartments_presence
     - test_metabolic_coverage
+    - test_unconserved_metabolites
+    - test_inconsistent_min_stoichiometry
   basic_met:
     title: "Metabolite Information"
     cases:

src/memote/suite/tests/test_consistency.py

@@ -41,42 +41,81 @@ def test_stoichiometric_consistency(model):
     either produce mass from nothing or consume mass from the model.
 
     Implementation:
-    This test first uses an implementation of the algorithm presented in
+    This test uses an implementation of the algorithm presented in
     section 3.1 by Gevorgyan, A., M. G Poolman, and D. A Fell.
     "Detection of Stoichiometric Inconsistencies in Biomolecular Models."
     Bioinformatics 24, no. 19 (2008): 2245.
     doi: 10.1093/bioinformatics/btn425
-    Should the model be inconsistent, then the list of unconserved metabolites
-    is computed using the algorithm described in section 3.2 of the same
-    publication. In addition, the list of min unconservable sets is computed
-    using the algorithm described in section 3.3.
 
     """
     ann = test_stoichiometric_consistency.annotation
     is_consistent = consistency.check_stoichiometric_consistency(model)
-    ann["data"] = {
-        "unconserved_metabolites": []
-        if is_consistent
-        else get_ids(consistency.find_unconserved_metabolites(model)),
-        "minimal_unconservable_sets": []
-        if is_consistent
-        else [
+    ann["data"] = is_consistent
+    ann["metric"] = 1.0 - float(is_consistent)
+    ann["message"] = wrapper.fill(
+        """This model's stoichiometry {}""".format(
+            "consistent" if is_consistent else "inconsistent"
+        )
+    )
+    assert is_consistent, ann["message"]
+
+
+@annotate(title="Uncoserved Metabolites", format_type="count")
+def test_unconserved_metabolites(model):
+    """
+    Report number all unconserved metabolites.
+
+    The list of unconserved metabolites is computed using the algorithm described
+    in section 3.2 of
+    "Detection of Stoichiometric Inconsistencies in Biomolecular Models."
+    Bioinformatics 24, no. 19 (2008): 2245.
+    doi: 10.1093/bioinformatics/btn425.
+    """
+    ann = test_unconserved_metabolites.annotation
+    is_consistent = consistency.check_stoichiometric_consistency(model)
+    ann["data"] = []
+    if not is_consistent:
+        ann["data"] = get_ids(consistency.find_unconserved_metabolites(model))
+    ann["metric"] = len(ann["data"])
+    ann["message"] = wrapper.fill(
+        """This model contains {} unconserved metabolites: {}""".format(
+            ann["metric"],
+            truncate(ann["data"]),
+        )
+    )
+    assert ann["metric"] == 0, ann["message"]
+
+
+@annotate(title="Minimal Inconsistent Net Stoichiometries", format_type="count")
+def test_inconsistent_min_stoichiometry(model):
+    """
+    Report inconsistent min stoichiometries.
+
+    Only 10 unconserved metabolites are reported and considered to
+    avoid computing for too long.
+
+    Implementation:
+    Algorithm described in section 3.3 of
+    "Detection of Stoichiometric Inconsistencies in Biomolecular Models."
+    Bioinformatics 24, no. 19 (2008): 2245.
+    doi: 10.1093/bioinformatics/btn425.
+    """
+    ann = test_inconsistent_min_stoichiometry.annotation
+    is_consistent = consistency.check_stoichiometric_consistency(model)
+    ann["data"] = []
+    if not is_consistent:
+        ann["data"] = [
             get_ids(mets)
             for mets in consistency.find_inconsistent_min_stoichiometry(model)
-        ],
-    }
-    ann["metric"] = len(ann["data"]["unconserved_metabolites"]) / len(model.metabolites)
+        ]
+    ann["metric"] = len(ann["data"])
     ann["message"] = wrapper.fill(
-        """This model contains {} ({:.2%}) unconserved
-        metabolites: {}; and {} minimal unconservable sets: {}""".format(
-            len(ann["data"]["unconserved_metabolites"]),
-            ann["metric"],
-            truncate(ann["data"]["unconserved_metabolites"]),
-            len(ann["data"]["minimal_unconservable_sets"]),
-            truncate(ann["data"]["minimal_unconservable_sets"]),
+        """This model contains {} minimal unconservable sets: {}""".format(
+            ann["metric"] if ann["metric"] > 10 else "more than 10",
+            truncate(ann["data"]),
         )
     )
-    assert is_consistent, ann["message"]
+    assert ann["metric"] == 0, ann["message"]
 
 
 @pytest.mark.parametrize("met", [x for x in consistency.ENERGY_COUPLES])
@@ -88,7 +127,7 @@ def test_stoichiometric_consistency(model):
     metric=dict(),
 )
 def test_detect_energy_generating_cycles(model, met):
-    u"""
+    """
     Expect that no energy metabolite can be produced out of nothing.
 
     When a model is not sufficiently constrained to account for the

src/memote/support/consistency.py

@@ -186,7 +186,7 @@ def find_unconserved_metabolites(model):
         )
 
 
-def find_inconsistent_min_stoichiometry(model, atol=1e-13):
+def find_inconsistent_min_stoichiometry(model, atol=1e-13, max_mets_computed=100):
     """
     Detect inconsistent minimal net stoichiometries.
 
@@ -197,6 +197,9 @@ def find_inconsistent_min_stoichiometry(model, atol=1e-13):
     atol : float, optional
         Values below the absolute tolerance are treated as zero. Expected to be
         very small but larger than zero.
+    max_mets_computed: int, optional
+        To avoid computing for too long, a soft cap is added to the number of
+        computed unconserved metabolites (trivial cases are ignored).
 
     Notes
     -----
@@ -239,6 +242,7 @@ def find_inconsistent_min_stoichiometry(model, atol=1e-13):
     LOGGER.debug("%s", str(problem))
     inc_minimal = set()
     cuts = list()
+    n_computed = 0
     for met in unconserved_mets:
         # always add the met as an uncoserved set if there is no left nullspace
         row = met_index[met]
@@ -247,6 +251,9 @@ def find_inconsistent_min_stoichiometry(model, atol=1e-13):
             # singleton set!
             inc_minimal.add((met,))
             continue
+        if n_computed >= max_mets_computed:
+            LOGGER.debug("max number of computed unconserved metabolites reached")
+            break
         # expect a positive mass for the unconserved metabolite
         problem.variables[met.id].lb = 1e-3
         status = problem.optimize()
@@ -264,6 +271,7 @@ def find_inconsistent_min_stoichiometry(model, atol=1e-13):
                 if var.primal > 0.2
             ]
             LOGGER.debug("%s: set size %d", met.id, len(solution))
+            n_computed += 1
             inc_minimal.add(tuple(solution))
             if len(solution) == 1:
                 break