Merge pull request #701 from SciML/PEtab_tutorial_2

Inverse problem doc remake and new petab tutorial
SciML · Nov 7, 2023 · 507ae86 · 507ae86
2 parents 7a53667 + c3550a9
commit 507ae86
Show file tree

Hide file tree

Showing 13 changed files with 739 additions and 15 deletions.
diff --git a/docs/Project.toml b/docs/Project.toml
@@ -1,17 +1,21 @@
 [deps]
 BifurcationKit = "0f109fa4-8a5d-4b75-95aa-f515264e7665"
 Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
+DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 HomotopyContinuation = "f213a82b-91d6-5c5d-acf7-10f1c761b327"
 Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
+Optim = "429524aa-4258-5aef-a3af-852621145aeb"
 Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"
 OptimizationOptimisers = "42dfb2eb-d2b4-4451-abcd-913932933ac1"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+PEtab = "48d54b35-e43e-4a66-a5a1-dde6b987cf69"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
+SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SciMLSensitivity = "1ed8b502-d754-442c-8d5d-10ac956f44a1"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
@@ -21,17 +25,21 @@ Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 [compat]
 BifurcationKit = "0.3"
 Catalyst = "13"
+DataFrames = "1"
 DifferentialEquations = "7.7"
 Distributions = "0.25"
 Documenter = "0.27"
 HomotopyContinuation = "2.6"
 Latexify = "0.15, 0.16"
 ModelingToolkit = "8.47"
 NonlinearSolve = "1.6.0, 2"
+Optim = "1"
 Optimization = "3.19"
 OptimizationOptimisers = "0.1.1"
 OrdinaryDiffEq = "6"
+PEtab = "2"
 Plots = "1.36"
+SciMLBase = "~2.5"
 SciMLSensitivity = "7.19"
 Setfield = "1.1"
 SpecialFunctions = "2.1"

diff --git a/docs/pages.jl b/docs/pages.jl
@@ -1,19 +1,20 @@
 pages = Any["Home" => "index.md",
             "Introduction to Catalyst" => Any["introduction_to_catalyst/catalyst_for_new_julia_users.md",
-                                              "introduction_to_catalyst/introduction_to_catalyst.md"],
+                                              "introduction_to_catalyst/introduction_to_catalyst.md" ],
             "Catalyst Functionality" => Any["catalyst_functionality/dsl_description.md",
                                             "catalyst_functionality/programmatic_CRN_construction.md",
                                             "catalyst_functionality/compositional_modeling.md",
                                             "catalyst_functionality/constraint_equations.md",
                                             "catalyst_functionality/parametric_stoichiometry.md",
-                                            "catalyst_functionality/network_analysis.md"],
+                                            "catalyst_functionality/network_analysis.md",
+                                            "Model creation examples" => Any["catalyst_functionality/example_networks/basic_CRN_examples.md",
+                                                                             "catalyst_functionality/example_networks/hodgkin_huxley_equation.md",
+                                                                             "catalyst_functionality/example_networks/smoluchowski_coagulation_equation.md"]],
             "Catalyst Applications" => Any["catalyst_applications/simulation_structure_interfacing.md",
                                            "catalyst_applications/advanced_simulations.md",
                                            "catalyst_applications/homotopy_continuation.md",
-                                           "catalyst_applications/bifurcation_diagrams.md",
-                                           "catalyst_applications/parameter_estimation.md"],
-            "Example Networks" => Any["example_networks/basic_CRN_examples.md",
-                                      "example_networks/hodgkin_huxley_equation.md",
-                                      "example_networks/smoluchowski_coagulation_equation.md"],
+                                           "catalyst_applications/bifurcation_diagrams.md"],
+            "Inverse Problems" => Any["inverse_problems/parameter_estimation.md",
+                                      "inverse_problems/petab_ode_param_fitting.md"],
             "FAQs" => "faqs.md",
