JumpProcesses.jl icon indicating copy to clipboard operation
JumpProcesses.jl copied to clipboard

Tracking issue for Documenter and design of docstrings and new tutorials.

Open 00krishna opened this issue 3 years ago • 23 comments

I am creating a new issue to help track the work on migrating some of the DiffEqJump.jl documentation to Documenter.jl. In talking with @isaacsas , he indicated that he wanted to revamp the introductory tutorials for DiffEqJump.

Originally I thought I would just migrate the tutorials to Documenter and incorporate the docstrings into the tutorials. But instead of that, perhaps it is better to review how others have implemented the docstrings for some other larger packages in the combined SciML docs.

This issue can just help us to coordinate efforts and agree on a design.

00krishna avatar Jun 14 '22 18:06 00krishna

Thanks! Yes, my plan, which hopefully I can get to this weekend or next week was to split the big tutorial into several pieces and clean them up. Something along the lines of pure jump models with ConstantRate and MassActionJumps (perhaps split into two pieces itself), followed by a tutorial involving VariableRateJumps.

The FAQ would be split out of the tutorial to a separate section. Basically I'm envisioning something like the Catalyst doc structure: https://catalyst.sciml.ai

We could also add and cleanup the jump diffusion tutorial -- I'd leave that decision to @ChrisRackauckas on where he feels that should live. (It would also make sense in StochasticDiffEq.)

Then we need to figure out how to reconcile the current JumpProblem page, and to a lesser extent the jump solve page, with the API doc strings to make an API page. They don't completely contain the same information, so it may be the doc strings need to be cleaned up and updated so the API docs are equivalent.

isaacsas avatar Jun 14 '22 18:06 isaacsas

I don't think docstrings should be in the tutorials -- see the Catalyst style. I'd prefer the tutorials to introduce how to use DiffEqJump in the most common use cases, and then leave the docstrings for more detailed info about the specific components.

isaacsas avatar Jun 14 '22 18:06 isaacsas

@00krishna if you want to take a first pass the simplest thing to do would be to initially split the tutorial in two parts, and have the second start with a note that one should be familiar with the first before looking at it. We could also convert it to use @example blocks so the code runs dynamically. Then put the FAQ on its own page/section, and just migrate the JumpProblem and jump solving pages (perhaps under a "Reference and Guidance" section). Finally, there can be an API section with docstrings. Initially that should probably be something like

JumpProblem
SSAStepper

# Jump Types
ConstantRateJump
VariableRateJump
MassActionJump
RegularJump

# Aggregators (see src/aggregators/aggregators.jl for a list)
Direct
FRM
NRM
...

That would be a good starting point to then build out from I think.

isaacsas avatar Jun 14 '22 18:06 isaacsas

@isaacsas Yeah, I can run with this. I will look at the Catalyst.jl docs as a model then, and try to emulate. I will post my design here so that you can make sure I am on the right track @isaacsas .

00krishna avatar Jun 14 '22 20:06 00krishna

TODO still:

  • [x] convert everything to @example blocks
  • [ ] add strict mode
  • [x] split the tutorials
  • [x] get all the API references in
  • [ ] merge/reconcile the jump and solver type documentation with the current doc strings
  • [ ] add RDME tutorial
  • [x] FAQ entry on getting past values of solution

isaacsas avatar Jun 19 '22 20:06 isaacsas

@isaacsas say, I was looking through the tutorials as a basis for the docstrings. Seems like a lot of the tutorial code is focused on using Catalyst.jl and those examples--the main one being the SIR model. For the docstrings, should I put so much emphasis on Catalyst, or should I setup the problems as regular DifferentialEquations.jl problems, like ODEProblem or I believe DiscreteProblem would also work in some cases.

00krishna avatar Jun 24 '22 16:06 00krishna

I'm not sure what you mean? The tutorial starts with Catalyst, but everything after it is pure DiffEqJump. Also, most of the core functions have docstrings showing now at https://jump.sciml.ai/dev/api/

isaacsas avatar Jun 24 '22 16:06 isaacsas

Ahh, I see. I was looking at the tutorial on https://jump.sciml.ai/dev/tutorials/discrete_stochastic_example/ . Right, so the original model is Catalyst for the SIR example, but then subsequent examples are regular DiffEq. Yeah, I think the DiffEq examples are one line functions, like rate1(u,p,t) = p[1]*u[1]*u[2] , so I just read over them too quickly.

Thanks for setting me straight.

00krishna avatar Jun 24 '22 17:06 00krishna

But maybe I should update that tutorial to move the Catalyst stuff to the end. Then it can start with DiffEqJump but present Catalyst at the end for people that are interested only in reactions systems and want an easier interface.

isaacsas avatar Jun 24 '22 17:06 isaacsas

I think that the SIR model example is a good one. You could define the SIR model in both DiffEq, as well as in Catayst, and present them together? I think that threw me off was that the tutorials starts with a detailed description of the SIR model in Catalyst form. So once that expectation was set, I did not notice the DiffEq examples below. So using an extended example is good, but perhaps showing the SIR example in both Catalyst and DiffEq forms? Like can I define a standard Diffeq model f(du, u, p, t) and then use that in the DiscreteProblem just like you used the Catalyst model setup for sir_model. I will try it out, but that was my thinking on what I want to check.

00krishna avatar Jun 24 '22 17:06 00krishna

Yeah, I think it makes sense to reverse the order. i.e. show SIR in DiffEqJump first, and then at the end point out it is easier to do these types of chemical systems with Catalyst (which has some other advantages). That way the information is there, but people who just want to learn DiffEqJump and don't care about chemical systems won't get hung up on it. I'm working on such an update now. Thanks for the suggestion!

isaacsas avatar Jun 24 '22 17:06 isaacsas

Say, @isaacsas I had an idea to adjust the formatting of the code blocks in the docstrings, but I wanted to check with you first. I thought an adjustment might make the DiffEqJump code blocks a bit more consistent with the DifferentialEquations codeblocks.

The current code blocks are formatted like this:

using DiffEqJump
β = 0.1 / 1000.0; ν = .01;
p = (β,ν)
rate1(u,p,t) = p[1]*u[1]*u[2]  # β*S*I
function affect1!(integrator)
  integrator.u[1] -= 1         # S -> S - 1
  integrator.u[2] += 1         # I -> I + 1
end
jump = ConstantRateJump(rate1,affect1!)

rate2(u,p,t) = p[2]*u[2]      # ν*I
function affect2!(integrator)
  integrator.u[2] -= 1        # I -> I - 1
  integrator.u[3] += 1        # R -> R + 1
end
jump2 = ConstantRateJump(rate2,affect2!)
u₀    = [999,1,0]
tspan = (0.0,250.0)
prob = DiscreteProblem(u₀, tspan, p)
jump_prob = JumpProblem(prob, Direct(), jump, jump2)
sol = solve(jump_prob, SSAStepper())

But most of the original DifferentialEquations.jl codeblocks have the function first and the parameters afterwards. This adjustment also gives a bit more prominence to the model function itself, so that a reader does not miss it. I was thinking to reformat the block like this.

using DiffEqJump

function rate1(u,p,t) 
    p[1]*u[1]*u[2]  # β*S*I
end

function affect1!(integrator)
  integrator.u[1] -= 1         # S -> S - 1
  integrator.u[2] += 1         # I -> I + 1
end

jump = ConstantRateJump(rate1,affect1!)

function rate2(u,p,t)
    p[2]*u[2]      # ν*I
end

function affect2!(integrator)
  integrator.u[2] -= 1        # I -> I - 1
  integrator.u[3] += 1        # R -> R + 1
end

jump2 = ConstantRateJump(rate2,affect2!)

β = 0.1 / 1000.0; ν = .01;
p = (β,ν)
u₀    = [999,1,0]
tspan = (0.0,250.0)
prob = DiscreteProblem(u₀, tspan, p)
jump_prob = JumpProblem(prob, Direct(), jump, jump2)
sol = solve(jump_prob, SSAStepper())

Let me know either way, and I will use that as a model. Like I said, it is just an idea.

00krishna avatar Jun 25 '22 23:06 00krishna

I don't think putting the parameters after is really needed, and maybe it is a bit confusing (since the earlier declaration makes clear the parameters are beta and nu, which are then referenced in subsequent comments -- moving it after means the parameter names are referenced before they are actually defined).

An alternative would be to move the p declaration after, but make p a named tuple and then access it via the named fields.

isaacsas avatar Jun 27 '22 15:06 isaacsas

Okay. So then something like this below? It adds a bit of space and makes sure that the rate and affect functions themselves are visible on quick glance.


β = 0.1 / 1000.0; ν = .01;
p = (β,ν)

function rate1(u,p,t) 
    p[1]*u[1]*u[2]  # β*S*I
end

function affect1!(integrator)
  integrator.u[1] -= 1         # S -> S - 1
  integrator.u[2] += 1         # I -> I + 1
end

jump = ConstantRateJump(rate1,affect1!)

function rate2(u,p,t)
    p[2]*u[2]      # ν*I
end

function affect2!(integrator)
  integrator.u[2] -= 1        # I -> I - 1
  integrator.u[3] += 1        # R -> R + 1
end

jump2 = ConstantRateJump(rate2,affect2!)

