sys-bio / tellurium

Python Environment for Modeling and Simulating Biological Systems
http://tellurium.analogmachine.org/
Apache License 2.0
110 stars 36 forks source link

Simulating stoichiometry with unknown coefficients #555

Closed janursa closed 1 year ago

janursa commented 2 years ago

In current version of tellerium, we can do,

model = """
        x = 1; 
        y = 10;
        k = .1;
        x_y = 0;
        x + 10  y -> x_y; k*x*y^10
    end
"""

what i want to do is to define 10 as a parameter (because it can be unknown) as,

model = """
        x = 1;
        y = 10;
        k = .1;
        x_y = 0;
        n=10;
        x + n  y -> x_y; k*x*y^n
    end
"""

This gives an error. I tried using n*y, but same. Is there a way to implement this? thanks

matthiaskoenig commented 2 years ago

Hi @janursa, variable stoichiometries are to my knowledge not supported (and you would run in many issues trying to integrate/simulate such a model not only with tellurium/roadrunner). I would try to create the various model versions with the different stoichiometries, i.e. and then compare the different models against each other (e.g. which model version fits best the data).

I.e. I would create multiple models along the line of

models = []
for n in range(11):

  model = f"""
          x = 1;
          y = 10;
          k = .1;
          x_y = 0;
          x + {n}  y -> x_y; k*x*y^{n}
      end
  """
  models.append(model)

Code not tested, but you get the general idea.

The second solution could be to have a model with all the reactions, i.e.

x + 1 y -> x_y; n_switch_1*k*x*y^1
...
x + 10  y -> x_y; n_switch_10*k*x*y^10

and an additional function which sets only one of the switch_n_* parameters to 1 and all other to 0, i.e. a piecwiese which is doing

if n==1:
   n_switch_1=1
   n_switch_2=0
   ...
   n_switch_10=0

if n==2:
   n_switch_1=0
   n_switch_2=1
   ...
   n_switch_10=0
...

The second version is most likely what you want because you can switch the stoichiometry based on the parameter n. If you need help with the piecwice let me know.

Best M

luciansmith commented 2 years ago

The limitation here is Antimony: antimony's syntax doesn't allow named stoichiometries. However, roadrunner itself will still simulate variable stoichiometries if you give it an SBML file with them, so if you want this, you'll need to create your SBML model with it, then load that in directly. (Support for this is limited, but should work in most situations.)

One option if you don't want to use a different model editor is to create the model in Antimony with a normal '10' as the stoichiometry, and then go into the SBML by hand and give that speciesReference an ID, i.e. change:

          <speciesReference species="S2" stoichiometry="10" constant="true"/>

to;

          <speciesReference id="n" species="S2" stoichiometry="10" constant="true"/>

Then you should be able to set the value of n inside tellurium with commands like 'r.n = 3' or what have you.

If you want 'n' to change over the course of the simulation itself, I would do essentially the same thing, but put in 'n' as a separate parameter. Then when you edit the SBML, change the ID of the speciesReference to 'n', and just delete the 'n' parameter.

If you would like this feature in Antimony directly (i.e. if you feel you're going to be doing this sort of thing more than once), please file an issue with that request at https://github.com/sys-bio/antimony/issues and I'll see what I can do!

janursa commented 2 years ago

@matthiaskoenig thanks for the suggestions. they might be practical for smaller n, but would be cumbersome if the n varies in a larger range

janursa commented 2 years ago

@luciansmith the proposed approach works but since i am in the development phase, where i would need to repeat this process over an over, editing the xml file each time would be difficult. i will create the issue and would be great if that gets integrated into the antimony directly

luciansmith commented 1 year ago

This is now implemented in Antimony, so we can close this. Do let me know if there are any other issues with this that need fixing.