Merge pull request #41 from COBREXA/mk-coupling

add support for loading the coupling bounds from FBC

.github/workflows/ci.yml

@@ -29,7 +29,7 @@ jobs:
         with:
           version: ${{ matrix.version }}
           arch: ${{ matrix.arch }}
-      - uses: julia-actions/cache@v1
+      - uses: julia-actions/cache@v2
       - uses: julia-actions/julia-buildpkg@latest
       - run: |
           git config --global user.name Tester

.github/workflows/docs.yml

@@ -19,7 +19,7 @@ jobs:
       - uses: julia-actions/setup-julia@latest
         with:
           version: '1.9'
-      - uses: julia-actions/cache@v1
+      - uses: julia-actions/cache@v2
       - name: Install dependencies
         run: julia --color=yes --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()'
       - name: Build and deploy

Project.toml

@@ -17,7 +17,7 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 
 [compat]
-AbstractFBCModels = "0.2.2"
+AbstractFBCModels = "0.3"
 Aqua = "0.7"
 Clarabel = "0.6"
 ConstraintTrees = "1.1"

README.md

@@ -11,7 +11,7 @@
 [docs-img-dev]: https://img.shields.io/badge/docs-latest-0af.svg
 [docs-url-dev]: https://cobrexa.github.io/COBREXA.jl/dev/
 [docs-url-examples]: https://cobrexa.github.io/COBREXA.jl/dev/examples/
-[docs-url-fba]: https://cobrexa.github.io/COBREXA.jl/dev/examples/02-flux-balance-analysis/
+[docs-url-fba]: https://cobrexa.github.io/COBREXA.jl/dev/examples/02a-flux-balance-analysis/
 
 [docker-url]: https://github.com/COBREXA/COBREXA.jl/pkgs/container/cobrexa.jl
 [docker-img]: https://ghcr-badge.egpl.dev/cobrexa/cobrexa.jl/size?color=%2362a0ea&tag=latest&label=docker&trim=

docs/src/examples/02b-model-modifications.jl → docs/src/examples/02c-model-modifications.jl

@@ -37,7 +37,7 @@
 # disadvantage is that the "common" FBC model interface does not easily express
 # various complicated constructions (communities, reaction coupling, enzyme
 # constraints, etc.) -- see the [example about modifying the
-# constraints](02c-constraint-modifications.md) for more details.
+# constraints](02d-constraint-modifications.md) for more details.
 #
 # ## Getting the base model
 
@@ -166,4 +166,28 @@ flux_changes =
 sort(collect(flux_changes), by = last)
 
 #md # !!! tip "Always use a uniquely defined flux solutions for flux comparisons"
-#md #     Since the usual flux balance allows a lot of freedom in the "solved" flux and the only value that is "reproducible" by the analysis is the objective, one should never compare the flux distributions directly. Typically, that may result in false-positive (and sometimes false-negative) differences. Use e.g. [parsimonious FBA](03-parsimonious-flux-balance.md) to obtain uniquely determined and safely comparable flux solutions.
+#md #     Since the usual flux balance allows a lot of freedom in the "solved" flux and the only value that is "reproducible" by the analysis is the objective, one should never compare the flux distributions directly. Typically, that may result in false-positive (and sometimes false-negative) differences. Use e.g. [parsimonious FBA](03b-parsimonious-flux-balance.md) to obtain uniquely determined and safely comparable flux solutions.
+
+# ## Coupling constraints
+#
+# Some model types support additional constraints over the reaction fluxes,
+# which are historically called "coupling". These allow to e.g. place a bound
+# on a total flux through several reactions.
+#
+# Canonical model supports these as "couplings":
+
+model.couplings["total_energy_intake"] = CM.Coupling(
+    lower_bound = 0,
+    upper_bound = 5,
+    reaction_weights = Dict("EX_glc__D_e" => -1.0, "EX_fru_e" => -1.0, "EX_pyr_e" => -1.0),
+)
+
+# The values of any coupling constraints can be inspected directly in the
+# solved model:
+
+solution_with_coupling = flux_balance_analysis(model, optimizer = GLPK.Optimizer)
+
+solution_with_coupling.coupling.total_energy_intake
+
+@test 0 <= solution_with_coupling.coupling.total_energy_intake <= 5 #src
+@test isapprox(solution_with_coupling.objective, 0.4155977750928965, atol = TEST_TOLERANCE) #src

docs/src/examples/02c-constraint-modifications.jl → docs/src/examples/02d-constraint-modifications.jl

@@ -17,7 +17,7 @@
 # # Making adjustments to the constraint system
 #
 # In the [previous example about model
-# adjustments](02b-model-modifications.md), we noted that some constraint
+# adjustments](02c-model-modifications.md), we noted that some constraint
 # systems may be too complex to be changed within the limits of the usual FBC
 # model view, and we may require a sharper tool to do the changes we need. This
 # example shows how to do that by modifying the constraint systems that are

src/frontend/balance.jl

@@ -52,6 +52,9 @@ function flux_balance_constraints(
     stoi = A.stoichiometry(model)
     bal = A.balance(model)
     obj = A.objective(model)
+    cpls = Symbol.(A.couplings(model))
+    cT = SparseArrays.SparseMatrixCSC(A.coupling(model)')
+    clbs, cubs = A.coupling_bounds(model)
 
     # The iteration through stoichiometry would be better done with eachrow(),
     # unfortunately it seems to be enormously inefficient on column-major
@@ -70,6 +73,14 @@ function flux_balance_constraints(
                 bound = C.EqualTo(b),
             ) for (met, row_idx, b) in zip(mets, 1:stoiT.n, bal)
         ) *
+        :coupling^C.ConstraintTree(
+            cpl => C.Constraint(
+                value = let i = cT.colptr[row_idx], e = cT.colptr[row_idx+1] - 1
+                    C.LinearValue(idxs = cT.rowval[i:e], weights = cT.nzval[i:e])
+                end,
+                bound = clb == cub ? C.EqualTo(clb) : C.Between(clb, cub),
+            ) for (cpl, row_idx, clb, cub) in zip(cpls, 1:cT.n, clbs, cubs)
+        ) *
         :objective^C.Constraint(C.LinearValue(SparseArrays.sparse(obj))),
     )