u₀    = [999,1,0]
tspan = (0.0,250.0)
prob = DiscreteProblem(u₀, tspan, p)
jump_prob = JumpProblem(prob, Direct(), jump, jump2)
sol = solve(jump_prob, SSAStepper())

00krishna avatar Jun 27 '22 16:06 00krishna

Sure

isaacsas avatar Jun 27 '22 16:06 isaacsas

@isaacsas I created the first set of doctests for the jump types. One thing to note is that the doctests are building, but they are indicating a failure a warning. I created the docstring for RegularJumps and included the code below. The warning I am getting is about

The RegularJump interface has changed to be matrix-free. See the documentation for more details.
│    └ @ DiffEqJump /media/hayagriva/Dropbox/sandbox/DiffEqJump.jl/src/jumps.jl:181

I used the code below.

using DiffEqJump

function rate(out,u,p,t)
    out[1] = (0.1/1000.0)*u[1]*u[2]
    out[2] = 0.01u[2]
end

function c(dc,u,p,t,mark)
    dc[1,1] = -1
    dc[2,1] = 1
    dc[2,2] = -1
    dc[3,2] = 1
end

dc = zeros(3,2)
rj = RegularJump(rate,c,dc;constant_c=true)

So the problem seems to be that dc is an AbstractMatrix and that is dispatching to the function with the warning. But I think it is okay to publish this code, because it will just produce a warning, but the code will otherwise work. If you tell me how to adjust this code to avoid the warning, then I can incorporate that as well.

00krishna avatar Jun 28 '22 02:06 00krishna

Someone changed the interface a while back but I don't think it has gotten much use since the change. You can look at the regular jump tests

https://github.com/SciML/DiffEqJump.jl/blob/master/test/regular_jumps.jl

and the definitions to see how it works now (notice it takes a counts argument now too)

https://github.com/SciML/DiffEqJump.jl/blob/277cb15fcd3be6cf80b6a3f1fc9fcd69219af831/src/jumps.jl#L94-L118

It is also used in StochasticDiffEq at

https://github.com/SciML/StochasticDiffEq.jl/blob/master/test/tau_leaping.jl

isaacsas avatar Jun 28 '22 12:06 isaacsas

Getting a doc string on RegularJumps with an example would be good.

isaacsas avatar Jun 28 '22 12:06 isaacsas

@isaacsas Okay good. I got a working example--no warnings. Note that in order to include this example, I will need to include StableRNGs, and LinearAlgebra in the Project.toml. Technically I could probably get away without StableRNGs, but keeping it seems helpful if people want to make replicable examples.

I will go ahead and put it in the PR with this code.

using DiffEqJump
using StableRNGs
using LinearAlgebra
rng = StableRNG(12345)

function regular_rate(out, u, p, t)
  out[1] = (0.1 / 1000.0) * u[1] * u[2]
  out[2] = 0.01u[2]
end

function regular_c(du, u, p, t, counts, mark)
    mul!(du, dc, counts)
end

dc = zeros(3, 2)
rj = RegularJump(regular_rate, regular_c, 2)
jumps = JumpSet(rj)
prob = DiscreteProblem([999, 1, 0], (0.0, 250.0))
jump_prob = JumpProblem(prob, Direct(), rj; rng = rng)
sol = solve(jump_prob, SimpleTauLeaping(); dt = 1.0)

00krishna avatar Jun 28 '22 17:06 00krishna

@isaacsas i am traveling the yesterday and today, but will get back to the PR in the next day or so. Just wanted to let you know.

00krishna avatar Jun 30 '22 17:06 00krishna

No worries. Have a good trip!

isaacsas avatar Jun 30 '22 19:06 isaacsas

@isaacsas Gradually getting back on the horse :). So I did more research and thinking and it seems like doctests in the docstrings only support test cases where some result is returned and compared to a desired result. Those @example blocks only work in markdown documents, and are not run when pulling docstrings from source files.

So we could just go with julia blocks for now, and then add tests to the tutorials. There are a couple of ways to do this as well. We could either see if the existing tutorials can be setup with @example blocks that are sufficient tests for CI purposes. Or, we can also create a set of hidden documents where we can write more comprehensive tests that are run at build time, but that users never see. This hidden option works if say the tutorials have code blocks that are not full tests or code runs, etc.

I will read over the new tutorials to figure out the options.

00krishna avatar Jul 11 '22 17:07 00krishna

I think julia blocks for docstrings with simple examples aimed at helping users see how to use the functions are fine. I don't really see the point of silent tests via the documentation -- at that point we might as well add the associated test as a normal package test instead. We have testing pretty well covered with the normal tests, and I've setup the new documentation to run all the tutorials as example blocks.

isaacsas avatar Jul 11 '22 19:07 isaacsas