Merge dcefeec into 3fd02b0

docs/Project.toml

@@ -1,19 +1,29 @@
 [deps]
 Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
+DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d"
 DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
+Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+GR_jll = "d2c73de3-f751-5644-a686-071e5b155ba9"
 Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
 JumpProcesses = "ccbc3e58-028d-4f4c-8cd5-9ae44345cda5"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
+PointProcesses = "af0b7596-9bb0-472a-a012-63904f2b4c55"
+StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
 
 [compat]
 Catalyst = "13.0"
-DifferentialEquations = "7.10"
-Documenter = "1"
-Graphs = "1.8.0"
+DensityInterface = "0.4"
+DifferentialEquations = "7.11"
+Distributions = "0.25"
+Documenter = "1.1"
+GR_jll = "0.72"
+Graphs = "1.9"
 JumpProcesses = "9.8"
-OrdinaryDiffEq = "6.56"
+OrdinaryDiffEq = "6.59"
 Plots = "1.39"
-StochasticDiffEq = "6.62"
+PointProcesses = "0.4"
+StableRNGs = "1.0"
+StochasticDiffEq = "6.63"

docs/make.jl

@@ -8,24 +8,24 @@ cp(joinpath(docpath, "Project.toml"), joinpath(assetpath, "Project.toml"), force
 include("pages.jl")
 
 mathengine = MathJax3(Dict(:loader => Dict("load" => ["[tex]/require", "[tex]/mathtools"]),
-                           :tex => Dict("inlineMath" => [["\$", "\$"], ["\\(", "\\)"]],
-                                        "packages" => [
-                                            "base",
-                                            "ams",
-                                            "autoload",
-                                            "mathtools",
-                                            "require",
-                                        ])))
+    :tex => Dict("inlineMath" => [["\$", "\$"], ["\\(", "\\)"]],
+        "packages" => [
+            "base",
+            "ams",
+            "autoload",
+            "mathtools",
+            "require",
+        ])))
 
 makedocs(sitename = "JumpProcesses.jl",
-         authors = "Chris Rackauckas",
-         modules = [JumpProcesses],
-         clean = true, doctest = false, linkcheck = true, warnonly = [:missing_docs],
-         format = Documenter.HTML(; assets = ["assets/favicon.ico"],
-                                  canonical = "https://docs.sciml.ai/JumpProcesses/",
-                                  prettyurls = (get(ENV, "CI", nothing) == "true"),
-                                  mathengine),
-         pages = pages)
+    authors = "Chris Rackauckas",
+    modules = [JumpProcesses],
+    clean = true, doctest = false, linkcheck = true, warnonly = [:missing_docs],
+    format = Documenter.HTML(; assets = ["assets/favicon.ico"],
+        canonical = "https://docs.sciml.ai/JumpProcesses/",
+        prettyurls = (get(ENV, "CI", nothing) == "true"),
+        mathengine),
+    pages = pages)
 
 deploydocs(repo = "github.com/SciML/JumpProcesses.jl.git";
-           push_preview = true)
+    push_preview = true)

docs/pages.jl

@@ -1,12 +1,14 @@
 # Put in a separate page so it can be used by SciMLDocs.jl
 
 pages = ["index.md",
+    "getting_started.md",
     "Tutorials" => Any["tutorials/simple_poisson_process.md",
-                       "tutorials/discrete_stochastic_example.md",
-                       "tutorials/jump_diffusion.md",
-                       "tutorials/spatial.md"],
+        "tutorials/discrete_stochastic_example.md",
+        "tutorials/jump_diffusion.md",
+        "tutorials/point_process.md",
+        "tutorials/spatial.md"],
+    "Type Documentation" => Any["Jumps, JumpProblem, and Aggregators" => "jump_types.md",
+        "Jump solvers" => "jump_solve.md"],
     "FAQ" => "faq.md",
-    "Type Documentation" => Any["Jump types and JumpProblem" => "jump_types.md",
-                                "Jump solvers" => "jump_solve.md"],
     "API" => "api.md",
 ]

docs/src/api.md

@@ -12,7 +12,7 @@ SSAStepper
 reset_aggregated_jumps!
 ```
 
-## Types of Jumps
+## Jump Types
 
 ```@docs
 ConstantRateJump
