diff --git a/CHANGELOG.md b/CHANGELOG.md index 228aa58..c91849a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,7 @@ ## Unreleased +* `pfba_enzymes` with `fraction_of_optimum` below 1 lets the objective fall by that fraction of the magnitude of its optimum, for a negative optimum and for a minimised objective as for a maximised positive one, as MATLAB's `getPfbaEnzymes` does; it raised `Infeasible` for these objectives. * `ec_fseof` ranks the targets by the score `|slope| * |correlation|` rounded to 9 significant digits, with scores that are equal after rounding in model order, as GECKO's `ecFSEOF` does; scores that differ by solver noise used to decide the order. * `geckopy.kcat_tuning.evotune` exports `best_parsimonious`, `n_changed`, `parsimony_frontier` and `isozyme_tie_map`, the helpers behind the parsimony sweep and the tying of isozyme copies, which GECKO has as `bestParsimonious`, `nChanged`, `parsimonyFrontier` and `isozymeTieMap`. * `read_kcat_source` and `load_kcat_source` skip a row whose fields are all empty quoted fields (`""""`), as they skip a row of tabs, and raise `ValueError` for a field that opens a double quote and does not close it. Before, the first was reported as an invalid row and the second was read as a row with a kcat that is not a number. diff --git a/src/geckopy/utilities/pfba_enzymes.py b/src/geckopy/utilities/pfba_enzymes.py index 66be1c9..44adbd9 100644 --- a/src/geckopy/utilities/pfba_enzymes.py +++ b/src/geckopy/utilities/pfba_enzymes.py @@ -48,9 +48,9 @@ def pfba_enzymes( ) -> cobra.Solution: """Parsimonious FBA minimising the L1 norm of ``usage_prot_*`` fluxes. - Sets the objective (if given), fixes it at - ``fraction_of_optimum * 0.999999 * optimal_value``, then minimises - ``sum(usage_prot_*.forward_variable)``. Usage reactions are + Sets the objective (if given), lets it fall from its optimum by + ``1 - fraction_of_optimum * 0.999999`` times the magnitude of the + optimum, then minimises ``sum(usage_prot_*.forward_variable)``. Usage reactions are forward-only in the GECKO 3 layout (``lb=0``); if a user has flipped one to allow reverse flux, its ``reverse_variable`` is included automatically. @@ -69,13 +69,16 @@ def pfba_enzymes( enzyme usage. If ``None``, the model's current objective is used as-is. fraction_of_optimum - Fraction of the objective's optimal value to fix as a + Fraction of the objective's optimal value to keep as a constraint before minimising enzyme usage (default 1.0). The fraction is multiplied by 0.999999, the margin of RAVEN's - ``solveLP`` that MATLAB's ``getPfbaEnzymes`` applies, so the - default fixes the objective at 0.999999 of its optimum and - 0.9 fixes it at 0.8999991. A lower value trades some - objective value for potentially lower enzyme usage. + ``solveLP`` that MATLAB's ``getPfbaEnzymes`` applies. The + objective may fall from its optimum by ``1 - fraction * 0.999999`` + times the magnitude of the optimum, so the default keeps it within + 1e-6 of the optimum and 0.9 within 10 %, for a negative optimum + (the exchange flux of a nutrient that is maximised) as for a + positive one. A lower value trades some objective value for + potentially lower enzyme usage. Returns ------- @@ -108,9 +111,16 @@ def pfba_enzymes( f"Objective is already {_PFBA_OBJ_NAME}; " "pfba_enzymes appears to be nested." ) - cobra.util.fix_objective_as_constraint( - m, fraction=fraction_of_optimum * _PFBA_MARGIN, - ) + optimum = m.slim_optimize(error_value=None) + if optimum is None: + raise cobra.exceptions.OptimizationError( + "pfba_enzymes: the model is infeasible for its objective." + ) + # The objective may fall by (1 - fraction * margin) of the magnitude + # of its optimum, for a positive and a negative optimum alike. + slack = (1.0 - fraction_of_optimum * _PFBA_MARGIN) * abs(optimum) + bound = optimum - slack if m.objective.direction == "max" else optimum + slack + cobra.util.fix_objective_as_constraint(m, bound=bound) usage_rxns = [r for r in m.reactions if r.id.startswith(USAGE_PREFIX)] variables = [] diff --git a/tests/test_pfba_enzymes.py b/tests/test_pfba_enzymes.py index 2259839..868be2c 100644 --- a/tests/test_pfba_enzymes.py +++ b/tests/test_pfba_enzymes.py @@ -112,3 +112,26 @@ def test_gecko_light_raises(): model.ec.gecko_light = True with pytest.raises(NotImplementedError): pfba_enzymes(model) + + +def test_fraction_of_a_negative_optimum_lets_the_objective_fall(): + model = _ectestgem_ec_model_with_kcats() + # Maximising -R3 with R3 >= 1 gives an optimum of -1. + model.reactions.get_by_id("R3").lower_bound = 1.0 + model.objective = {model.reactions.get_by_id("R3"): -1.0} + optimum = model.slim_optimize() + assert optimum == pytest.approx(-1.0) + sol = pfba_enzymes(model, fraction_of_optimum=0.9) + assert sol.status == "optimal" + assert -sol.fluxes["R3"] >= optimum * 1.1 - 1e-6 + + +def test_fraction_of_a_minimised_objective_lets_it_rise(): + model = _ectestgem_ec_model_with_kcats() + model.reactions.get_by_id("R3").lower_bound = 1.0 + model.objective = "R3" + model.objective_direction = "min" + optimum = model.slim_optimize() + sol = pfba_enzymes(model, fraction_of_optimum=0.9) + assert sol.status == "optimal" + assert sol.fluxes["R3"] <= optimum * 1.1 + 1e-6