-            "API" => "api/catalyst_api.md"]
+            "API" => "api/catalyst_api.md"]
diff --git a/docs/src/assets/petab_best_objective.svg b/docs/src/assets/petab_best_objective.svg
diff --git a/docs/src/assets/petab_waterfall.svg b/docs/src/assets/petab_waterfall.svg
diff --git a/docs/src/catalyst_applications/simulation_structure_interfacing.md b/docs/src/catalyst_applications/simulation_structure_interfacing.md
@@ -99,7 +99,7 @@ sol[:k1]
 ```
 Finally, we note that we cannot change the values of solution states or parameters (i.e. both `sol[:X1] = 0.0` and `sol[:k1] = 0.0` generate errors).
 
-## Interfacing using symbolic representation
+## [Interfacing using symbolic representation](@id simulation_structure_interfacing_symbolic_representation)
 
 Catalyst is built on an *intermediary representation* implemented by (ModelingToolkit.jl)[https://github.com/SciML/ModelingToolkit.jl]. ModelingToolkit is a modelling framework where one first declares a set of symbolic variables and parameters using e.g.
 ```@example ex2

diff --git a/docs/src/catalyst_functionality/constraint_equations.md b/docs/src/catalyst_functionality/constraint_equations.md
@@ -85,7 +85,7 @@ sol = solve(oprob, Tsit5())
 plot(sol)
 ```
 
-## Adding events
+## [Adding events](@id constraint_equations_events)
 Our current model is unrealistic in assuming the cell will grow exponentially
 forever. Let's modify it such that the cell divides in half every time its
 volume reaches a size of `2`. We also assume we lose half of the protein upon

diff --git a/docs/src/catalyst_functionality/dsl_description.md b/docs/src/catalyst_functionality/dsl_description.md
@@ -412,7 +412,7 @@ sol = solve(oprob)
 plot(sol)
 ```
 
-## Setting initial conditions that depend on parameters
+## [Setting initial conditions that depend on parameters](@id dsl_description_parametric_initial_conditions)
 It is possible to set the initial condition of one (or several) species so that they depend on some system parameter. This is done in a similar way as default initial conditions, but giving the parameter instead of a value. When doing this, we also need to ensure that the initial condition parameter is a variable of the system:
 ```@example tut2
 rn = @reaction_network begin

diff --git a/...rc/example_networks/basic_CRN_examples.md → ...ty/example_networks/basic_CRN_examples.md b/...rc/example_networks/basic_CRN_examples.md → ...ty/example_networks/basic_CRN_examples.md
@@ -89,7 +89,7 @@ jprob = JumpProblem(rs, dprob, Direct())
 jsol = solve(jprob, SSAStepper())
 
 plot(plot(osol; title = "Reaction Rate Equation ODEs"),
-     plot(jsol; title = "Gillespie Jump Simulation");
+     plot(jsol; title = "Stochastic Chemical Kinetics Jump Processes");
      layout = (2, 1))
 ```
 

diff --git a/...ample_networks/hodgkin_huxley_equation.md → ...ample_networks/hodgkin_huxley_equation.md b/...ample_networks/hodgkin_huxley_equation.md → ...ample_networks/hodgkin_huxley_equation.md
diff --git a/...orks/smoluchowski_coagulation_equation.md → ...orks/smoluchowski_coagulation_equation.md b/...orks/smoluchowski_coagulation_equation.md → ...orks/smoluchowski_coagulation_equation.md
@@ -127,7 +127,7 @@ plot!(ϕ, sol[3,:]/Nₒ, line = (:dot,4,:purple), label="Analytical sol--X3",
 ```
 For the **additive kernel** we find
 
-![additive_kernel](../assets/additive_kernel.svg)
+![additive_kernel](../../assets/additive_kernel.svg)
 
 ---
 ## References

diff --git a/docs/src/catalyst_functionality/parametric_stoichiometry.md b/docs/src/catalyst_functionality/parametric_stoichiometry.md
@@ -72,12 +72,12 @@ oprob = ODEProblem(osys, u₀, (0.0, 1.0), p)
 nothing # hide
 ```
 We can now solve and plot the system
-```@example s1
+```@julia
 sol = solve(oprob, Tsit5())
 plot(sol)
 ```
 *If we had used a vector to store parameters, `m` and `n` would be converted to
-floating point giving an error when solving the system.*
+floating point giving an error when solving the system.* **Note, currently a [bug](https://github.com/SciML/ModelingToolkit.jl/issues/2296) in ModelingToolkit has broken this example by converting to floating point when using tuple parameters, see the alternative approach below for a workaround.**
 
 An alternative approach to avoid the issues of using mixed floating point and
 integer variables is to disable the rescaling of rate laws as described in

diff --git a/...lyst_applications/parameter_estimation.md → .../inverse_problems/parameter_estimation.md b/...lyst_applications/parameter_estimation.md → .../inverse_problems/parameter_estimation.md