@@ -22,7 +22,7 @@ RegularJump
 JumpSet
 ```
 
-## Aggregators
+## Aggregator Types
 
 Aggregators are the underlying algorithms used for sampling
 [`ConstantRateJump`](@ref)s, [`MassActionJump`](@ref)s, and

docs/src/assets/Project.toml

@@ -1,19 +1,29 @@
 [deps]
 Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
+DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d"
 DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
+Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+GR_jll = "d2c73de3-f751-5644-a686-071e5b155ba9"
 Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
 JumpProcesses = "ccbc3e58-028d-4f4c-8cd5-9ae44345cda5"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
+PointProcesses = "af0b7596-9bb0-472a-a012-63904f2b4c55"
+StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
 
 [compat]
 Catalyst = "13.0"
-DifferentialEquations = "7.10"
-Documenter = "1"
-Graphs = "1.8.0"
+DensityInterface = "0.4"
+DifferentialEquations = "7.11"
+Distributions = "0.25"
+Documenter = "1.1"
+GR_jll = "0.72"
+Graphs = "1.9"
 JumpProcesses = "9.8"
-OrdinaryDiffEq = "6.56"
+OrdinaryDiffEq = "6.59"
 Plots = "1.39"
-StochasticDiffEq = "6.62"
+PointProcesses = "0.4"
+StableRNGs = "1.0"
+StochasticDiffEq = "6.63"

docs/src/faq.md

@@ -50,7 +50,7 @@ vj2 = VariableRateJump(rate4, affect4!)
 vjtuple = (vj1, vj2)
 
 jset = JumpSet(; constant_jumps = cjvec, variable_jumps = vjtuple,
-               massaction_jumps = mass_act_jump)
+    massaction_jumps = mass_act_jump)
 ```
 
 ## How can I set the random number generator used in the jump process sampling algorithms (SSAs)?
@@ -62,7 +62,7 @@ argument. Continuing the previous example:
 #] add RandomNumbers
 using RandomNumbers
 jprob = JumpProblem(dprob, Direct(), maj,
-                    rng = Xorshifts.Xoroshiro128Star(rand(UInt64)))
+    rng = Xorshifts.Xoroshiro128Star(rand(UInt64)))
 ```
 
 uses the `Xoroshiro128Star` generator from
@@ -148,7 +148,7 @@ function paffect!(integrator)
     nothing
 end
 sol = solve(jprob, SSAStepper(), tstops = [20.0],
-            callback = DiscreteCallback(pcondit, paffect!))
+    callback = DiscreteCallback(pcondit, paffect!))
 ```
 
 Here at time `20.0` we turn off production of `u[2]` while activating production

docs/src/getting_started.md

