[jplfaria]

Credit. Andrew Freiburger found and fixed this bug in his fork:
freiburgermsu/ModelSEEDpy#19, in build_and_gapfill, his
fork's name for dev's build(). This PR carries his fix onto dev. It uses his helper, apply_gapfill_reversals,
with the same semantics. His PR also gives the effect on a 200-genome benchmark: 107 of 2,362 gap-filled models needed
a reversal, and in an earlier run that dropped reversals, 136 models could not grow.

Why this matters in our workflow

We build models with MSBuilder.build(). When its gap-fill needs a reaction already in the model to run in a
direction the model blocks, build() returns a model without that change, and the model may not grow on the medium
it was gap-filled on. Our glucose builds have not needed a reversal so far, but in Andrew's benchmark (gap-filling on
complete medium) 107 of 2,362 gap-filled models, about 5%, needed a reversal.

What happens

msbuilder.py L1041-L1053 at 062714b:

gapfill_res = gapfill.run_gapfilling(gapfill_media)
...
gap_sol = {}
for rxn_id, d in gapfill_res["new"].items():
    if rxn_id[:-1] in self.template.reactions:
        gap_sol[rxn_id[:-1]] = get_reaction_constraints_from_direction(d)
model_gapfilled = model_base.copy()
MSBuilder.integrate_gapfill_solution(self.template, model_gapfilled, gap_sol)
return model_gapfilled, atp_correction, tests, gap_sol

run_gapfilling returns "new" (template reactions to add) and "reversed" (model reactions whose needed direction
the model blocks). For a model reaction whose template reaction has GapfillDirection =, the gap-fill package
widens the model's bounds to the template reaction's bounds
(gapfillingpkg.py L434-L449).
So a reversal happens when the model's bounds are narrower than the template's: for example, a reaction narrowed in the
model after it was built, or a model built from another template. build() reads only "new", and
integrate_gapfill_solution (L962-L976)
only adds reactions, so every reversal is dropped. The package's other integration paths do apply reversals:
MSGapfill.integrate_gapfill_solution and MSBuilder.gapfill_model, both to ±100.

How to reproduce

The reproduction uses the package's e_coli_core test model. ATPM_c0's lower bound is set to 0 (as
MSBuilder.add_atpm makes it), and GAPD_c0 is made irreversible the wrong way, (-1000, 0). In
template_core_bigg.json, GAPD_c is reversible, with GapfillDirection =. Only build_base_model, which needs
an annotated genome, and load_default_medias, whose media use ModelSEED compound ids, are replaced. Gap-fill runs
for real; ATP correction runs with no ATP media.

model = get_model()                       # e_coli_core fixture from tests/core/test_msgapfill.py
model.reactions.BIOMASS_Ecoli_core_w_GAM_c0.id = "bio1"; model.repair()
model.reactions.ATPM_c0.lower_bound = 0
model.reactions.GAPD_c0.bounds = (-1000, 0)
MSBuilder.build_base_model = lambda self, *args: model
modelseedpy.core.msatpcorrection.load_default_medias = lambda: []
model_gapfilled, _, _, gap_sol = MSBuilder(MSGenome(), template, template_core=template).build(
    "test", media_glucose_aerobic, annotate_with_rast=False)

	gap-fill result	gap_sol	GAPD_c0 in returned model	growth on the gap-fill medium
dev	reversed {'GAPD_c0': '>'}	{}	(-1000, 0)	0.0
this branch	reversed {'GAPD_c0': '>'}	{}	(-1000, 1000)	0.9166

The fix

• New static MSBuilder.apply_gapfill_reversals(model, reversed_reactions) (Andrew's helper). For each reversed
  reaction, it opens the gap-filled direction to the bound get_reaction_constraints_from_direction gives, which is
  ±1000, the bound build() already gives new reactions. The other bound is kept, and so is a bound that is already
  wider. A reaction that is not in the model is skipped with a warning. The helper returns the ids it opened.
• build() calls it right after integrate_gapfill_solution. The return value of build() is unchanged.
• Differences from Andrew's version. It is called from dev's build(). It iterates the dict as given instead of
  sorted(...). It does not accept None: on dev, gapfill_res["reversed"] is always a dict, because when
  run_gapfilling returns None, build() already fails on gapfill_res["new"]. It has no debug log line.
• ATP safety. run_test_conditions already evaluates the solution with the reversed direction open, up to the
  template's bound. So opening it in the final model does not admit anything the tests did not see.

Two package functions call build(): build_from_annotaton_ontology (L791) and build_metabolic_model (L1235). Both
use the pre-refactor signature and fail before reaching this code (separate PR). External callers get the same return
value, and the returned model now includes the reversals.

Tests

• tests/core/test_msbuilder.py::test_apply_gapfill_reversals covers the helper on a small model: a reverse opening,
  a forward opening with the other bound kept, a wider bound kept, and an absent reaction skipped.
• tests/core/test_msgapfill.py::test_build_applies_gapfill_reversals is the reproduction above, using that file's
  fixtures.

Against dev's code, both fail:

E       assert (-1000, 0) == (-1000, 1000)
tests/core/test_msgapfill.py:118: AssertionError
E       AttributeError: type object 'MSBuilder' has no attribute 'apply_gapfill_reversals'

With the fix, both pass. Full suite (Python 3.10, cobra 0.32.1, CPLEX 22.1.1): tests/core goes from 56 to 58
passed, and tests/ from 58 to 60. black 22.8.0 and the other pre-commit hooks pass.

Real-model check

We ran MSBuilder.build() on two genome-scale builds (E. coli K-12 and a Pantoea strain, from RAST annotations,
GramNeg template v7.1 with Core-V6.1, minimal glucose), once on dev and once on this branch. Neither gap-fill returned
a reversal, and the two results are identical: final models, every gap-fill result, ATP correction outputs and growth.

Related

• Issue: MSBuilder.build() drops the gap-filler's reversed reactions, so the model can fail to grow on its own gap-fill medium jplfaria/ModelSEEDpy#3
• Andrew Freiburger's fix: build_and_gapfill: apply the gapfiller's reversals, not only its new reactions freiburgermsu/ModelSEEDpy#19

🤖 Generated with Claude Code