cobrapy
cobrapy copied to clipboard
Exchanges and external compartment don't match in models with SBO terms
Problem description
If a model annotates reactions as exchanges those are picked up by model.exchanges, however, the external compartment detection does not take this into account. So the model will return exchanges that are not in the external compartment.
Code Sample
This uses Recon 2.2 from https://github.com/u003f/recon2/tree/master/models.
In [23]: from cobra.io.sbml import read_sbml_model, LOGGER
In [24]: LOGGER.setLevel("ERROR")
In [25]: recon2 = read_sbml_model("recon_2.2.xml")
In [26]: recon2.exchanges[0]
Out[26]: <Reaction EX_10fthf_b at 0x7f2233f1e5b0>
In [27]: recon2.exchanges[0].compartments
Out[27]: {'b'}
In [28]: from cobra.medium import find_external_compartment
In [29]: find_external_compartment(recon2)
Out[29]: 'e'
It's okay that we can sometimes not identify the external compartment, especially since Recon2 uses a non-standard identifier but it has to be consistent. Either return exchanges from the compartment "e" or report the external compartment as "b".
Context
System Information
==================
OS Linux
OS-release 5.10.0-4-amd64
Python 3.9.2
Package Versions
================
appdirs~=1.4 not installed
black 20.8b1
bumpversion 0.6.0
cobra 0.21.0
depinfo 1.6.0
diskcache~=5.0 not installed
future 0.18.2
httpx~=0.14 not installed
importlib_resources 3.3.0
isort; extra == "development" not installed
numpy~=1.13 not installed
optlang<1.4.6 not installed
pandas~=1.0 not installed
pip 21.0.1
pydantic~=1.6 not installed
python-libsbml==5.19.0 not installed
rich~=6.0 not installed
ruamel.yaml~=0.16 not installed
scipy 1.6.1
setuptools 54.0.0
six 1.15.0
swiglpk 5.0.0
tox; extra == "development" not installed
wheel 0.36.2
Hmm, this is a bit complicated and is related to boundary metabolites yet again. So the problem is reactions of the form R_e: A_e <->
where A_e
is marked at a boundary metabolite in SBML but the R_e
reaction is annotated with BO:0000627
which is an exchange reaction. Parsing the model will yield to adding in the exchange for the boundary metabolite thus yielding a new reaction R_b: A_b <->
and modifying R_e: A_e <-> A_b
. We still detect R_e
as the exchange due to the SBO term which has priority in the code. However, R_e
is not a valid exchange reaction anymore and R_b
has no SBO term. This can create all kinds of problems. @matthiaskoenig would love to hear your opinion. In my view, this is a slightly incorrect model specification but it won't be caught by the validators. There are several potential ways to fix that each with its own advantages/disadvantages:
- remove the SBO term on the "wrong" exchanges
- add the SBO term to the newly added boundary reactions
R_b
- ignore SBO terms for exchanges and only use the heuristic
I think another good helper function for cobrapy would be a reaction to prune linear chains of boundary reactions to a single reaction so <-> A_b <-> A_c <-> A_e
becomes <-> A_e
. This is favorable for model interpretation but also for computation as it simplifies the stoichiometry.,
I think another good helper function for cobrapy would be a reaction to prune linear chains of boundary reactions to a single reaction so
<-> A_b <-> A_c <-> A_e
becomes<-> A_e
. This is favorable for model interpretation but also for computation as it simplifies the stoichiometry.,
Side note: I would only do that if there are no gene rules associated with such reactions.
Good point, this could be an argument to that function.