@@ -0,0 +1,189 @@
+# Getting Started with JumpProcesses in Julia
+
+Jumps (or point) processes are stochastic processes with discrete state changes
+driven by a `rate` function. Historically, jump processes have been developed in
+the context of dynamical systems to describe dynamics with sudden changes —
+the jumps — in a system's value at random times. In contrast, the development
+of point processes has been more focused on describing the occurrence of random
+events — the points — over a support. In reality, jump and point processes
+share many things in common which make JumpProcesses ideal for both.
+
+Processes involving multiple jumps are known as compound jump (or point)
+processes.
+
+```math
+du = \sum_{j} h_j(u,p,t) dN_j(t)
+```
+
+where ``N_j`` is a jump process with rate ``\lambda_i(u,p,t)``.
+
+The homogeneous Poisson process is the canonical point process with a constant
+rate of change. A compound Poisson process is a continuous-time Markov Chain
+where the time to the next jump is exponentially distributed as determined by
+the aggregate rate. In the statistics literature, the composition of Poisson
+processes is described by the superposition theorem. Simulation algorithms for
+these types of processes are known in biology and chemistry as Gillespie methods
+or Stochastic Simulation Algorithms (SSA), with the time evolution that the
+probability these processes are in a given state at a given time satisfying the
+Chemical Master Equation (CME).
+
+Any differential equation can be extended by jumps. For example, we have an ODE
+with jumps, denoted by
+
+```math
+du = f(u,p,t)dt + \sum_{j} h_j(u,p,t) dN_j(t)
+```
+
+Extending a stochastic differential equation (SDE) to have jumps is commonly known as a jump-
+diffusion, and is denoted by
+
+```math
+du = f(u,p,t)dt + \sum_{i}g_i(u,t)dW_i(t) + \sum_{j}h_i(u,p,t)dN_i(t)
+```
+
+The general workflow with any of the jump processes above is to define the base
+and jump problem, solve the jump problem and then analyze the solution. The full
+code for a jump process with no other dynamics apart from jumps is:
+
+```@example ex0
+using JumpProcesses
+u0 = [0]
+tspan = (0.0, 10.0)
+dprob = DiscreteProblem(u0, tspan)
+rate(u, p, t) = 2.0
+affect!(integrator) = (integrator.u[1] += 1)
+jump = ConstantRateJump(rate, affect!)
+jprob = JumpProblem(dprob, Direct(), jump)
+sol = solve(jprob, SSAStepper())
+
+using Plots
+plot(sol, title = "Sample path from a jump process with constant rate",
+    label = "N(t)", xlabel = "t", legend = :bottomright)
+```
+
+## Step 1: Defining a problem
+
+The first thing you want to do is to define your base problem from the many
+options available. For dynamics that involve only jumps we employ
+a `DiscreteProblem` as our base problem.
+
+```@example ex0
+using JumpProcesses
+u0 = [0]
+tspan = (0.0, 10.0)
+dprob = DiscreteProblem(u0, tspan)
+```
+
+For our example, notice that `u0` is a `Int[]` which is an appropriate choice in
+this case because we will only be working with jumps. Since we will be modifying
+`u0` only when a jump occurs via its `affect!`, we initialize `DiscreteProblem`
+without any function mapping which means `DiscreteProblem` use the default
+identity mapping.
+
+JumpProcesses exports
+[`DiscreteProblem`](https://docs.sciml.ai/DiffEqDocs/stable/types/discrete_types/)
+from OrdinaryDiffEq. In case you want
+to model the base dynamics as an
+[`ODEProblem`](https://docs.sciml.ai/DiffEqDocs/stable/types/ode_types/) or as
+an [`SDEProblem`](https://docs.sciml.ai/DiffEqDocs/stable/types/sde_types/)
+you will need to import
+OrdinaryDiffEq.
+
+## Step 2: Defining a jump
+
+Jumps are implemented as [callbacks
+functions](https://docs.sciml.ai/DiffEqDocs/stable/features/callback_functions/)
+from DifferentialEquations. Our library takes care of adding the appropriate
+callbacks to the base problem to ensure that jumps are properly initialized and
+executed.
+
+Here we add a `ConstantRateJump` to the base problem. As such the user must
+define the `rate` and `affect!` which determines the frequency and the effect of
+the jumps on the base problem.
+
+```@example ex0
+rate(u, p, t) = 2.0
+affect!(integrator) = (integrator.u[1] += 1)
+jump = ConstantRateJump(rate, affect!)
+```
+
+There are many types of jumps. Some of them allow for variable rates. You can
+read about [different jump types in this page](@ref jump_types).
+
+Once the jump is defined, we initialize the `JumpProblem`. _Aggegators_ are
+algorithms that determines jump times. We call them _aggregators_ because they
+aggregate all jump callbacks into a single callback. Alternatively, we can think
+of aggregators as the jump simulation algorithms. In this case we use the
+`Direct` _aggregator_.
+
+```@example ex0
+jprob = JumpProblem(dprob, Direct(), jump)
+```
+
+JumpProblem can be initialized with different combination of jumps and
+aggregators. If you want to understand more about the options available read
+about [JumpProblem initialization in here](@ref defining_jump_problem).
+
+In addition to that, aggregators come with different trade-offs. Not all
+aggregators accept all jump types. To learn [more about _aggregators_ check this
+section](@ref Jump-Aggregators-for-Exact-Simulation).
+
+## Step 3: Solving a problem
+
+After defining a problem, we solve it using `solve` with an appropriate stepper.
+
+```@example ex0
+sol = solve(jprob, SSAStepper())
+```
+
+Steppers tell how we evolve time in our base problem and required for any
+numerical simulator. JumpProcesses offers the `SSAStepper` which steps through
+time one jump candidate at a time.
+
+If you are modelling other types of base problems like `ODEProblem` or
+`SDEProblem`, you will not be able to use `SSAStepper` since these problems
+require more fine-grained time evolution.
+
+Apart from time-stepping, you might also be interested in controlling the saving
+frequency of the state variable `u`. This control can be thought as orthogonal
+to how the stepper evolves time. To avoid saving at every jump, we can
+initialize jumps as following.
+
+```@example ex0
+jprob = JumpProblem(dprob, Direct(), jump; save_positions = (false, false))
+```
+
+Finally, to solve the problem at regular intervals we can use `saveat`.
+
+```@example ex0
+sol = solve(jprob, SSAStepper(); saveat = 1.0)
+```
+
+When you do not save the jump events, be careful when analysing interpolated
+values as they will not be an accurate representation of the sampled path. This
+can be particularly problematic plotting the data.
+
+## Going Beyond the Poisson process: How to Use the Documentation
+
+This tutorial covered only the basics of the Poisson process with constant rate.
+We mostly focused on the basic pattern to define and solve a JumpProblem.
+
+In case you want to go further, JumpProcesses is a component package in the
+[SciML](https://sciml.ai/) ecosystem, and one of the core solver libraries
+included in
+[DifferentialEquations.jl](https://docs.sciml.ai/DiffEqDocs/stable/).
+
+The documentation includes tutorials and examples:
+
+  - [Simulating basic Poisson processes](@ref poisson_proc_tutorial)
+  - [Simulating jump processes via SSAs (i.e., Gillespie methods)](@ref ssa_tutorial)
+  - [Simulating jump-diffusion processes](@ref jump_diffusion_tutorial)
+  - [Temporal point processes (TPP)](@ref tpp_tutorial)
+  - [Spatial SSAs](@ref Spatial-SSAs-with-JumpProcesses.jl)
+
+In addition to that the document contains references to guide you through:
+
+  - [References on the types of jumps and available simulation methods](@ref jump_problem_type)
+  - [References on jump time stepping methods](@ref jump_solve)
+  - [FAQ with information on changing parameters between simulations and using callbacks](@ref FAQ)
+  - [API documentation](@ref JumpProcesses.jl-API)