From 2da2af0ab76fc8befd8e6dbb1a0a685d4ea67558 Mon Sep 17 00:00:00 2001 From: Eduard Kerkhoven Date: Thu, 5 Nov 2026 12:44:27 +0100 Subject: [PATCH] Drop the complexes without an ortholog under the drop policy The policy drop removes each complex that has a subunit without an ortholog and removes the reaction when no complex remains. It removed the reaction as soon as one subunit had no ortholog. --- CHANGELOG.md | 4 +++ .../reconstruction/homology/homology.py | 34 +++++++++---------- tests/test_reconstruction_homology.py | 25 ++++++++++++-- 3 files changed, 43 insertions(+), 20 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index a975bb1..448ee4c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,10 @@ ## Unreleased +### Breaking changes + +* `get_model_from_homology` with `complex_policy="drop"` removes each complex that has a subunit without an ortholog, and removes the reaction only when no complex remains (was: the reaction was removed as soon as one subunit had no ortholog). For the rule `(A and B) or (C and D)` without an ortholog of `B`, the draft has `C and D`. + ### Data * The KEGG assets are the `kegg120` release of `toolbox-assets` (was: `kegg118`): the reference model and tables, the taxonomy and the profile HMM libraries of 4246 prokaryotic and 2724 eukaryotic KOs. The baked registry and `data/manifest.json` list them. diff --git a/src/raven_toolbox/reconstruction/homology/homology.py b/src/raven_toolbox/reconstruction/homology/homology.py index 2239f16..6e631ae 100644 --- a/src/raven_toolbox/reconstruction/homology/homology.py +++ b/src/raven_toolbox/reconstruction/homology/homology.py @@ -91,27 +91,23 @@ def _rewrite_node(node, ortho: dict): return _group("and", children) -def _contains_unmapped(expr) -> bool: - if isinstance(expr, _Unmapped): - return True - return isinstance(expr, tuple) and any(_contains_unmapped(c) for c in expr[1]) - - def _apply_policy(expr, policy: str, model_id: str): """Resolve the unmapped leaves of a rewritten tree under ``complex_policy``. - ``flag`` turns each into an ``OLD__`` placeholder gene, ``keep`` removes - them (a group left empty disappears with them) and ``drop`` removes the whole reaction - as soon as one is left. Returns ``None`` when nothing remains. + ``flag`` turns each into an ``OLD__`` placeholder gene and ``keep`` removes + them (a group left empty disappears with them). ``drop`` removes each complex that has + one: an AND group with an unmapped subunit is removed as a whole, and an OR group keeps + its remaining alternatives. Returns ``None`` when nothing remains. """ - if policy == "drop": - return None if _contains_unmapped(expr) else expr if isinstance(expr, _Unmapped): return f"OLD_{model_id}_{expr.gene}" if policy == "flag" else None if isinstance(expr, str): return expr op, children = expr - return _group(op, [c for c in (_apply_policy(c, policy, model_id) for c in children) if c is not None]) + resolved = [_apply_policy(c, policy, model_id) for c in children] + if policy == "drop" and op == "and" and any(c is None for c in resolved): + return None + return _group(op, [c for c in resolved if c is not None]) def _prune_genes(node, removed: set): @@ -334,11 +330,15 @@ def get_model_from_homology( complex_policy : str What happens to an AND-linked subunit with no ortholog. Default ``"flag"``. ``"flag"`` keeps the reaction with a placeholder gene (``OLD__``) for - later curator review, ``"keep"`` removes that subunit, ``"drop"`` drops the whole - reaction. An isozyme group counts as a subunit without an ortholog only when none of - its members has one, and then one placeholder stands for it. ``"drop"`` removes the - reaction whenever a placeholder would remain, also when another isozyme branch is - complete. Any other value raises :class:`ValueError`. + later curator review, ``"keep"`` removes that subunit, ``"drop"`` removes each complex + that has such a subunit and removes the reaction when no complex remains. An isozyme + group counts as a subunit without an ortholog only when none of its members has one, + and then one placeholder stands for it. Under ``"drop"``, the rule + ``(A and B) or (C and D)`` with no ortholog of ``B`` becomes ``C and D``, and the rule + ``(A or B) and C`` with no ortholog of ``A`` or ``B`` removes the reaction. A complex + that contains another complete complex is kept, and a gene shared by several complexes + removes only the complexes in which it has no ortholog. Any other value raises + :class:`ValueError`. only_genes_in_models : bool Restricts the hits table to hits whose template gene appears in the given ``models`` before mapping orthologs, so a gene the templates lack cannot be a best hit. It raises diff --git a/tests/test_reconstruction_homology.py b/tests/test_reconstruction_homology.py index d2d9d0e..fef710d 100644 --- a/tests/test_reconstruction_homology.py +++ b/tests/test_reconstruction_homology.py @@ -130,9 +130,28 @@ def test_flagged_subunit_stays_beside_a_complete_isozyme_branch(): assert _policy_draft("A or (B and C)", ["A", "C"], "keep") == {"R1": "nA or nC"} -def test_complex_policy_drop_removes_reaction_with_any_flagged_subunit(): - assert _policy_draft("A or (B and C)", ["A", "C"], "drop") == {} - assert _policy_draft("(A and B) or C", ["C"], "drop") == {} +def test_complex_policy_drop_removes_each_complex_with_a_flagged_subunit(): + assert _policy_draft("A or (B and C)", ["A", "C"], "drop") == {"R1": "nA"} + assert _policy_draft("(A and B) or C", ["C"], "drop") == {"R1": "nC"} + assert _policy_draft("A and B", ["A"], "drop") == {} + assert _policy_draft("(A and B) or (C and D)", ["A", "C", "D"], "drop") == {"R1": "nC and nD"} + + +def test_complex_policy_drop_keeps_complete_complexes_of_a_shared_subunit(): + rule = "(A and B) or (A and C) or (D and E) or (A and B and F)" + assert _policy_draft(rule, ["A", "B", "C", "D", "F"], "drop") == { + "R1": "(nA and nB) or (nA and nC) or (nA and nB and nF)"} + assert _policy_draft("(A and B) or B", ["B"], "drop") == {"R1": "nB"} + assert _policy_draft("(A and B) or B", ["A", "B"], "drop") == {"R1": "(nA and nB) or nB"} + + +def test_complex_policy_drop_keeps_isozyme_groups_inside_a_complex(): + hits = make_ortholog_hits([("A", "n1"), ("A", "n2"), ("B", "n3")], "T", "N") + draft = get_model_from_homology([_one_rule_template("T", "R1", "A and B")], hits, "N", complex_policy="drop") + assert _rules(draft) == {"R1": "(n1 or n2) and n3"} + + +def test_complex_policy_flag_and_keep_on_an_incomplete_complex_beside_an_isozyme(): assert _policy_draft("(A and B) or C", ["C"], "flag") == {"R1": "(OLD_T_A and OLD_T_B) or nC"} assert _policy_draft("(A and B) or C", ["C"], "keep") == {"R1": "nC"}