Proper syntax for adding reactions to a reaction network using addreaction!() function in Catalyst.jl

Specific Domains Modelling & Simulations
, ,
1

I recently started to thinking about switching my simulations to Julia and got to know about Catalyst.jl package and its convenient macro to create and modify reaction networks.

I wanted to add a reaction to an already defined reaction network using the addreaction!() function but am unable to do so because I am unable to understand the syntax mentioned.
This is my network:

@reaction_func Hsn(thr,N,lbd,cop) = ((thr^cop)/((thr^cop) + (N^cop))) - lbd*(1-((thr^cop)/((thr^cop) + (N^cop))))
@reaction_func Hsp(thr,N,lbd,cop) = (((thr^cop)/((thr^cop) + (N^cop))) - lbd*(1-((thr^cop)/((thr^cop) + (N^cop)))))/lbd
coupled_switch = @reaction_network begin
    (kA,kB,kC,kD), (A,B,C,D) → ∅
    gA*(Hsn(tBA,B,fBA,nBA)*Hsn(tCA,C,fCA,nCA)), ∅ → A
    gB*(Hsn(tAB,A,fAB,nAB)*Hsn(tDB,B,fDB,nDB)), ∅ → B
    gC*(Hsn(tAC,A,fAC,nAC)*Hsn(tDC,D,fDC,nDC)), ∅ → C
    # This is the reaction I am trying to add as a test.
    #gD*(Hsn(tBD,B,fBD,nBD)*Hsn(tCD,C,fCD,nCD)), ∅ → D 
end

And to add the reaction I am doing:

addreaction!(coupled_switch, gD*(Hsn(tBD,B,fBD,nBD)*Hsn(tCD,C,fCD,nCD)), ∅ → D )

and I am getting the following error:

UndefVarError: Hsn not defined
in top-level scope at r.jl:45

I realise this is completely wrong, but I am unable to find the right method to do so as I am unable to understand the documentation.
The documentation does not mention anything about the format of the the reaction to be passed to the addreaction!() function. Does anyone know the right way to do it?

1 Like
2

Yes, this should be doable. Sorry we haven’t yet gotten a nice tutorial up on how to use the network extension methods – we do need better documentation here.

In this case it is probably easiest, and most flexible, to use the underlying ModelingToolkit representations to do this. The @reaction_network macro produces a ModelingToolkit.ReactionSystem, which stores a collection of ModelingToolkit.Reactions. addreaction! modifies this system by adding in a Reaction that you just create yourself. So for your example, try something like:

# these commands create symbolic variables in ModelingToolkit for the rates and species
@parameters t kA kB kC kD gA tBA fBA nBA tCA fCA nCA
@variables A(t) B(t) C(t)

# this creates a normal Julia function that will work with the symbolic variables 
Hsn(thr,N,lbd,cop) = ((thr^cop)/((thr^cop) + (N^cop))) - lbd*(1-((thr^cop)/((thr^cop) + (N^cop))))

# this creates a new Reaction and adds it to coupled_switch, see
# https://catalyst.sciml.ai/dev/api/catalyst_api/#ModelingToolkit.Reaction
# for more examples of creating Reactions
addreaction!(coupled_switch, Reaction(gA*(Hsn(tBA,B,fBA,nBA)*Hsn(tCA,C,fCA,nCA)), nothing, [A]))
3

Just to add, we also have a @add_reactions macro that uses the DSL syntax. It doesn’t work on your example, as it requires any species that appear in a rate to also be used as a substrate or product (whereas in your example species B only appears in the rate). A simpler example where it would work is:

@add_reactions coupled_switch begin
       gD*D, ∅ → D 
end gD

This lets you use the DSL syntax to extend coupled_switch. I opened an issue, so hopefully we can get it working more generally on your example.

4

I guess I just need to do some minor adjustments to my code to get it working. That’s a very convenient way to go about it.
Thank You very much for clearing my doubt!!!