Merge pull request #1341 from cdiener/fix/cycle_free

fix the objective in CycleFreeFlux

documentation_builder/plot_helper.py

@@ -1,6 +1,6 @@
 # -*- coding: utf-8 -*-
 from matplotlib.pyplot import figure, xlim, ylim, gca, arrow, text, scatter
-from mpl_toolkits.axes_grid.axislines import SubplotZero
+from mpl_toolkits.axisartist.axislines import SubplotZero
 from numpy import linspace, arange, sqrt, pi, sin, cos, sign
 from IPython.display import set_matplotlib_formats
 

documentation_builder/requirements.txt

@@ -3,3 +3,4 @@ sphinxcontrib-napoleon
 sphinx-autoapi
 nbsphinx>=0.2.4
 ipykernel
+matplotlib>=3.6

release-notes/next-release.md

@@ -9,6 +9,9 @@ COBRApy is now compatible with pydantic version 2.
 COBRApy is now compatible with the latest numpy versions and pin to a lower version
 has been removed.
 
+`loopless_solution` now fixes the objective to its optimum as in the
+originally published method and returns the objective value in the solution object.
+
 ## Other
 
 Backwards compatibility for pickled models has been improved.

src/cobra/flux_analysis/loopless.py

@@ -154,14 +154,18 @@ def loopless_solution(
     The returned flux solution has the following properties:
 
     - It contains the minimal number of loops possible and no loops at all
-      if all flux bounds include zero.
-    - It has an objective value close to the original one and the same
-      objective value id the objective expression can not form a cycle
+      if all flux bounds include zero and the objective is not in a cycle.
+    - It has the same objective value as the original flux solution and assumes
+      that the objective does not participate in a cycle
       (which is usually true since it consumes metabolites).
     - It has the same exact exchange fluxes as the previous solution.
     - All fluxes have the same sign (flow in the same direction) as the
       previous solution.
 
+    When providing fluxes to the method, please note that those have to come from the
+    exact same model that you provided, meaning that no bounds or coefficients have
+    been changed, and the optimum has remained the same.
+
     References
     ----------
     .. [1] CycleFreeFlux: efficient removal of thermodynamically infeasible
@@ -176,13 +180,16 @@ def loopless_solution(
     if fluxes is None:
         sol = model.optimize(objective_sense=None)
         fluxes = sol.fluxes
+        opt = sol.objective_value
+    else:
+        opt = model.slim_optimize()
 
     with model:
         prob = model.problem
-        # Needs one fixed bound for cplex...
+        # Fix the objective
         loopless_obj_constraint = prob.Constraint(
             model.objective.expression,
-            lb=-1e32,
+            lb=opt,
             name="loopless_obj_constraint",
         )
         model.add_cons_vars([loopless_obj_constraint])

tests/test_core/test_core_reaction.py

@@ -1036,7 +1036,7 @@ def test_compartment_changes(model: Model) -> None:
 def test_gpr_serialization(model: Model) -> None:
     """Verify that reactions GPRs are serialized compactly as str."""
     state = model.reactions[0].__getstate__()
-    assert type(state["_gpr"]) == str
+    assert isinstance(state["_gpr"], str)
 
 
 def test_gpr_serialization_backwards_compatibility(model: Model) -> None:

tests/test_flux_analysis/test_loopless.py

@@ -62,19 +62,23 @@ def test_loopless_benchmark_after(benchmark: Callable) -> None:
 
 def test_loopless_solution(ll_test_model: Model) -> None:
     """Test loopless_solution()."""
+    opt_feasible = ll_test_model.slim_optimize()
     solution_feasible = loopless_solution(ll_test_model)
     ll_test_model.reactions.v3.lower_bound = 1
-    ll_test_model.optimize()
+    opt_infeasible = ll_test_model.slim_optimize()
     solution_infeasible = loopless_solution(ll_test_model)
     assert solution_feasible.fluxes["v3"] == 0.0
+    assert solution_feasible.objective_value == pytest.approx(opt_feasible)
     assert solution_infeasible.fluxes["v3"] == 1.0
+    assert solution_infeasible.objective_value == pytest.approx(opt_infeasible)
 
 
 def test_loopless_solution_fluxes(model: Model) -> None:
     """Test fluxes of loopless_solution()."""
-    fluxes = model.optimize().fluxes
-    ll_solution = loopless_solution(model, fluxes=fluxes)
+    sol = model.optimize()
+    ll_solution = loopless_solution(model, fluxes=sol.fluxes)
     assert len(ll_solution.fluxes) == len(model.reactions)
+    assert ll_solution.objective_value == pytest.approx(sol.objective_value)
 
 
 def test_add_loopless(ll_test_model: Model) -> None: