Align task lists, GPR parsing and YAML warnings with RAVEN

src/raven_toolbox/init/genes.py

@@ -79,7 +79,8 @@ def remove_low_score_genes(
         negative alternatives was cut to one of several alternatives with the same
         best score are appended to it, once per reaction. In such a rule the first of
         the tied alternatives in the rule is kept, so the result depends on the order
-        of the rule; RAVEN picks one of them at random.
+        of the rule. RAVEN keeps the first one too (CompMetEng/RAVEN#11); before that
+        it kept the alphabetically first.
 
     Returns
     -------

src/raven_toolbox/io/yaml.py

@@ -779,6 +779,17 @@ def write_yaml_model(
         version = stored_meta.get("version")
     foreign = model_notes.pop("_yaml_sections", None) or {}
 
+    unnamed = [met.id for met in model.metabolites if not met.name]
+    if unnamed:
+        # RAVEN's writeYAMLmodel stops on these (checkModelStruct: "metabolites have
+        # empty names"), so a file written here with `name: ""` may not load there.
+        shown = ", ".join(unnamed[:5]) + (", ..." if len(unnamed) > 5 else "")
+        warnings.warn(
+            f"{len(unnamed)} metabolite(s) have an empty name ({shown}); the YAML is written "
+            "anyway, but RAVEN's writeYAMLmodel refuses such a model.",
+            stacklevel=2,
+        )
+
     doc = OrderedDict(_to_plain(model_to_dict(model)))
 
     # cobra's model_to_dict serialises model.notes verbatim into doc["notes"],

src/raven_toolbox/tasks/check.py

@@ -20,6 +20,7 @@
 
 import concurrent.futures as cf
 import pickle
+import warnings
 from collections.abc import Iterable
 from dataclasses import dataclass
 from pathlib import Path
@@ -116,6 +117,13 @@ def touch(mid: str) -> list[float]:
                 if not ids:
                     missing.append(token)
                     continue
+                if len(ids) > 1:
+                    warnings.warn(
+                        f"task {task.id!r}: {token!r} matches {len(ids)} model metabolites "
+                        f"({', '.join(ids)}); all of them are constrained (RAVEN constrains "
+                        "only the first).",
+                        stacklevel=2,
+                    )
             for mid in ids:
                 b = touch(mid)
                 if is_input:
@@ -138,8 +146,9 @@ def task_name_maps(model: cobra.Model) -> tuple[dict[str, list[str]], dict[str,
     """Build ``name[comp]→[ids]`` and ``comp→[ids]`` lookups for a model's metabolites.
 
     ``name[comp]`` maps to a *list* because a model can carry several metabolites with
-    the same name and compartment; a task referencing it constrains all of them (as
-    RAVEN does), rather than an arbitrary one.
+    the same name and compartment. A task referencing it constrains all of them. RAVEN's
+    ``checkTasks`` constrains only the first match; both sides warn when a task hits
+    duplicates, so the difference is visible instead of silent.
     """
     name_to_ids: dict[str, list[str]] = {}
     comp_to_ids: dict[str, list[str]] = {}
@@ -177,10 +186,15 @@ def apply_task_constraints(
             tmp = model.reactions.get_by_id(f"TASK_TMP_{i}")
             task_mets |= {m.id for m in tmp.metabolites if m.id in existing}
 
-    for rxn_id, lb, ub in task.changed:
+    for rxn_id, new_lb, new_ub in task.changed:
         if rxn_id not in model.reactions:
             return set(), f"CHANGED RXN not in model: {rxn_id!r}"
-        model.reactions.get_by_id(rxn_id).bounds = (lb, ub)
+        rxn = model.reactions.get_by_id(rxn_id)
+        # An empty bound in the task list means "leave this bound as it is".
+        rxn.bounds = (
+            rxn.lower_bound if new_lb is None else new_lb,
+            rxn.upper_bound if new_ub is None else new_ub,
+        )
 
     model.objective = model.problem.Objective(Zero, direction="max")  # feasibility only
     return task_mets, None

src/raven_toolbox/tasks/tasklist.py

@@ -14,7 +14,9 @@
 A task spans consecutive rows; only its first row carries an ID. Metabolites are
 written ``name[compartment]``; several in one cell are separated by ``;`` (sharing
 that row's bounds). ``IN``/``OUT`` default LB 0, UB 1000; ``EQU`` defaults LB
--1000 if reversible (``<=>``) else 0, UB 1000. The special tokens ``ALLMETS`` and
+-1000 if reversible (``<=>``) else 0, UB 1000. An empty ``CHANGED LB``/``CHANGED UB``
+leaves that bound of the reaction unchanged. ``SHOULD FAIL``/``PRINT FLUX`` are false
+when empty or ``0``/``false``/``no`` (any case), true otherwise. The special tokens ``ALLMETS`` and
 ``ALLMETSIN[comp]`` allow free uptake/excretion of all metabolites (only the upper
 bound is used).
 """
@@ -34,7 +36,11 @@
 
 @dataclass
 class Task:
-    """One metabolic task. Bounds are ``(metabolite_or_reaction, lb, ub)`` triples."""
+    """One metabolic task. Bounds are ``(metabolite_or_reaction, lb, ub)`` triples.
+
+    In ``changed``, a bound is ``None`` when the cell was empty, meaning "leave that
+    bound of the reaction as it is" (RAVEN reads it as NaN).
+    """
 
     id: str
     description: str = ""
@@ -44,7 +50,7 @@ class Task:
     inputs: list[tuple[str, float, float]] = field(default_factory=list)
     outputs: list[tuple[str, float, float]] = field(default_factory=list)
     equations: list[tuple[str, float, float]] = field(default_factory=list)
-    changed: list[tuple[str, float, float]] = field(default_factory=list)
+    changed: list[tuple[str, float | None, float | None]] = field(default_factory=list)
 
 
 def _truthy(value: str) -> bool:
@@ -56,8 +62,18 @@ def _num(value: str, default: float) -> float:
     return float(value) if value else default
 
 
+def _opt_num(value: str) -> float | None:
+    value = value.strip()
+    return float(value) if value else None
+
+
 def _read_rows(path: str | Path) -> list[list[str]]:
     path = Path(path)
+    if path.suffix.lower() == ".xls":
+        raise ValueError(
+            f"{path}: the old binary .xls format is not supported. Save the task list as "
+            ".xlsx (sheet 'TASKS') or as a tab-delimited .txt file."
+        )
     if path.suffix.lower() in (".xlsx", ".xlsm"):
         try:
             from openpyxl import load_workbook
@@ -145,5 +161,5 @@ def _add_row(task: Task, row: list[str], cell) -> None:
         ub = _num(cell(row, "EQU UB"), 1000.0)
         task.equations.append((equ.strip(), lb, ub))
     if chg := cell(row, "CHANGED RXN"):
-        lb, ub = _num(cell(row, "CHANGED LB"), -1000.0), _num(cell(row, "CHANGED UB"), 1000.0)
-        task.changed += [(r.strip(), lb, ub) for r in chg.split(";") if r.strip()]
+        chg_lb, chg_ub = _opt_num(cell(row, "CHANGED LB")), _opt_num(cell(row, "CHANGED UB"))
+        task.changed += [(r.strip(), chg_lb, chg_ub) for r in chg.split(";") if r.strip()]

src/raven_toolbox/utils/gpr.py

@@ -17,7 +17,9 @@
 from __future__ import annotations
 
 import ast
+import logging
 import statistics
+import warnings
 from collections.abc import Callable
 from dataclasses import dataclass
 
@@ -88,19 +90,48 @@ def is_dnf(gpr: GPR | str | None) -> bool:
         A cobra :class:`~cobra.core.gene.GPR`, a grRule string, or ``None``.
         An empty/``None`` rule is trivially DNF.
 
+    Raises
+    ------
+    ValueError
+        If a rule string is not valid (unbalanced parentheses, a dangling
+        operator, two genes without an operator between them, ``()``).
+
     Examples
     --------
     >>> is_dnf("(G1 and G2) or G3")
     True
     >>> is_dnf("(G1 or G2) and G3")
     False
     """
-    parsed = GPR.from_string(gpr) if isinstance(gpr, str) else gpr
+    parsed = _parse_rule(gpr) if isinstance(gpr, str) else gpr
     if parsed is None:
         return True
     return _is_dnf_node(parsed.body)
 
 
+def _parse_rule(rule: str) -> GPR | None:
+    """Parse a grRule string; ``None`` for a blank rule, ``ValueError`` if malformed.
+
+    ``GPR.from_string`` swallows a syntax error, logs it and returns an empty GPR, which
+    would make any malformed rule look trivially DNF. RAVEN's ``isDnfGrRule`` raises
+    ``RAVEN:badGrRule`` instead, so this does too.
+    """
+    if not rule.strip():
+        return None
+    gene_logger = logging.getLogger("cobra.core.gene")
+    was_disabled = gene_logger.disabled
+    gene_logger.disabled = True
+    try:
+        with warnings.catch_warnings():
+            warnings.simplefilter("ignore")
+            parsed = GPR.from_string(rule)
+    finally:
+        gene_logger.disabled = was_disabled
+    if parsed.body is None:
+        raise ValueError(f"Malformed grRule: {rule!r}")
+    return parsed
+
+
 @dataclass(frozen=True)
 class GPRIssue:
     """A reaction whose GPR is flagged by the linter.

tests/test_init_genes.py

@@ -1,7 +1,8 @@
 """remove_low_score_genes — the three RAVEN docstring oracle cases.
 
-Scores use distinct values to avoid the random tie-break RAVEN mentions when all
-isozyme alternatives are negative.
+Scores use distinct values in the oracle cases. Ties are tested separately: when
+all isozyme alternatives are negative and equal, the one listed first in the rule is
+kept, on both sides (RAVEN's removeLowScoreGenes does the same).
 """
 import cobra
 
@@ -71,6 +72,22 @@ def test_removed_genes_reported_and_pruned():
     assert "G1" not in {g.id for g in out.genes}
 
 
+def test_tie_between_negative_isozymes_keeps_the_first_in_the_rule():
+    # G2 comes first in the rule although G1 comes first alphabetically.
+    assert _result("G2 or G1", {"G1": -1.0, "G2": -1.0}) == "G2"
+    assert _result("G4 or G3 or G2", {"G2": -1.0, "G3": 1.0, "G4": -1.0}) == "G3"
+
+
+def test_tie_between_negative_complexes_keeps_the_first_in_the_rule():
+    scores = {"G1": -1.0, "G2": -1.0, "G3": -1.0, "G4": -1.0}
+    assert _result("(G3 and G4) or (G1 and G2)", scores) == _norm("G3 and G4")
+
+
+def test_positive_isozymes_keep_rule_order():
+    scores = {"G1": 1.0, "G2": 1.0, "G3": 1.0}
+    assert _result("G3 or G1 or G2", scores) == _norm("G3 or G1 or G2")
+
+
 def test_tie_among_negative_isozymes_reported():
     """All four isozymes negative with the same score: the first is kept and the
     reaction is reported."""

tests/test_io_yaml.py

@@ -392,6 +392,27 @@ def test_sbo_term_written_as_a_list_is_read_as_a_string(tmp_path):
     assert model.reactions.get_by_id("r1").annotation["sbo"] == "SBO:0000627"
 
 
+def test_write_warns_on_empty_metabolite_names(tmp_path):
+    """RAVEN's writeYAMLmodel refuses empty metabolite names; Python writes but warns."""
+    m = cobra.Model("m")
+    m.add_metabolites(
+        [
+            cobra.Metabolite("a_c", name="A", compartment="c"),
+            cobra.Metabolite("b_c", name="", compartment="c"),
+        ]
+    )
+    with pytest.warns(UserWarning, match=r"1 metabolite\(s\) have an empty name \(b_c\)"):
+        write_yaml_model(m, tmp_path / "out.yml")
+    assert (tmp_path / "out.yml").exists()
+
+
+def test_write_is_silent_when_all_metabolites_are_named(tmp_path, recwarn):
+    m = cobra.Model("m")
+    m.add_metabolites([cobra.Metabolite("a_c", name="A", compartment="c")])
+    write_yaml_model(m, tmp_path / "out.yml")
+    assert not [w for w in recwarn if "empty name" in str(w.message)]
+
+
 def test_write_quotes_yaml_1_1_booleans(tmp_path):
     """Nitric oxide (name and formula NO), formula Y and compartment n are
     YAML 1.1 booleans, which PyYAML reads as True or False when unquoted."""

tests/test_tasks.py

@@ -211,3 +211,82 @@ def test_check_tasks_accepts_a_file_path(model, tmp_path):
     )
     results = check_tasks(model, p)  # path, parsed internally
     assert _by_id(results)["make_B"].passed
+
+
+# --------------------------------------------------------------------------- #
+# Parity with RAVEN parseTaskList / checkTasks (flags, lists, CHANGED, duplicates)
+# --------------------------------------------------------------------------- #
+def test_parse_flags_false_texts(tmp_path):
+    p = tmp_path / "flags.txt"
+    p.write_text(
+        "ID\tSHOULD FAIL\tPRINT FLUX\n"
+        "A\tfalse\t0\n"
+        "B\tNo\tFALSE\n"
+        "C\tyes\t1\n"
+        "D\t\tTRUE\n"
+    )
+    flags = {t.id: (t.should_fail, t.print_fluxes) for t in parse_task_list(p)}
+    assert flags == {"A": (False, False), "B": (False, False), "C": (True, True), "D": (False, True)}
+
+
+def test_parse_semicolon_lists_are_trimmed_and_drop_empty_items(tmp_path):
+    p = tmp_path / "lists.txt"
+    p.write_text(
+        "ID\tIN\tIN UB\tCHANGED RXN\tCHANGED UB\n"
+        "T1\tA[c]; B[c];\t5\t;r1 ; r2\t10\n"
+    )
+    (t,) = parse_task_list(p)
+    assert t.inputs == [("A[c]", 0.0, 5.0), ("B[c]", 0.0, 5.0)]
+    assert [c[0] for c in t.changed] == ["r1", "r2"]
+
+
+def test_parse_empty_changed_bounds_are_none(tmp_path):
+    p = tmp_path / "changed.txt"
+    p.write_text("ID\tCHANGED RXN\tCHANGED LB\tCHANGED UB\nT1\tr1\t\t10\nT2\tr1\t-5\t\nT3\tr1\t1\t2\n")
+    t1, t2, t3 = parse_task_list(p)
+    assert t1.changed == [("r1", None, 10.0)]
+    assert t2.changed == [("r1", -5.0, None)]
+    assert t3.changed == [("r1", 1.0, 2.0)]
+
+
+def test_parse_headers_ignore_case_and_whitespace(tmp_path):
+    p = tmp_path / "hdr.txt"
+    p.write_text(" id \tIn\t in ub \nT1\tA[c]\t7\n")
+    (t,) = parse_task_list(p)
+    assert t.inputs == [("A[c]", 0.0, 7.0)]
+
+
+def test_parse_xls_gives_a_clear_error(tmp_path):
+    p = tmp_path / "old.xls"
+    p.write_bytes(b"")
+    with pytest.raises(ValueError, match=r"\.xls.*not supported.*\.xlsx"):
+        parse_task_list(p)
+
+
+def test_empty_changed_bound_leaves_reaction_bound_alone(model):
+    """A None bound keeps the reaction's own bound; r1 is [0, 1000]."""
+    inputs, outputs = [("A[c]", 0, 1000)], [("B[c]", 1, 1000)]
+    keep_lb = Task("keep_lb", inputs=inputs, outputs=outputs, changed=[("r1", None, 1000.0)])
+    (res,) = check_tasks(model, [keep_lb])
+    assert res.passed
+    # UB 0 with the LB left at 0 blocks r1, so B cannot be made.
+    ub_zero = Task("ub_zero", inputs=inputs, outputs=outputs, changed=[("r1", None, 0.0)])
+    (res2,) = check_tasks(model, [ub_zero])
+    assert not res2.passed
+    # The model passed in is not modified.
+    assert model.reactions.r1.bounds == (0, 1000)
+
+
+def test_duplicate_name_comp_warns_and_constrains_all():
+    m = cobra.Model("dup")
+    a1, a2, b = _met("A1_c", "A"), _met("A2_c", "A"), _met("B_c", "B")
+    m.add_metabolites([a1, a2, b])
+    r1 = cobra.Reaction("r1", lower_bound=0, upper_bound=1000)
+    r1.add_metabolites({a1: -1, b: 1})
+    r2 = cobra.Reaction("r2", lower_bound=0, upper_bound=1000)
+    r2.add_metabolites({a2: -1, b: 1})
+    m.add_reactions([r1, r2])
+    task = Task("dup", inputs=[("A[c]", 0, 1000)], outputs=[("B[c]", 1, 1000)])
+    with pytest.warns(UserWarning, match=r"matches 2 model metabolites.*RAVEN constrains only the first"):
+        (res,) = check_tasks(m, [task])
+    assert res.passed

tests/test_tasks_essential.py

@@ -5,6 +5,7 @@
 e[c]->e[s] producer); the alternative internal paths make nothing else essential.
 """
 import cobra
+import pytest
 from tinit_oracles import (
     TEST_MODEL_TASK_ESSENTIAL_PREMERGE,
     make_test_model,
@@ -195,6 +196,7 @@ def test_duplicate_name_comp_metabolites_both_constrained():
     m.objective = "R1"
     # Output b from input a -> 'a[s]' matches both a1 and a2; must not be reported missing.
     task = Task(id="t", inputs=[("a[s]", 0.0, 1000.0)], outputs=[("b[s]", 1.0, 1.0)])
-    res = find_task_essential_reactions(m, [task])
+    with pytest.warns(UserWarning, match="matches 2 model metabolites"):
+        res = find_task_essential_reactions(m, [task])
     assert res.failed_tasks == []  # 'a[s]' resolved (to both a1 and a2), task feasible
     assert "R1" in res.reactions

tests/test_utils_gpr.py

@@ -82,3 +82,18 @@ def test_find_non_dnf_grrules_flags_only_offenders():
 def test_find_non_dnf_grrules_empty_when_all_clean():
     model = _model_with_rules({"R1": "G1 or G2", "R2": "(G1 and G2) or G3"})
     assert find_non_dnf_grrules(model) == []
+
+
+@pytest.mark.parametrize(
+    "rule",
+    ["G1 and", "(G1 or G2", "G1 or or G2", "G1 G2", "and G1", "(G1 and G2)) or G3", "()"],
+)
+def test_is_dnf_raises_on_malformed_rule(rule):
+    """RAVEN's isDnfGrRule errors on these (RAVEN:badGrRule); cobra alone returns an empty GPR."""
+    with pytest.raises(ValueError, match="Malformed grRule"):
+        is_dnf(rule)
+
+
+@pytest.mark.parametrize("rule", ["", "   ", None])
+def test_is_dnf_blank_rule_is_still_true(rule):
+    assert is_dnf(rule) is True