Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enable functional/time dependent parameters #1191

Merged
merged 22 commits into from
Apr 3, 2025
Merged

Conversation

TorkelE
Copy link
Member

@TorkelE TorkelE commented Mar 7, 2025

Adds docs and tests. Also modifies the main code to enable these to be used within Catalyst.

Example:

using Catalyst, DataInterpolations, OrdinaryDiffEqDefault, Plots

# Create input function/data
tend = 120.0
ts = collect(0.0:1.0:tend)
light = sin.(ts/6) .+ 1
light = [max(0.0, l - rand()) for l in light]
interpolated_light = LinearInterpolation(light, ts)
plot(interpolated_light)

image

# Create model
t = default_t()
in_type = typeof(interpolated_light)
@parameters kA kD (light_in::in_type)(..)
@species Pₐ(t) Pᵢ(t)
rxs = [
    Reaction(kA*light_in(t), [Pᵢ], [Pₐ]),
    Reaction(kD, [Pₐ], [Pᵢ])
]
@named rs = ReactionSystem(rxs, t)
rs = complete(rs)

# Simulate the model.
u0 = [Pᵢ => 1.0, Pₐ => 0.0]
ps = [kA => 1.5, kD => 1.0, light_in => interpolated_light]
oprob = ODEProblem(rs, u0, tend, ps)
sol = solve(oprob)
plot(sol)

image

@TorkelE TorkelE changed the title [WIP] Enable functional/time depe3ndent parameters [WIP] Enable functional/time dependent parameters Mar 9, 2025
@TorkelE TorkelE changed the title [WIP] Enable functional/time dependent parameters Enable functional/time dependent parameters Mar 9, 2025
@TorkelE
Copy link
Member Author

TorkelE commented Mar 9, 2025

Builds on the change in #1192, however, but otherwise this is now ready.

@species I(default_t())
inf_rate = I_rate(I)
sir = @reaction_network rs begin
k1*$inf_rate, S --> I
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It would be nicer if we could make something like

Suggested change
k1*$inf_rate, S --> I
k1*($I_rate)(I), S --> I

work.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Having to declare t and/or I(t) symbolically is a bit clunky.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed. had wanted to do this but get a

ERROR: LoadError: TypeError: in isdefined, expected Symbol, got a value of type Expr
Stacktrace:
  [1] recursive_escape_functions!(expr::Expr, syms_skip::Vector{Any})
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:980
  [2] (::Catalyst.var"#232#233"{Expr, Vector{Any}})(i::Int64)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:978
  [3] foreach
    @ ./abstractarray.jl:3187 [inlined]
  [4] recursive_escape_functions!(expr::Expr, syms_skip::Vector{Any})
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:978
  [5] recursive_escape_functions!
    @ ~/.julia/dev/Catalyst/src/dsl.jl:977 [inlined]
  [6] get_rxexpr(rx::Catalyst.DSLReaction)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:555
  [7] #219
    @ ~/.julia/dev/Catalyst/src/dsl.jl:546 [inlined]
  [8] foreach(f::Catalyst.var"#219#221"{Expr}, itr::Vector{Catalyst.DSLReaction})
    @ Base ./abstractarray.jl:3187
  [9] get_rxexprs(reactions::Vector{Catalyst.DSLReaction}, equations::Vector{Expr}, all_syms::Vector{Union{Expr, Symbol}})
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:546
 [10] make_reaction_system(ex::Expr, name::QuoteNode)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:323
 [11] make_rs_expr(name::QuoteNode, network_expr::Expr; complete::Bool)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:153
 [12] make_rs_expr(name::QuoteNode, network_expr::Expr)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:152
 [13] var"@reaction_network"(__source__::LineNumberNode, __module__::Module, name::Symbol, network_expr::Expr)
    @ Catalyst ~/.julia/dev/Catalyst/src/dsl.jl:92
in expression starting at /home/loman/Desktop/Julia Playground/Environment - Catalyst Docs (full)/catalyst_docs_playground.jl:71

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let’s get an issue to follow up and fix this then. (Post V15.)

Copy link
Member

@isaacsas isaacsas left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I still think the initial example is likely to just confuse readers as to what a functional parameter is, and would be better off just using a basic function instead of bringing in DataInterpolations (or else add an explanation that you want to show how one can use functional parameters to represent time-dependent parameters that are only measured at a discrete set of points via building an interpolant and setting it to be the functional parameter). But that said, feel free to merge if you want.

@TorkelE
Copy link
Member Author

TorkelE commented Apr 3, 2025

I added a note at the end of that example to clarify things a bit, but I have now expanded further to

!!! note
    For this simple example, $(2 + t)/(1 + t)$ could have been used directly as a reaction rate (or written as a normal function), technically making the functional parameter approach unnecessary. However, here we used this function as a simple example of how discrete data can be made continuous using DataInterpolations, and then have its values inserted using a (functional) parameter.

Do you think that is helpful?

@isaacsas
Copy link
Member

isaacsas commented Apr 3, 2025

Sure, seems helpful. In any case, finish it as you see fit and then merge. (I should move on to other V15 blocker PRs.)

@TorkelE
Copy link
Member Author

TorkelE commented Apr 3, 2025

Sounds good, I will merge like this (once things pass). Happy to get back to this though. It might even make sense to have a dedicated example tutorial to the various ways of modelling functions.

@TorkelE TorkelE merged commit d71d11a into master Apr 3, 2025
13 checks passed
@TorkelE TorkelE deleted the functional_parameters branch April 3, 2025 20:46
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants