Subsampling an ODESolution at Integer Times (and then plotting)

New to Julia
1

Hi,

I’m currently playing around with DifferentialEquations.jl and Catalyst.jl for simulating the paths of Markov Jump Processes, with an end goal of eventually using it to prototype some Bayesian inference algorithms.

Having followed the documentation here, I’m able to simulate paths of my MJP of interest, and plot their solution with Plots.jl, as below.

using Catalyst, DifferentialEquations

lotvol_model = @reaction_network begin
    c1, x1 --> 2x1
    c2, x1 + x2 --> 2x2
    c3, x2 --> 0
end c1 c2 c3

rates = (1.0, 0.005, 0.6)
prob = DiscreteProblem(lotvol_model, [70, 80], (0.0, 50.0), rates)
jump_prob = JumpProblem(lotvol_model, prob, Direct())

sol = solve(jump_prob, SSAStepper())
using Plots; plot(sol)

and everything is looking great so far.

What I’d like to do next is take my solution, subsample it at integer times, and be able to plot that in the same format (as this is the sort of data upon which I will eventually be running inference algorithms).

Thus far, I’ve tried both

sol_int_times = [ sol(t) for t in 1:50]

and

sol_int_times = [(t, sol(t)) for t in 1:50]

With the former attempt, I can plot it, but it forgets the temporal structure of the solution, and so is not really what I’m looking for. The issue with the latter attempt is that it won’t plot at all. Upon inspecting sol_int_times, the time and position components seem to be nested unevenly, or something like this - I don’t completely understand, but clearly I’m transforming things wrongly somehow.

Anyways, I would like to be able to subsample the solution to the MJP at the specified times, and have it stored in a format which will plot cleanly (and hopefully, be easy to work with downstream). What’s the correct syntax for doing so?

1 Like
2

sol = solve(jump_prob, SSAStepper(), saveat=1.0)

2 Likes
3

Great, thanks. Am I correct in understanding that this re-runs the solver from scratch? If so, is there a way of extracting the original path without re-solving?

4

If you solve like this it solves once and saves at the integers. Isn’t that what you want?

5

I think it ought to be fine for the most part — I was thinking that, further down the track, it might be useful to be able to subsample the solution more generally, if manipulating the solutions in this way is not too messy. It’s not super important though, I think that I should be able to do what I want with this.

6

I mean, from the way you solved it you could do sol(1:50).

1 Like
7

Perfect - thank you!