[API] How should SBML reaction compartments be handled?

[cdiener]

Moved here from the PR discussion:

MK:
For the transport reactions you just add a compartment such as id="pm" name="plasma membrane" spatialDimensions="2" and then set transport_reaction.compartment=pm. This also allows nicely to define at which interface compartments transporters are located. The spatialDimensions information is especially important in scaling of systems, because volume and area scale differently, e.g. if a cell growths. So you would scale proteins differently if a cell growths.

I agree, copying of the complete compartment is nice to have, but should not be the default behavior.

Dev meeting would be great. Would be happy to participate.

Sorry forgot. There are also many membrane bound reactions which catalyze on one side of the compartment. I.e. the enymes are associated with the membrane, but all the metabolites involved are on one side of the membrane. Inferring compartments from metabolites gives incorrect compartments here.

MB:

Sorry forgot. There are also many membrane bound reactions which catalyze on one side of the compartment. I.e. the enymes are associated with the membrane, but all the metabolites involved are on one side of the membrane. Inferring compartments from metabolites gives incorrect compartments here.

I don't see yet how this gives the wrong information. If an enzyme is membrane bound but carries out the reaction only in one compartment, that is where the reaction happens, right?

CD:
Agree with @Midnighter here. Reaction is a bit of an abstract concept but if you go by physical chemistry it is a collision that passes an energy threshold so it requires physical contact. So the actual reaction would still take place where the substrate is located. For channels the initial adsorption/entry takes place on one side of the membrane and secretion/export on the other side. The substrate and product would have different compartments and the reaction would have two compartments due to the multi-step nature. Which is how it is handled in cobrapy right now. If you want to specifically model what happens with the molecule inside the channel you have to add a species in the membrane compartment as well. I totally get that enzymes can have a single location in that case but reactions are not exactly the same thing.

MK:
@cdiener and @Midnighter: For scaling of reaction bounds based on RNA or protein data the incorrect compartment will give the incorrect scaling. I.e. proteins (reactions) in the membrane should be scaled based on Area, not volume. For this to work the reactions should be assigned the correct location, i.e. it happens at the membrane.
I.e.
See the following paper for discussion in the kinetic context: https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-020-03913-8

However, many processes happen at membranes and thus not in 3D, but in 2D. The
scaling of the rates of these processes poses a special problem, if volumes of compart-
ments are changed. They will typically scale with an area, but not with the volume of
the involved compartment. However, commonly, this is neglected when setting up
models and/or volume scaling also sometimes automatically happens when using
modelling software in the field.

This is a big issue when using Data to set bounds for reactions, or also to just account for dilution with cell growth (e.g. in approaches such as RBA). Because cobrapy is meant as a core library on which other packages can build more complex analysis the correct handling of this cases would be important. If cobrapy is just a simple tool for FBA, never mind.

Similar issues exist for visualization algorithms using the information in SBML. The important information about the localization of the reaction at the membrane is lost and the network is visualized incorrectly.
Yes, it makes sense to infer localization from metabolites if no localization information is available for the reaction. But if the user provides information this should not be discarded. I.e. it must be possible to set the compartment of a reaction and existing information should not be overwritten by heuristics.

CD:
@matthiaskoenig if the reaction happens in the membrane I think everybody agrees. But the argument from the paper seems incomplete as only one component is bound to the membrane but the other may still diffuse freely in the cytosol. If you scale the cell volume up this would still enter as a volumetric measure for that second freely diffusible component. And regarding the scope of cobrapy: as the name implies we provide infrastructure for constraint-based modeling which is admittedly not very well defined but I would say that kinetic models and (general) visualization would generally not fall within that scope. However, I do agree that information should not be dropped and we could save that either in the annotations or a separate attribute. Though it would be nice if SBML could support multiple compartments for an enzyme since membrane-bound components or channels do have parts embedded in several compartments. I am also not super eager of allowing that reactions may involve metabolites that have no overlap with the enzyme compartment now. Having a reaction assigned to the nucleus using metabolites from the extracellular space does not look right to me.

UA:
@cdiener @matthiaskoenig
I looked at the SBML standard, and it seems that only species have compartments. I think that means that per SBML, only metabolites have compartments, but I could be wrong. Compartments can have any number of dimensions between 0 and 3, so you could set up the membrane as a compartment with two dimensions (and then it would be scaled by area).
There also seems to be an ability to set up compartment type and compartment relationships using the multi plugin in SBML.
Some questions

Did I understand the SBML standard correctly?
What would you like to have compartments able to be assigned in cobrapy? Metabolites? Genes? Reactions? All of them? Would a compartment be unique to one type of object, or can the same one have genes and metabolites, for example.

MK:
@akaviaLab you are wrong. Please just look at the SBML specification. You definitely did not look in the document. Reactions have compartments.

So the main question is what to do if the compartment attribute is already set on a reaction or directly assigned. This should have precedence over any heuristics. That is all I am advocating here.

