Add SBML support for non-mass action rates

[WilliamIX]

Minimally support all the rate-types in bioscrape. Any thoughts on how to best deal with general rate expressions? I am not super happy with using sympy to parse expressions as is done in bioscrape.

[ayush9pandey]

SBML can parse almost all rate expressions using parseL3Formula, for example,

parseL3Formula('k1f * X^2 / (X^2 + Y^2)')

will look for species with id X and Y, and create the rate law accordingly. For more complicated rate expressions, there are FunctionDefinition, RateRule, etc. available in SBML.

Where in the code currently do we feel that there is a need for more complicated rate expressions?

[WilliamIX]

Currently the reaction class in chemical_reaction_network.py has two attributes:

1: mass_action (bool) I think this can be expanded/improved so that instead of using mass_action as a boolean we have a string "reaction type". Then we can have some automatic support for certain kinds of reactions (eg hill functions) and also the option to use custom rate equations. I would ideally like automatic support (without parsing) for all the kinds of propensities supported by bioscrape.

2: rate_formula The rate formula is supposed to take formulas of the kind you wrote down which can then be converted into SBML. That is not super straight forward though because you need to know which species are involved in the rate formula so you can add them as modifier species to the kinetic_rate_law sbml object. Also, rate formulas will have to be generated when saving sbml objects from the innately supported propensity types.

As a note, the reason I want to match bioSCRAPES propensities innately is because the cython implementations are much faster to simulate than the parsed "custom" rate laws in the current bioscrape implementation. So, if BioCRNpyler simulates directly through bioscrape without converting to SBML it should have considerably better performance than it would get if all non-mass action rate laws have to be parsed.

[ayush9pandey]

Yeah, having an automatic support for certain kinds of reactions and propensities sounds like a good feature to have.

For the rate_formula of the kind that I wrote down, SBML can automatically detect the species given in the string and add them as modifier species to the KineticLaw object under the reaction to which that string will be passed.

Also, regarding fast simulations, a highly developed and sophisticated simulator called libRoadRunner is available for SBML models. libRoadRunner has been in development for many years and their focus has been optimal and fast simulation of models. They have algorithms for optimal generation of C models (much like bioscrape Cython) and similar other features for fast simulations. In my opinion, there is a two fold advantage in going on this route -

1. You'll always generate SBML models, which is the standard, and SBML supports every kind of rate law (from the most simplest to the most complicated)
2. You'll be able to use libRoadRunner which is a highly developed and fast simulation tool. (and also currently supported/improved by the developers).

[murrayrm]

Side comment: If libRoadRunner is faster than the cython simulator in bioscrape, might be worth replacing the bioscrape simulator with libRoadRunner and leaving the system identification functions in place (which should run faster).

[WilliamIX]

I just looked into libRoadRunner - without doing any speed tests I can already see that there are quite a few differences in terms of functionality. Roughly speaking:

LibRoadRunner has many more features for doing analytics / theoretical computations about models (eigenvalues, frequency domain simulation, etc). LibRoadRunner is also written in C++ so it should be at least as fast as bioSCRAPE. However, because libRoadRunner is written in C++, it will be considerably harder to modify/add features. Additionally, accessing model parameters / attributes via libRoadRunner is only accomplished via an api that directly connects to SBML attributes - therefore in order to understand how to programatically access libRoadRunner models (say from python), some knowledge of SBML is required.

BioSCRAPE has more kinds of simulators - in particular it supports stochastic delay differential equations and cell lineage tracking for arbitrary stochastic simulations. Additionally BioSCRAPE supports algebraic assignment rules, which are a very flexible way to implement anything from conservation laws to steady-state approximations. Finally, the structure of BioSCRAPE makes it easy to add new simulators of different kinds. Once my bioSCRAPE model API is complete, it should also be easier to access bioSCRAPE objects in a programatic way without any knowledge of SBML.

It would be very easy to add bioCRNpyler functionality to automatically simulate with libRoadRunner by first writing and SBML file and then loading it. If libRoadRunner is better for parameter inference (no idea if this is true), this would be a straightforward pipeline to implement. On the other hand, bioSCRAPE will run more quickly if bioCRNpyler talks to it directly compared to if bioSCRAPE first loads sbml. Additionally, it should be possible to add libRoadRunner as a simulator inside bioSCRAPE - for ODE and SSA simulations, this could potentially result in somewhat improved performance.

[WilliamIX]

This works and has been tested for everything but general propensity types

[WilliamIX]

@ayush9pandey can you test this, including with your newest bioscrape changes and as soon as it is officially ready close the issue? I think we are 95% there.

[ayush9pandey]

I checked this code and it works fine. I would say we are 99% there now :), let's just wait until sbml_fix gets merged to master before closing this issue.

[zoltuz]

PR #83 has been merged, have a look and close this issue if it works now.

[WilliamIX]

I think this works now.