CD:
@matthiaskoenig we can't remove the reaction compartments logic in cobrapy right now due to backwards compatibility. Also the interpretation is different and more correct in cobrapy IMO as it actually ensures that reaction compartments are compatible with metabolite compartments (but we can ask the community to decide). But we can definitely store the SBML-assigned compartment and have that assigned as a different attribute like rxn.enzyme_compartment (which is what SBML seems to describe) or in the annotations. SBML writing would than use this compartment.

But either way it's unrelated to the PR and can be addressed in a separate discussion.

UA:
@matthiaskoenig @cdiener
Sorry, I got the standard wrong.
SBML requires that species have a compartment, but it seems optional for reactions (see the Reaction section 4.11.1 p 70 of the standard).
Could we set the code that if reaction compartment was set implicitly, it will be returned (and written to SBML). Otherwise, reaction.compartments keeps returning the current, which is a set of the compartments the metabolites have?

If both reactions and metabolites have compartments, I don't think it would make sense for genes to have compartments. Maybe in the future, if we want to do ME models (or similar stuff), there can be a protein that derives from gene and has compartments. Section 7.12 (p 137 of the standard). @cdiener and @Midnighter, do you agree genes shouldn't have compartments now?

[akaviaLab]

I put a very early draft of what I'd like to do at #1193
This is just to see if the API I'm trying makes sense to people before I flesh it out. model.py hasn't been changed yet.

[matthiaskoenig]

• Don't support compartment_type. This is only part of the old SBML Level 2 and not used for anything. This is why it was dropped in SBML level 3. I.e. drop the compartment_type.
• Reaction should have a compartment attribute! and behave as the compartment of the metabolite. I would probably split the members in metabolites and reactions on the compartment. I.e. it will create many issues if the members are metabolites + reactions.

[akaviaLab]

Wasn't this added back in the multi package of SBML 3 (section 3.5.2)? If it was added to multi, do you want to support it, or not now.

[matthiaskoenig]

I don't see any use for it in the context of constraint-based modeling, so I would not add it.

[matthiaskoenig]

Just to comment. I don't see how to a compartment on reaction will break anything ?! @cdiener could you add some information here. At the moment we only have a property

    @property
    def compartments(self) -> Set:
        """Return set of compartments the metabolites are in.

which will remain with the same behavior. The new compartment would allow to store the compartment for the reaction. As I understand the code base there is no issue at all to support reaction.compartment

[Midnighter]

How about Reaction.localization or .location? Or is that attribute name too far from its intended value, i.e., a compartment? I agree that compartment and compartments is too close and I would not want to rename to metabolite_compartments because that would break existing logic.

[matthiaskoenig]

location is perfect. I agree rxn.compartments is something different and will just stay what it is for backwards compatibility.

[cdiener]

Ok, perfect, I just misunderstood and we are all agreeing already 😅

[akaviaLab]

Okay. Updated my draft with the comments on the API.
In compartment, there is
direct_reactions, which are reactions that were assigned (have reaction.location be this compartment)
inferred_reactions, which are the reactions that have metabolites with this compartment
reactions, which includes both if you don't distinguish between them.

There is also metabolites

Since it inherits from group, you can access both of them with the members property, but you probably shouldn't.

[akaviaLab]

Replying to #1193 (comment)
Should reactions_from_metabolites include

1. Only reactions that have ALL metabolites in that compartment
2. Reactions that have ANY metabolites in that compartment
3. Something else, like reactions that have all products or reagents in the same compartment
4. Something completely different?

[matthiaskoenig]

Also GeneProducts could have compartments, but somehow indirectly.
I.e.

A GeneProduct may, optionally, refer to a Species so as to provide compatibility with the Manchester style encoding of gene-protein associations. In this case the attribute associatedSpecies is of type SIdRef and, if defined,
must point to an existing Species in the model.

I.e. a GeneProduct can reference a species (often used to encode the actual protein which is the gene product) which can have a localization. But this would be a bit indirect. I would suggest to not set a compartment on Gene or GeneProduct, but only have a compartment on metabolite and reaction.

[Midnighter]

I agree on this. Makes complete sense to me.

[akaviaLab]

Okay. This is for the various developers, including @Midnighter @matthiaskoenig @cdiener
What do you think Compartment should inherit from?
@cdiener suggested Group, but there seems to be some disagreement. I think it could be Object or maybe even Species, but I'm neutral on whatever makes sense.

Related question - does Optlang (should) have any compartment features? Right now, it doesn't, and I don't think it should, but I wanted to ask as part of this discussion.

[Midnighter]

Related question - does Optlang (should) have any compartment features? Right now, it doesn't, and I don't think it should, but I wanted to ask as part of this discussion.

Compartments are modeled implicitly in optlang due to metabolite constraint names having a compartment suffix so that there is no name collision either. Optlang is still a general purpose optimization package so I don't think a compartment is a concept that belongs there.

[Midnighter]

What do you think Compartment should inherit from?
@cdiener suggested Group, but there seems to be some disagreement. I think it could be Object or maybe even Species, but I'm neutral on whatever makes sense.

I think Object with an interface similar to Group makes the most sense to me.