From 0aa6921b8d508f21257a2e6577bcbb7a89d6dfc1 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Mon, 19 Oct 2026 16:45:21 +0200 Subject: [PATCH 1/8] add GlobalSensitivity to `Project.toml` --- exercises/solved_notebooks/Project.toml | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/exercises/solved_notebooks/Project.toml b/exercises/solved_notebooks/Project.toml index e3eaa3c..d9953e8 100644 --- a/exercises/solved_notebooks/Project.toml +++ b/exercises/solved_notebooks/Project.toml @@ -1,6 +1,7 @@ [deps] Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" +GlobalSensitivity = "af5da776-676b-467e-8baf-acd8249e4f0f" InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240" JumpProcesses = "ccbc3e58-028d-4f4c-8cd5-9ae44345cda5" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" @@ -18,6 +19,6 @@ StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" Turing = "fce5fe82-541a-59a6-adf8-730c64b5f9a0" [compat] -Turing = "0.41" Catalyst = "15" -julia = "1.12" \ No newline at end of file +Turing = "0.41" +julia = "1.12" From 6465a283caac08a7e4de73ac69be25307eb047a0 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Mon, 19 Oct 2026 16:45:32 +0200 Subject: [PATCH 2/8] add GSA to intro notebook --- .../P7_sens_uncert/sens_intro.jl | 397 +++++++++--------- 1 file changed, 199 insertions(+), 198 deletions(-) diff --git a/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl b/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl index 3fa7fac..f89be20 100644 --- a/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl +++ b/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl @@ -1,5 +1,5 @@ ### A Pluto.jl notebook ### -# v0.20.21 +# v1.0.1 using Markdown using InteractiveUtils @@ -7,8 +7,10 @@ using InteractiveUtils # ╔═╡ e1e7bc8e-7264-4cbc-98d2-aa73679fa2df begin # add this cell if you want the notebook to use the environment from where the Pluto server is launched - using Pkg - Pkg.activate("..") + using Pkg; Pkg.activate("..") + using Catalyst + using OrdinaryDiffEq, StatsPlots + using PlutoUI; TableOfContents() end # ╔═╡ 5f4fea06-0632-11ef-102e-21f5606d2056 @@ -17,15 +19,12 @@ using Markdown # ╔═╡ 489b5399-fe4c-481a-834f-0101bbe28cea using InteractiveUtils -# ╔═╡ 22489bd4-ab64-4bf9-ad03-5372ea273935 -using Catalyst - -# ╔═╡ cdcebbb1-40e0-457f-a6ec-b769f6b1f2e9 -using OrdinaryDiffEq, StatsPlots - # ╔═╡ 65885bbe-7b73-4efd-a8e7-984a34b80548 using ForwardDiff +# ╔═╡ a493c9e6-e250-4b7b-8cf2-03e3ac40cf19 +using GlobalSensitivity + # ╔═╡ 427b509f-d08e-4d93-99ba-a79f9c244b28 md""" # Introduction to sensitivity analysis @@ -33,9 +32,12 @@ md""" # ╔═╡ 49e7085c-3691-4182-9630-68dc9371ad18 md""" -## Goal of this practicum +## Goal of this practical """ +# ╔═╡ faf316bb-0ccc-4471-a53b-2476c3646705 +md"### Sensitivity functions" + # ╔═╡ c41894dd-0f6c-483f-b19f-dbf8148f776f md""" Sensitivity functions indicate how sensitive the model output is to a change in parameter values. When a model output is very sensitive to a certain parameter, a small change in the value of this parameter will have a large influence on the value of the model output. Sensitivity functions thus provide important information about the model and are implicitely used to estimate parameters and explicitely in the context of optimal experimental design. @@ -48,19 +50,36 @@ The sensitivity function that measures how sensitive output $y_j$ is to changes $$\cfrac{\partial \hat{y}_i(\theta)}{\partial \theta_j}\tag{1}$$ """ -# ╔═╡ abd70f04-1bf8-454f-869a-0d8082d68453 +# ╔═╡ ec1ee750-1f35-4f51-a8c7-d42a18a7b233 md""" Try to understand why the above expression $(1)$ does indeed give us the information we were promised in the first paragraph. How will we be able to see from the value determined by the above expression $(1)$ that whether or not output $y_j$ is sensitive to a change in $\theta_i$? """ -# ╔═╡ e76aedec-69f2-4301-b015-6960e4503c42 +# ╔═╡ f4ee83d1-259e-4dea-a86c-8c909ae4ec61 md""" -Answer: missing +!!! example + Consider the function $f(x) = a \, e^{b \, x}$, with parameters $\vec{θ} = (a, b)$. + + The sensitivity to $a$ is + + ```math + \dfrac{\partial f}{\partial a} = e^{b \, x} \, , + ``` + and for $b$: + + ```math + \dfrac{\partial f}{\partial b} = a \, x \, e^{b \, x} \, . + ``` + + If we now take a set of values for the input and parameters, such as $(x = 2, a = 1, b = 1)$, we can calculate the sensitivities at this point: + - For $a$: $\frac{\partial f}{\partial a} = e^2$ + - For $b$: $\frac{\partial f}{\partial b} = 2 e^2$ + from which we learn the function is twice as sensitive to $b$ compared to $a$ here. """ # ╔═╡ a6a97279-0b64-45c3-8312-c22b1a8425d0 md""" -Sometimes expression $(1)$ can be evaluated analytically. Usually, however, we will have to approximate the partial derivative numerically. Expression $(1)$ can be made more specific: +For simple functions, $(1)$ can be evaluated analytically. Usually, however, we will have to approximate the partial derivative numerically. Expression $(1)$ can be made more specific: $$\cfrac{\partial \hat{y}_i(\theta)}{\partial \theta_j} \approx \cfrac{\hat{y}_i(\theta_j+\Delta\theta_j)-\hat{y}_i(\theta_j)}{\Delta\theta_j}$$ """ @@ -70,6 +89,9 @@ md""" Thus, to calculate the sensitivity function numerically, the model is evaluated for the parameter values $\theta_i$ and $\theta_i+\Delta\theta_i$ and the difference between these evaluations is taken. """ +# ╔═╡ 33a71f15-dc39-45ed-b26c-855f42c55cae +md"### Normalized sensitivity" + # ╔═╡ 2a544bef-d7fe-4300-bc59-69e6f7429304 md""" Since quantity $(1)$ is dependent on the units, a normalized variant is often used: @@ -81,9 +103,14 @@ The interpretation of $(2)$ is how much the output changes per cent if the param Using the normalized variant allows you to compare all possible sensitivity functions with each other. """ -# ╔═╡ 7107a3a9-15ef-488e-8709-52a6444d0e1c +# ╔═╡ ec590df1-91be-40a2-844e-4ec3ea5b1077 +md"### Global sensitivity analysis" + +# ╔═╡ 3e2b5d93-4385-4261-9486-d93511ee7c26 md""" -We now calculate and interpret sensitivity functions for some given models. To illustrate the concepts, we first consider three simple models modelling the growth of grass. +The sensitivity function evaluates the partial derivative of your function **given some values for your parameters** (see also the given example). Therefore, the results of your sensitivity analysis will change based on the value of θ. As these methods are dependent on the current location of θ, they are called **local sensitivity analysis** methods. + +**Global sensitivity methods** try to circumvent this problem by quantifying sensitivity using **multiple** values of the output $y$. We consider two methods in this course: the Sobol method and the Morris method. For more information regarding these methods, we refer to the theory. """ # ╔═╡ 04e95855-7c4a-4d2c-b836-c5dde291adad @@ -107,9 +134,17 @@ with output $W$ the grass yield, and $W_f$, $\mu$ and $D$ parameters. The table | Exponential | 0.02 | 10.0 | | | Gompertz | 0.09 | | 0.04 | -We will use an initial condition of $W_0 = 2.0$ for each and a simulation time of $100$ days. +For all simulations (regardless of model), we will use an initial condition of $W_0 = 2.0$, a simulation time of $100$ days and save our solution every $0.5$ days. """ +# ╔═╡ 0b3d35bb-c5b0-44b7-94b3-06fa571d339e +begin + const u0 = [:W => 2.0] + const tspan = (0.0, 100.0) + const Δt = 0.5 + const t_vals = tspan[1]:Δt:tspan[2] # time values of the simulation +end; + # ╔═╡ 8d55bc42-1f8f-4bc1-b9fa-e8152fc125ce md""" We will illustrate how to compute the local sensitivity functions for the logistic model. The same will be left as exercises below for the exponential and gompertz models. @@ -145,32 +180,21 @@ Convert the *reaction model* to check that we work with the correct differential """ # ╔═╡ ba181db8-d176-4b83-9168-d5939ffe9661 -osys_log = convert(ODESystem, growth_log) - -# ╔═╡ cd8a7ba1-194a-4b42-8c87-2c9b1fe6b475 -md""" -Setting initial conditions, timespan and parameter values: -""" - -# ╔═╡ cbb2ca49-b019-495b-9310-83fcc00cad26 -u0_log = [:W => 2.0] - -# ╔═╡ 0b3d35bb-c5b0-44b7-94b3-06fa571d339e -tspan = (0.0, 100.0) # this will be the same for the three models +osys_log = convert(ODESystem, growth_log) # ╔═╡ be565a3c-31b6-4df1-b73b-08f308a8c09b md""" For the sake of clarity, we will use the variables `μ_log` and `Wf_log` to store the parameter values. """ -# ╔═╡ f62806d1-77e1-470b-9711-33a924c788cc -μ_log = 0.07 +# ╔═╡ dca01fd2-28f9-4d36-9b8f-b209ab97d3db +μ_log = 0.07; -# ╔═╡ 546ed163-26a7-4235-982f-7568ed609488 -Wf_log = 10.0 +# ╔═╡ 292aff06-7dfa-424b-8a90-173c37b9e67b +Wf_log = 10.0; # ╔═╡ 0da53fa2-5a42-46e6-8bd8-45d6aa903d46 -params_log = [:μ => μ_log, :Wf => Wf_log] +params_log = [:μ => μ_log, :Wf => Wf_log]; # ╔═╡ b1298f40-4696-49d0-ac94-896e0cdbc996 md""" @@ -178,12 +202,12 @@ Creating and solving the ODEProblem and plotting results: """ # ╔═╡ 270647d2-1371-4272-8bc1-3a6ad77bc716 -oprob_log = ODEProblem(growth_log, u0_log, tspan, params_log) +oprob_log = ODEProblem(growth_log, u0, tspan, params_log); # Also possible: # oprob_log = ODEProblem(growth_mod_log, [], tspan, []) # ╔═╡ ac235d86-1d93-4944-aa89-1b4fd38f0e6e -osol_log = solve(oprob_log, Tsit5(), saveat=0.5) +osol_log = solve(oprob_log, Tsit5(), saveat = Δt); # ╔═╡ 3e13efa1-9bc6-456f-8e62-ecd3165e2a65 plot(osol_log) @@ -206,30 +230,18 @@ We need to write a solution function with as argument a vector of the parameters # ╔═╡ bc201b61-0f58-49d2-a20b-4f18c42fcc96 function growth_sim_log(params) μ, Wf = params - u0_log = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_log = ODEProblem(growth_log, u0_log, tspan, [:μ=>μ, :Wf=>Wf]) - osol_log = solve(oprob_log, Tsit5(), saveat=0.5) + oprob_log = ODEProblem(growth_log, u0, tspan, [:μ=>μ, :Wf=>Wf]) + osol_log = solve(oprob_log, Tsit5(), saveat = Δt) return osol_log end # ╔═╡ 11f15a0e-2958-4075-928c-5a24bcc00c69 md""" -Next we will need to make a function based on the solution function that returns a single output. +Next we will need to make a function based on the solution function that returns the output of a single variable. """ # ╔═╡ d0f3197f-3094-4f6e-b33c-e5b74e0947f6 -growth_sim_W_log(params) = growth_sim_log(params)[:W] - -# ╔═╡ ba4eaa90-b97a-4b71-b3e1-968340c6def5 -md""" -Now make a time vector that is the same as the time vector from the solution. -""" - -# ╔═╡ 71f897a6-1e91-412b-ad58-c5f1e7cd1adb -t_vals_log = 0:0.5:100.0 -# Alternative: -# t_vals_log = tspan[1]:0.5:tspan[2] +growth_sim_W_log(params) = growth_sim_log(params)[:W]; # ╔═╡ b3b97b57-a4f7-4cb6-80c8-6ad893ade75d md""" @@ -237,7 +249,7 @@ Compute the single output with the given parameter values. """ # ╔═╡ 99864b21-3db0-4d60-a84e-8e96db4de4ae -W_log = growth_sim_W_log([μ_log, Wf_log]) +W_log = growth_sim_W_log([μ_log, Wf_log]); # ╔═╡ 344eb3a7-8614-4349-a536-62d9acef6bda md""" @@ -245,7 +257,7 @@ Use the function `ForwardDiff.jacobian` to compute the sensitivities. This funct """ # ╔═╡ 670ccfa4-24fa-4167-a46a-2a48dc19538b -sens_W_log = ForwardDiff.jacobian(growth_sim_W_log, [μ_log, Wf_log]) +sens_W_log = ForwardDiff.jacobian(growth_sim_W_log, [μ_log, Wf_log]); # ╔═╡ c262d3ed-9f62-45ad-a7ee-ec78ab16f35c md""" @@ -255,10 +267,10 @@ To get the sensitivities of $W$ on $\mu$, and of $W$ on $W_f$, you need to use i """ # ╔═╡ 6600d28e-2522-4069-a5e0-643be43f6117 -sens_W_on_μ_log = sens_W_log[:,1] # sensitivity of W on μ +sens_W_on_μ_log = sens_W_log[:,1]; # sensitivity of W on μ # ╔═╡ d9e3a4ac-c138-4b93-a19c-737166a0f0ea -sens_W_on_Wf_log = sens_W_log[:,2] # sensitivity of W on Wf +sens_W_on_Wf_log = sens_W_log[:,2]; # sensitivity of W on Wf # ╔═╡ f07d8205-e172-45fa-b955-792bd95f3023 md""" @@ -266,50 +278,81 @@ We now calculate the normalized sensitivities. For that we need to multiply by t """ # ╔═╡ 1c142797-ada6-4d77-a22f-b981ffd38956 -sens_W_on_μ_rel_log = sens_W_on_μ_log .* μ_log ./ W_log +sens_W_on_μ_rel_log = sens_W_on_μ_log .* μ_log ./ W_log; # ╔═╡ cb397e0b-56f9-420f-a86e-bfab35286b44 -sens_W_on_Wf_rel_log = sens_W_on_Wf_log .* Wf_log ./ W_log +sens_W_on_Wf_rel_log = sens_W_on_Wf_log .* Wf_log ./ W_log; # ╔═╡ 9e6946a4-f207-4bfb-9e19-aa7b77c2a05b md""" -We are now ready to plot the two sensitivity functions. We provide the time vector (first argument) and a vector of the two sensitivity functions (second argument). Additionally, you can provide a title, (legend) labels and a x- and/or y-label. +We can now plot the two sensitivity functions as follows: """ -# ╔═╡ f58ce914-f366-4087-a7d1-8cfe69ac623b -plot(t_vals_log, [sens_W_on_μ_rel_log, sens_W_on_Wf_rel_log], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") +# ╔═╡ dc0507ee-0f67-4556-bdc3-1294177c86b7 +begin + plot( + t_vals, sens_W_on_μ_rel_log, title = "Normalized sensitivities", + label = "W on μ", xlabel = "Time (day)" + ) + plot!(t_vals, sens_W_on_Wf_rel_log, label="W on Wf") +end -# ╔═╡ 05972c7f-f64f-4865-b0a7-f33029d0a6fa +# ╔═╡ 04a1d1ad-53a9-4146-be09-65aa422e3730 md""" -Notice that in the `label` option there is no comma separating the labels. +Conclusions: +- From the sensitivity plot of $W$ on $\mu$ it can be seen that $W$ is most sensitive to $\mu$ in the time region $[0, 30]$. The latter corresponds to the region where the yield rate is largest (i.e., when the growth is largest). This makes sense because when looking at the differential equation, $\mu$ is approximately the growth rate for relatively small $W$ values. +- From the sensitivity plot of $W$ on $W_f$ it can be seen that $W$ is most sensitive to $W_f$ in the region where time values are large (cf. operating point). The latter corresponds to the region where the yield rate stagnates (i.e., when the yield reaches a steady value). This makes sense because when looking at the differential equation, $W_f$ is the steady state value. """ -# ╔═╡ 88aa13f6-e16c-40cd-ac10-3ede8b9fb429 -# ╠═╡ disabled = true -#=╠═╡ +# ╔═╡ 5152abfc-98fa-4e10-aa28-6455b996a256 +md"### Global sensitivity analysis (GSA)" + +# ╔═╡ 8af97e64-6a14-4c7d-b935-002af339c94a md""" -#### Alternative +To perform global sensitivity analysis, we use the `gsa` function from the `GlobalSensitivity` package. This function requires 4 inputs: +1. Similarly to local sensitivity analysis, we require a function that takes the model parameters as its input and returns a model output of interest. We can just reuse `growth_sim_log` here to get the model outputs at all timepoints. +1. A GSA method. The methods considered in this course are the Sobol method `Sobol()` and the Morris method `Morris()`. We will use Sobol here. +1. A range for every parameter from where values will be sampled. For our grass growth models, we will assume realistic values for $μ$ are $[0, 0.2]$ and for $W_f$ are $[1, 20]$. +1. The amount of samples. As the Sobol method is computationally expensive, you typically don't want to use too many for an exercise. For the Morris method on the other hand, you typically want a large amount of samples for stable results (> `1_000_000`). +""" + +# ╔═╡ b86d3fe4-c3e5-4c39-8118-7970d609c17a +p_range = [[0.0, 0.2], [1.0, 20.0]]; -For each plot command we need to provide the time vector (first argument) and the column vector with the local sensitivity (second argument). Additionally, you can provide a title, a (legend) label and a x- and/or y-label. +# ╔═╡ 66083965-8402-47b6-9ccd-9ced9c7412e7 +gsa_res = gsa(growth_sim_log, Sobol(), p_range, samples = 100) -Beware that if you want to execute multiple commands (here: `plot`) in a single cell, you need to put them in a `begin`-and-`end` block. Also if you want to visualize subsequent graphs in the same plot, the forthcoming plot command names should be followed by a `!`, as in `plot!(...)`. +# ╔═╡ 97eed321-4dfe-468b-8b4e-6748608a781d +md""" +!!! tip + Does a function return some weird type you don't know what to do with? You can find out the fields of a variable `x` using `fieldnames(typeof(x))`. """ - ╠═╡ =# -# ╔═╡ dc0507ee-0f67-4556-bdc3-1294177c86b7 -# ╠═╡ disabled = true -#=╠═╡ +# ╔═╡ 8c923fea-2e75-42c1-ae73-0bdf8d1bee91 +fieldnames(typeof(gsa_res)) + +# ╔═╡ 3016e5f5-2f13-4877-9ed8-e6a355e63d5e +md""" +The resulting variable contains the first, second and total order Sobol indices, which we can access using respectively `gsa_res.S1`, `gsa_res.S2` and `gsa_res.ST`. We will plot the total effects here: +""" + +# ╔═╡ e2ec4060-342e-41c6-8bb2-1a20192ceeff begin -plot(t_vals_log, sens_W_on_μ_rel_log, title="Normalized sensitivities", label="W on μ", xlabel="Time (day)") -plot!(t_vals_log, sens_W_on_Wf_rel_log, label="W on Wf") + plot(t_vals, gsa_res.ST[1][1, :], label = "μ", xlabel = "Time (day)", + title = "Global sensitivities", ylims = (0, 1) + ) + plot!(t_vals, gsa_res.ST[2][1, :], label = "Wf") end - ╠═╡ =# -# ╔═╡ 04a1d1ad-53a9-4146-be09-65aa422e3730 +# ╔═╡ 73990cc1-77d5-4f09-ab2e-fb67c7ff81a4 md""" -Conclusions: -- From the sensitivity plot of $W$ on $\mu$ it can be seen that $W$ is most sensitive to $\mu$ in the time region $[0, 30]$. The latter corresponds to the region where the yield rate is largest (i.e., when the growth is largest). This makes sense because when looking at the differential equation, $\mu$ is approximately the growth rate for relatively small $W$ values. -- From the sensitivity plot of $W$ on $W_f$ it can be seen that $W$ is most sensitive to $W_f$ in the region where time values are large (cf. operating point). The latter corresponds to the region where the yield rate stagnates (i.e., when the yield reaches a steady value). This makes sense because when looking at the differential equation, $W_f$ is the steady state value. +!!! question + GSA returns an appreciably positive value for the effect of growth rate μ on final grass yield, while it was approximately $0$ in LSA. Why is this effect larger when we perform a global analysis where we consider multiple possible values for $μ$ and $W_f$? +""" + +# ╔═╡ 51488d93-48ad-4a02-9b51-00ef26f64699 +md""" +**Answer**: missing """ # ╔═╡ 79f0f1dd-850f-4dfd-b895-ff41a2d8adb8 @@ -352,30 +395,16 @@ Convert the system to a symbolic differential equation model (name it: `osys_exp # osys_exp = missing # Uncomment and complete the instruction osys_exp = convert(ODESystem, growth_exp) -# ╔═╡ 63952c54-3304-4500-9536-b375c5c8f280 -md""" -Initialize a vector `u0_exp` with the initial condition: -""" - -# ╔═╡ 7dc013f7-7b35-42e4-aa77-0dca3755f389 -# u0_exp = missing # Uncomment and complete the instruction -u0_exp = [:W => 2.0] - -# ╔═╡ 86785c55-b567-449c-ae76-bc15a16223bc -md""" -We will use the same timespan as before, so no need to redefine it. -""" - # ╔═╡ 05e96232-ce95-452b-98e4-e79817d45ae2 md""" For the sake of clarity, we will use the variables `μ_exp` and `Wf_exp` to store the parameter values. """ # ╔═╡ 9a705b3f-3c45-4896-9b89-e9c37d47fd89 -μ_exp = 0.02 +μ_exp = 0.02; # ╔═╡ 55a13ebc-09e7-433c-83a8-9b6805e2fbec -Wf_exp = 10.0 +Wf_exp = 10.0; # ╔═╡ f1293675-b4d0-4827-8714-65a593ecc00e md""" @@ -384,16 +413,16 @@ Initialize a vector `params_exp` with the parameter values: # ╔═╡ 7005d88f-9e28-47fe-9a32-dec023903dc3 # params_exp = missing # Uncomment and complete the instruction -params_exp = [:μ => μ_exp, :Wf => Wf_exp] +params_exp = [:μ => μ_exp, :Wf => Wf_exp]; # ╔═╡ c8853554-26c5-491c-a536-4c07e8c6a986 md""" -Create the ODE problem and store it in `oprob_exp`: +Create the ODE problem and store it in `oprob_exp`, using the same initial condition and timespan as before. """ # ╔═╡ 025c1154-b2ae-4e1c-af1b-277b24d648a4 # oprob_exp = missing # Uncomment and complete the instruction -oprob_exp = ODEProblem(growth_exp, u0_exp, tspan, params_exp) +oprob_exp = ODEProblem(growth_exp, u0, tspan, params_exp); # ╔═╡ f58825d2-55fc-43f1-b164-9555bf9f5b84 md""" @@ -402,7 +431,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os # ╔═╡ 8044dd5f-9996-4062-a54a-8fb0eee30b46 # osol_exp = missing # Uncomment and complete the instruction -osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5) +osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5); # ╔═╡ a1752ad4-2c11-416c-8321-9da058f9aaea md""" @@ -426,30 +455,19 @@ Write a solution function with as argument a vector of the parameters (that you # end function growth_sim_exp(params) μ, Wf = params - u0_exp = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_exp = ODEProblem(growth_exp, u0_exp, tspan, [:μ=>μ, :Wf=>Wf]) + oprob_exp = ODEProblem(growth_exp, u0, tspan, [:μ=>μ, :Wf=>Wf]) osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5) return osol_exp end # ╔═╡ b31b4d24-ae30-4cff-b6b7-348200feaa5d md""" -Make a function based on the solution function that returns a single output. +Make a function based on the solution function that returns the output of a single variable. """ # ╔═╡ f2fa716f-d7b5-4c95-a8c6-6b16fbfe1499 # growth_sim_W_exp(params) = missing # Uncomment and complete the instruction -growth_sim_W_exp(params) = growth_sim_exp(params)[:W] - -# ╔═╡ 0b413b98-a800-44e9-a917-ba61467bd613 -md""" -Make the time vector. -""" - -# ╔═╡ 6be6fea7-561b-4eec-b249-fa4522a5b039 -# t_vals = missing # Uncomment and complete the instruction -t_vals_exp = 0:0.5:100.0 +growth_sim_W_exp(params) = growth_sim_exp(params)[:W]; # ╔═╡ 98070d41-e8c5-49e6-9eea-04c718cbff65 md""" @@ -458,7 +476,7 @@ Compute the output for the given parameter values. # ╔═╡ cef751aa-f9c8-46da-9436-9cc2e5d7515c # W_exp = missing # Uncomment and complete the instruction -W_exp = growth_sim_W_exp([μ_exp, Wf_exp]) +W_exp = growth_sim_W_exp([μ_exp, Wf_exp]); # ╔═╡ 9e7bb7ec-8fe3-422e-920f-13c2ef055feb md""" @@ -467,7 +485,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s # ╔═╡ b13b5e77-ca60-465b-9f63-be9e5da0482f # sens_W_exp = missing # Uncomment and complete the instruction -sens_W_exp = ForwardDiff.jacobian(growth_sim_W_exp, [μ_exp, Wf_exp]) +sens_W_exp = ForwardDiff.jacobian(growth_sim_W_exp, [μ_exp, Wf_exp]); # ╔═╡ 735f5191-259e-440f-a379-20a2a70ec72c md""" @@ -476,11 +494,11 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. # ╔═╡ 689556c9-a690-4e9e-9052-7ac66e999d4d # sens_W_on_μ_exp = missing # Uncomment and complete the instruction -sens_W_on_μ_exp = sens_W_exp[:,1] +sens_W_on_μ_exp = sens_W_exp[:, 1]; # ╔═╡ b85f0d09-a102-435b-a898-149ba4b29caa # sens_W_on_Wf_exp = missing # Uncomment and complete the instruction -sens_W_on_Wf_exp = sens_W_exp[:,2] +sens_W_on_Wf_exp = sens_W_exp[:, 2]; # ╔═╡ d067f814-1dd5-4e59-8cfd-3bc0b2d97612 md""" @@ -489,11 +507,11 @@ Compute the normalized sensitivities. # ╔═╡ e01657f3-33ac-4eb8-b39e-e2441a165a8d # sens_W_on_μ_rel_exp = missing # Uncomment and complete the instruction -sens_W_on_μ_rel_exp = sens_W_on_μ_exp .* μ_exp ./ W_exp +sens_W_on_μ_rel_exp = sens_W_on_μ_exp .* μ_exp ./ W_exp; # ╔═╡ 7231e54e-b293-48fc-b1f0-273f45f51539 # sens_W_on_Wf_rel_exp = missing # Uncomment and complete the instruction -sens_W_on_Wf_rel_exp = sens_W_on_Wf_exp .* Wf_exp ./ W_exp +sens_W_on_Wf_rel_exp = sens_W_on_Wf_exp .* Wf_exp ./ W_exp; # ╔═╡ 0524a976-9d4d-4ef9-b088-bdd9102415e6 md""" @@ -502,7 +520,7 @@ Plot both normalized sensitivity functions (with appropriate title and labels): # ╔═╡ 33031b53-6f2d-4256-a0aa-0ba97beaa7af # missing # Uncomment and complete the instruction -plot(t_vals_exp, [sens_W_on_μ_rel_exp, sens_W_on_Wf_rel_exp], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") +plot(t_vals, [sens_W_on_μ_rel_exp, sens_W_on_Wf_rel_exp], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") # ╔═╡ 22b58c97-1aa2-49c4-8a0a-488c63014a90 md" @@ -512,6 +530,22 @@ Draw your conclusions: - missing " +# ╔═╡ ee17e112-a422-4b0e-be4c-c9e6abb040cb +md""" +Perform a global sensitivity analysis on this model, now using the **Morris method** and plot the results. You can use the same parameter range as before. +""" + +# ╔═╡ a8c930a2-a119-4681-a063-7a6c99615c5b +gsa_res_exp = gsa(growth_sim_exp, Morris(num_trajectory = 50), p_range, samples = 1_000_000); + +# ╔═╡ 1f9d7f84-fae7-4f0e-8783-7a8ce44051a3 +begin + plot(t_vals, gsa_res_exp.means[1][1, :], label = "μ", xlabel = "Time (day)", + title = "Global sensitivities" + ) + plot!(t_vals, gsa_res_exp.means[2][1, :], label = "Wf") +end + # ╔═╡ 41abf8c4-e67e-4f66-a38a-7204a878d98d md""" ### Exercise 2 - Sensitivity analysis of the Gompertz growth model @@ -547,48 +581,29 @@ Convert the system to a symbolic differential equation model (name it: `osys_gom # osys_gom = missing # Uncomment and complete the instruction osys_gom = convert(ODESystem, growth_gom) -# ╔═╡ ed1150f1-9c75-4867-b7f5-535d605810f4 -md""" -Initialize a vector `u0_gom` with the initial condition: -""" - -# ╔═╡ d35ef155-61aa-4b19-872f-d6e621f96572 -# u0_gom = missing # Uncomment and complete the instruction -u0_gom = [:W => 2.0] - -# ╔═╡ 3f546d58-15ff-442d-a4e8-7c0783ed2fe0 -md""" -We will use the same timespan as before, so no need to redefine it. -""" - -# ╔═╡ b778f786-87bc-4103-a792-2ae5cf004d50 -md""" -For the sake of clarity, we will use the variables `μ_gom` and `D_gom` to store the parameter values. -""" - -# ╔═╡ edd76a4e-e541-4a83-923f-6f593d8fadad -μ_gom = 0.09 - -# ╔═╡ 05f1ae6c-1229-4924-9c41-9554d6ca4261 -D_gom = 0.04 - # ╔═╡ 70df05e7-4527-4087-a762-942a33e89f74 md""" Initialize a vector `params_gom` with the parameter values: """ +# ╔═╡ e9408c0f-cbcb-4703-8b85-e349e37b2025 +μ_gom = 0.09; + +# ╔═╡ 24c4e348-7322-4c70-b8da-ec70d25356a9 +D_gom = 0.04; + # ╔═╡ bf9291d5-6fc5-4359-874e-167237a22147 # params_gom = missing # Uncomment and complete the instruction -params_gom = [:μ => μ_gom, :D => D_gom] +params_gom = [:μ => μ_gom, :D => D_gom]; # ╔═╡ 9e3158c5-d4b2-456c-b920-8431e2047770 md""" -Create the ODE problem and store it in `oprob_gom`: +Create the ODE problem and store it in `oprob_gom`, using the same initial conditions and timespan as before: """ # ╔═╡ 3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 # oprob_gom = missing # Uncomment and complete the instruction -oprob_gom = ODEProblem(growth_gom, u0_gom, tspan, params_gom) +oprob_gom = ODEProblem(growth_gom, u0, tspan, params_gom); # ╔═╡ 30c35371-d106-4673-8569-d07b39edbdfe md""" @@ -597,7 +612,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os # ╔═╡ c01fe052-3701-4e22-9a64-f681d52a445a # osol_gom = missing # Uncomment and complete the instruction -osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5) +osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5); # ╔═╡ 4f959d21-954b-4116-a498-bafc1aada47a md""" @@ -621,9 +636,7 @@ Write a solution function with as argument a vector of the parameters (that you # end function growth_sim_gom(params) μ, D = params - u0_gom = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_gom = ODEProblem(growth_gom, u0_gom, tspan, [:μ=>μ, :D=>D]) + oprob_gom = ODEProblem(growth_gom, u0, tspan, [:μ=>μ, :D=>D]) osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5) return osol_gom end @@ -635,16 +648,7 @@ Make a function based on the solution function that returns a single output. # ╔═╡ 69f96b57-7c8d-4e64-9e6a-9804d001b0b1 # growth_sim_W_gom(params) = missing # Uncomment and complete the instruction -growth_sim_W_gom(params) = growth_sim_gom(params)[:W] - -# ╔═╡ f7f6c56b-4910-4641-a462-eacfa4b4d034 -md""" -Make the time vector. -""" - -# ╔═╡ 84c1bb38-23c3-4f70-be67-7da2300a737b -# t_vals_gom = missing # Uncomment and complete the instruction -t_vals_gom = 0:0.5:100.0 +growth_sim_W_gom(params) = growth_sim_gom(params)[:W]; # ╔═╡ 09eeab39-79a1-420a-8298-72b7002ec168 md""" @@ -653,7 +657,7 @@ Compute the output for the given parameter values. # ╔═╡ 9e15601f-4ad2-41be-8168-db30548a1c3b # W_gom = missing # Uncomment and complete the instruction -W_gom = growth_sim_W_gom([μ_gom, D_gom]) +W_gom = growth_sim_W_gom([μ_gom, D_gom]); # ╔═╡ 2b731e4a-d6eb-48f2-91f0-27e916d85683 md""" @@ -662,7 +666,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s # ╔═╡ 68ad7e05-7f79-4e8d-9031-fefb5f5c0897 # sens_W_gom = missing # Uncomment and complete the instruction -sens_W_gom = ForwardDiff.jacobian(growth_sim_W_gom, [μ_gom, D_gom]) +sens_W_gom = ForwardDiff.jacobian(growth_sim_W_gom, [μ_gom, D_gom]); # ╔═╡ 362a8e24-8fac-4978-bc0e-13596a05d39e md""" @@ -671,11 +675,11 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. # ╔═╡ e9114d64-fe3e-420d-9fde-4b0fcbbc8527 # sens_W_on_μ_gom = missing # Uncomment and complete the instruction -sens_W_on_μ_gom = sens_W_gom[:,1] +sens_W_on_μ_gom = sens_W_gom[:,1]; # ╔═╡ 7c9e1112-0a63-417f-a65a-9c2bdfa8bc47 # sens_W_on_D_gom = missing # Uncomment and complete the instruction -sens_W_on_D_gom = sens_W_gom[:,2] +sens_W_on_D_gom = sens_W_gom[:,2]; # ╔═╡ 94e66042-b352-4bf2-a24b-b15789e10fe3 md""" @@ -684,11 +688,11 @@ Compute the normalized sensitivities. # ╔═╡ 1711346a-c161-415a-96f9-1235d786a584 # sens_W_on_μ_rel_gom = missing # Uncomment and complete the instruction -sens_W_on_μ_rel_gom = sens_W_on_μ_gom .* μ_gom ./ W_gom +sens_W_on_μ_rel_gom = sens_W_on_μ_gom .* μ_gom ./ W_gom; # ╔═╡ 7a87bfbb-a8c0-4071-8f86-66513ae40968 # sens_W_on_D_rel_gom = missing # Uncomment and complete the instruction -sens_W_on_D_rel_gom = sens_W_on_D_gom .* D_gom ./ W_gom +sens_W_on_D_rel_gom = sens_W_on_D_gom .* D_gom ./ W_gom; # ╔═╡ 7593a4bf-ad0f-4fd3-aa3c-8b6a5073cb02 md" @@ -697,7 +701,7 @@ Plot both sensitivity functions (with appropriate title and labels): # ╔═╡ ad70a1d5-9b68-4cd6-9259-a3c25cad706b # missing # Uncomment and complete the instruction -plot(t_vals_gom, [sens_W_on_μ_rel_gom, sens_W_on_D_rel_gom], title="Normalized sensitivities", label=["W on μ" "W on D"], xlabel="Time (day)") +plot(t_vals, [sens_W_on_μ_rel_gom, sens_W_on_D_rel_gom], title="Normalized sensitivities", label=["W on μ" "W on D"], xlabel="Time (day)") # ╔═╡ 12f64333-36d0-4e2a-9061-e3dcc7a4ae96 md" @@ -713,30 +717,29 @@ Draw your conclusions: # ╠═e1e7bc8e-7264-4cbc-98d2-aa73679fa2df # ╟─427b509f-d08e-4d93-99ba-a79f9c244b28 # ╟─49e7085c-3691-4182-9630-68dc9371ad18 +# ╟─faf316bb-0ccc-4471-a53b-2476c3646705 # ╟─c41894dd-0f6c-483f-b19f-dbf8148f776f # ╟─7b58f2f9-cb93-4a20-bf90-62a9006b57d6 -# ╟─abd70f04-1bf8-454f-869a-0d8082d68453 -# ╟─e76aedec-69f2-4301-b015-6960e4503c42 +# ╟─ec1ee750-1f35-4f51-a8c7-d42a18a7b233 +# ╟─f4ee83d1-259e-4dea-a86c-8c909ae4ec61 # ╟─a6a97279-0b64-45c3-8312-c22b1a8425d0 # ╟─25d3ef62-382c-455c-92bb-fadfc650c5a6 +# ╟─33a71f15-dc39-45ed-b26c-855f42c55cae # ╟─2a544bef-d7fe-4300-bc59-69e6f7429304 -# ╟─7107a3a9-15ef-488e-8709-52a6444d0e1c +# ╟─ec590df1-91be-40a2-844e-4ec3ea5b1077 +# ╟─3e2b5d93-4385-4261-9486-d93511ee7c26 # ╟─04e95855-7c4a-4d2c-b836-c5dde291adad # ╟─4673bfdf-8f4a-42bd-a026-21a66d800f2b +# ╠═0b3d35bb-c5b0-44b7-94b3-06fa571d339e # ╟─8d55bc42-1f8f-4bc1-b9fa-e8152fc125ce -# ╠═22489bd4-ab64-4bf9-ad03-5372ea273935 -# ╠═cdcebbb1-40e0-457f-a6ec-b769f6b1f2e9 # ╟─676119ce-f4fe-41f9-8121-b2c21f0dd28c # ╟─5492210b-6ee6-4baf-b37e-d21358cdeb60 # ╠═18eece66-46fa-4458-aaed-f4c8fa002c20 # ╟─8bb0b24b-6cce-49cd-a625-5f375b92d9b7 # ╠═ba181db8-d176-4b83-9168-d5939ffe9661 -# ╟─cd8a7ba1-194a-4b42-8c87-2c9b1fe6b475 -# ╠═cbb2ca49-b019-495b-9310-83fcc00cad26 -# ╠═0b3d35bb-c5b0-44b7-94b3-06fa571d339e # ╟─be565a3c-31b6-4df1-b73b-08f308a8c09b -# ╠═f62806d1-77e1-470b-9711-33a924c788cc -# ╠═546ed163-26a7-4235-982f-7568ed609488 +# ╠═dca01fd2-28f9-4d36-9b8f-b209ab97d3db +# ╠═292aff06-7dfa-424b-8a90-173c37b9e67b # ╠═0da53fa2-5a42-46e6-8bd8-45d6aa903d46 # ╟─b1298f40-4696-49d0-ac94-896e0cdbc996 # ╠═270647d2-1371-4272-8bc1-3a6ad77bc716 @@ -749,8 +752,6 @@ Draw your conclusions: # ╠═bc201b61-0f58-49d2-a20b-4f18c42fcc96 # ╟─11f15a0e-2958-4075-928c-5a24bcc00c69 # ╠═d0f3197f-3094-4f6e-b33c-e5b74e0947f6 -# ╟─ba4eaa90-b97a-4b71-b3e1-968340c6def5 -# ╠═71f897a6-1e91-412b-ad58-c5f1e7cd1adb # ╟─b3b97b57-a4f7-4cb6-80c8-6ad893ade75d # ╠═99864b21-3db0-4d60-a84e-8e96db4de4ae # ╟─344eb3a7-8614-4349-a536-62d9acef6bda @@ -762,20 +763,25 @@ Draw your conclusions: # ╠═1c142797-ada6-4d77-a22f-b981ffd38956 # ╠═cb397e0b-56f9-420f-a86e-bfab35286b44 # ╟─9e6946a4-f207-4bfb-9e19-aa7b77c2a05b -# ╠═f58ce914-f366-4087-a7d1-8cfe69ac623b -# ╟─05972c7f-f64f-4865-b0a7-f33029d0a6fa -# ╠═88aa13f6-e16c-40cd-ac10-3ede8b9fb429 # ╠═dc0507ee-0f67-4556-bdc3-1294177c86b7 # ╟─04a1d1ad-53a9-4146-be09-65aa422e3730 +# ╟─5152abfc-98fa-4e10-aa28-6455b996a256 +# ╠═a493c9e6-e250-4b7b-8cf2-03e3ac40cf19 +# ╟─8af97e64-6a14-4c7d-b935-002af339c94a +# ╠═b86d3fe4-c3e5-4c39-8118-7970d609c17a +# ╠═66083965-8402-47b6-9ccd-9ced9c7412e7 +# ╟─97eed321-4dfe-468b-8b4e-6748608a781d +# ╠═8c923fea-2e75-42c1-ae73-0bdf8d1bee91 +# ╟─3016e5f5-2f13-4877-9ed8-e6a355e63d5e +# ╠═e2ec4060-342e-41c6-8bb2-1a20192ceeff +# ╟─73990cc1-77d5-4f09-ab2e-fb67c7ff81a4 +# ╠═51488d93-48ad-4a02-9b51-00ef26f64699 # ╟─79f0f1dd-850f-4dfd-b895-ff41a2d8adb8 # ╟─fa2270d0-2548-409c-a23f-4369d8bce8ec # ╟─dca50b37-b06c-4efe-881e-cf966ebc8fd7 # ╠═2836f231-3b90-4829-9012-3ed9a09239ff # ╟─e8db84be-31a0-415d-a194-064c4c87a293 # ╠═78d175b4-a9b2-49b5-bbc8-eb348799985b -# ╟─63952c54-3304-4500-9536-b375c5c8f280 -# ╠═7dc013f7-7b35-42e4-aa77-0dca3755f389 -# ╟─86785c55-b567-449c-ae76-bc15a16223bc # ╟─05e96232-ce95-452b-98e4-e79817d45ae2 # ╠═9a705b3f-3c45-4896-9b89-e9c37d47fd89 # ╠═55a13ebc-09e7-433c-83a8-9b6805e2fbec @@ -791,8 +797,6 @@ Draw your conclusions: # ╠═eb1e5f97-3e09-43f0-b6ce-cc6cbc42ec2f # ╟─b31b4d24-ae30-4cff-b6b7-348200feaa5d # ╠═f2fa716f-d7b5-4c95-a8c6-6b16fbfe1499 -# ╟─0b413b98-a800-44e9-a917-ba61467bd613 -# ╠═6be6fea7-561b-4eec-b249-fa4522a5b039 # ╟─98070d41-e8c5-49e6-9eea-04c718cbff65 # ╠═cef751aa-f9c8-46da-9436-9cc2e5d7515c # ╟─9e7bb7ec-8fe3-422e-920f-13c2ef055feb @@ -806,18 +810,17 @@ Draw your conclusions: # ╟─0524a976-9d4d-4ef9-b088-bdd9102415e6 # ╠═33031b53-6f2d-4256-a0aa-0ba97beaa7af # ╟─22b58c97-1aa2-49c4-8a0a-488c63014a90 +# ╟─ee17e112-a422-4b0e-be4c-c9e6abb040cb +# ╠═a8c930a2-a119-4681-a063-7a6c99615c5b +# ╠═1f9d7f84-fae7-4f0e-8783-7a8ce44051a3 # ╟─41abf8c4-e67e-4f66-a38a-7204a878d98d # ╟─f311e943-c96b-4af3-8757-c662fd302a88 # ╠═1ae76036-42ef-46e9-88bf-d66d4267addc # ╟─e190ceca-30b1-49e2-baf5-ceb62929f4c0 # ╠═acaf5b91-accc-4dfa-9371-6495a47c8736 -# ╟─ed1150f1-9c75-4867-b7f5-535d605810f4 -# ╠═d35ef155-61aa-4b19-872f-d6e621f96572 -# ╟─3f546d58-15ff-442d-a4e8-7c0783ed2fe0 -# ╟─b778f786-87bc-4103-a792-2ae5cf004d50 -# ╠═edd76a4e-e541-4a83-923f-6f593d8fadad -# ╠═05f1ae6c-1229-4924-9c41-9554d6ca4261 # ╟─70df05e7-4527-4087-a762-942a33e89f74 +# ╠═e9408c0f-cbcb-4703-8b85-e349e37b2025 +# ╠═24c4e348-7322-4c70-b8da-ec70d25356a9 # ╠═bf9291d5-6fc5-4359-874e-167237a22147 # ╟─9e3158c5-d4b2-456c-b920-8431e2047770 # ╠═3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 @@ -829,8 +832,6 @@ Draw your conclusions: # ╠═be21e269-16c4-42b0-8930-2d009bd91161 # ╟─e57f5916-b753-429b-8434-8a85a67bd1fb # ╠═69f96b57-7c8d-4e64-9e6a-9804d001b0b1 -# ╟─f7f6c56b-4910-4641-a462-eacfa4b4d034 -# ╠═84c1bb38-23c3-4f70-be67-7da2300a737b # ╟─09eeab39-79a1-420a-8298-72b7002ec168 # ╠═9e15601f-4ad2-41be-8168-db30548a1c3b # ╟─2b731e4a-d6eb-48f2-91f0-27e916d85683 From 1bc4b06b796775f4ab9c86e32a8a8c4a13636bf7 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Mon, 19 Oct 2026 17:12:51 +0200 Subject: [PATCH 3/8] add more `# missing` to intro --- .../P7_sens_uncert/sens_intro.jl | 20 +++++++++++++------ 1 file changed, 14 insertions(+), 6 deletions(-) diff --git a/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl b/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl index f89be20..e11bcd3 100644 --- a/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl +++ b/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl @@ -371,7 +371,7 @@ Create a *reaction network object* for the exponential growth model. Name it `gr """ # ╔═╡ 2836f231-3b90-4829-9012-3ed9a09239ff -# Uncomment and complete the instruction +# # Uncomment and complete the instruction # growth_exp = @reaction_network begin # @species missing # @parameters missing @@ -448,7 +448,7 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ eb1e5f97-3e09-43f0-b6ce-cc6cbc42ec2f -# Uncomment and complete the instruction +# # Uncomment and complete the instruction # function growth_sim_exp(params) # missing # ... @@ -520,7 +520,11 @@ Plot both normalized sensitivity functions (with appropriate title and labels): # ╔═╡ 33031b53-6f2d-4256-a0aa-0ba97beaa7af # missing # Uncomment and complete the instruction -plot(t_vals, [sens_W_on_μ_rel_exp, sens_W_on_Wf_rel_exp], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") +begin + plot(t_vals, sens_W_on_μ_rel_exp, title = "Normalized sensitivities", + label = "W on μ", xlabel = "Time (day)") + plot!(t_vals, sens_W_on_Wf_rel_exp, label = "W on Wf") +end # ╔═╡ 22b58c97-1aa2-49c4-8a0a-488c63014a90 md" @@ -536,9 +540,11 @@ Perform a global sensitivity analysis on this model, now using the **Morris meth """ # ╔═╡ a8c930a2-a119-4681-a063-7a6c99615c5b +# gsa_res_exp = missing gsa_res_exp = gsa(growth_sim_exp, Morris(num_trajectory = 50), p_range, samples = 1_000_000); # ╔═╡ 1f9d7f84-fae7-4f0e-8783-7a8ce44051a3 +# missing begin plot(t_vals, gsa_res_exp.means[1][1, :], label = "μ", xlabel = "Time (day)", title = "Global sensitivities" @@ -557,7 +563,7 @@ Create a *reaction network object* for the Gompertz growth model. Name it `growt """ # ╔═╡ 1ae76036-42ef-46e9-88bf-d66d4267addc -# Uncomment and complete the instruction +# # Uncomment and complete the instruction # growth_gom = @reaction_network begin # @species missing # @parameters missing @@ -587,9 +593,11 @@ Initialize a vector `params_gom` with the parameter values: """ # ╔═╡ e9408c0f-cbcb-4703-8b85-e349e37b2025 +# μ_gom = missing μ_gom = 0.09; # ╔═╡ 24c4e348-7322-4c70-b8da-ec70d25356a9 +# D_gom = missing D_gom = 0.04; # ╔═╡ bf9291d5-6fc5-4359-874e-167237a22147 @@ -602,7 +610,7 @@ Create the ODE problem and store it in `oprob_gom`, using the same initial condi """ # ╔═╡ 3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 -# oprob_gom = missing # Uncomment and complete the instruction +# oprob_gom = missing # Uncomment and complete the instruction oprob_gom = ODEProblem(growth_gom, u0, tspan, params_gom); # ╔═╡ 30c35371-d106-4673-8569-d07b39edbdfe @@ -629,7 +637,7 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ be21e269-16c4-42b0-8930-2d009bd91161 -# Uncomment and complete the instruction +# # Uncomment and complete the instruction # function growth_sim_gom(params) # missing # ... From 785dc77b939e2df33fc404964700c2d51199319a Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Mon, 19 Oct 2026 17:13:21 +0200 Subject: [PATCH 4/8] suffix name of solved intro with `_sol` --- .../P7_sens_uncert/{sens_intro.jl => sens_intro_sol.jl} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename exercises/solved_notebooks/P7_sens_uncert/{sens_intro.jl => sens_intro_sol.jl} (100%) diff --git a/exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl b/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl similarity index 100% rename from exercises/solved_notebooks/P7_sens_uncert/sens_intro.jl rename to exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl From bab72cabac04048efcc31266487e8e8521cfe497 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Tue, 20 Oct 2026 10:15:46 +0200 Subject: [PATCH 5/8] remove standalone `...` from intro solution --- exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl b/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl index e11bcd3..c098a28 100644 --- a/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl +++ b/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl @@ -451,7 +451,6 @@ Write a solution function with as argument a vector of the parameters (that you # # Uncomment and complete the instruction # function growth_sim_exp(params) # missing -# ... # end function growth_sim_exp(params) μ, Wf = params @@ -640,7 +639,6 @@ Write a solution function with as argument a vector of the parameters (that you # # Uncomment and complete the instruction # function growth_sim_gom(params) # missing -# ... # end function growth_sim_gom(params) μ, D = params From 0281db6335144f028532b356261c7b7216d11858 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Tue, 20 Oct 2026 10:17:11 +0200 Subject: [PATCH 6/8] remove `# Uncomment and complete the instruction` --- .../P7_sens_uncert/sens_intro_sol.jl | 56 +++++++++---------- 1 file changed, 26 insertions(+), 30 deletions(-) diff --git a/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl b/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl index c098a28..26fea06 100644 --- a/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl +++ b/exercises/solved_notebooks/P7_sens_uncert/sens_intro_sol.jl @@ -371,7 +371,6 @@ Create a *reaction network object* for the exponential growth model. Name it `gr """ # ╔═╡ 2836f231-3b90-4829-9012-3ed9a09239ff -# # Uncomment and complete the instruction # growth_exp = @reaction_network begin # @species missing # @parameters missing @@ -392,7 +391,7 @@ Convert the system to a symbolic differential equation model (name it: `osys_exp """ # ╔═╡ 78d175b4-a9b2-49b5-bbc8-eb348799985b -# osys_exp = missing # Uncomment and complete the instruction +# osys_exp = missing osys_exp = convert(ODESystem, growth_exp) # ╔═╡ 05e96232-ce95-452b-98e4-e79817d45ae2 @@ -412,7 +411,7 @@ Initialize a vector `params_exp` with the parameter values: """ # ╔═╡ 7005d88f-9e28-47fe-9a32-dec023903dc3 -# params_exp = missing # Uncomment and complete the instruction +# params_exp = missing params_exp = [:μ => μ_exp, :Wf => Wf_exp]; # ╔═╡ c8853554-26c5-491c-a536-4c07e8c6a986 @@ -421,7 +420,7 @@ Create the ODE problem and store it in `oprob_exp`, using the same initial condi """ # ╔═╡ 025c1154-b2ae-4e1c-af1b-277b24d648a4 -# oprob_exp = missing # Uncomment and complete the instruction +# oprob_exp = missing oprob_exp = ODEProblem(growth_exp, u0, tspan, params_exp); # ╔═╡ f58825d2-55fc-43f1-b164-9555bf9f5b84 @@ -430,7 +429,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os """ # ╔═╡ 8044dd5f-9996-4062-a54a-8fb0eee30b46 -# osol_exp = missing # Uncomment and complete the instruction +# osol_exp = missing osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5); # ╔═╡ a1752ad4-2c11-416c-8321-9da058f9aaea @@ -439,7 +438,7 @@ Plot the result: """ # ╔═╡ 0f58fd3d-1495-46c7-9d92-5539395cf12e -# missing # Uncomment and complete the instruction +# missing plot(osol_exp) # ╔═╡ 6321774c-5bbe-42a7-bdad-0168c891b5ce @@ -448,7 +447,6 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ eb1e5f97-3e09-43f0-b6ce-cc6cbc42ec2f -# # Uncomment and complete the instruction # function growth_sim_exp(params) # missing # end @@ -465,7 +463,7 @@ Make a function based on the solution function that returns the output of a sing """ # ╔═╡ f2fa716f-d7b5-4c95-a8c6-6b16fbfe1499 -# growth_sim_W_exp(params) = missing # Uncomment and complete the instruction +# growth_sim_W_exp(params) = missing growth_sim_W_exp(params) = growth_sim_exp(params)[:W]; # ╔═╡ 98070d41-e8c5-49e6-9eea-04c718cbff65 @@ -474,7 +472,7 @@ Compute the output for the given parameter values. """ # ╔═╡ cef751aa-f9c8-46da-9436-9cc2e5d7515c -# W_exp = missing # Uncomment and complete the instruction +# W_exp = missing W_exp = growth_sim_W_exp([μ_exp, Wf_exp]); # ╔═╡ 9e7bb7ec-8fe3-422e-920f-13c2ef055feb @@ -483,7 +481,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s """ # ╔═╡ b13b5e77-ca60-465b-9f63-be9e5da0482f -# sens_W_exp = missing # Uncomment and complete the instruction +# sens_W_exp = missing sens_W_exp = ForwardDiff.jacobian(growth_sim_W_exp, [μ_exp, Wf_exp]); # ╔═╡ 735f5191-259e-440f-a379-20a2a70ec72c @@ -492,11 +490,11 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. """ # ╔═╡ 689556c9-a690-4e9e-9052-7ac66e999d4d -# sens_W_on_μ_exp = missing # Uncomment and complete the instruction +# sens_W_on_μ_exp = missing sens_W_on_μ_exp = sens_W_exp[:, 1]; # ╔═╡ b85f0d09-a102-435b-a898-149ba4b29caa -# sens_W_on_Wf_exp = missing # Uncomment and complete the instruction +# sens_W_on_Wf_exp = missing sens_W_on_Wf_exp = sens_W_exp[:, 2]; # ╔═╡ d067f814-1dd5-4e59-8cfd-3bc0b2d97612 @@ -505,11 +503,11 @@ Compute the normalized sensitivities. """ # ╔═╡ e01657f3-33ac-4eb8-b39e-e2441a165a8d -# sens_W_on_μ_rel_exp = missing # Uncomment and complete the instruction +# sens_W_on_μ_rel_exp = missing sens_W_on_μ_rel_exp = sens_W_on_μ_exp .* μ_exp ./ W_exp; # ╔═╡ 7231e54e-b293-48fc-b1f0-273f45f51539 -# sens_W_on_Wf_rel_exp = missing # Uncomment and complete the instruction +# sens_W_on_Wf_rel_exp = missing sens_W_on_Wf_rel_exp = sens_W_on_Wf_exp .* Wf_exp ./ W_exp; # ╔═╡ 0524a976-9d4d-4ef9-b088-bdd9102415e6 @@ -518,7 +516,7 @@ Plot both normalized sensitivity functions (with appropriate title and labels): """ # ╔═╡ 33031b53-6f2d-4256-a0aa-0ba97beaa7af -# missing # Uncomment and complete the instruction +# missing begin plot(t_vals, sens_W_on_μ_rel_exp, title = "Normalized sensitivities", label = "W on μ", xlabel = "Time (day)") @@ -562,7 +560,6 @@ Create a *reaction network object* for the Gompertz growth model. Name it `growt """ # ╔═╡ 1ae76036-42ef-46e9-88bf-d66d4267addc -# # Uncomment and complete the instruction # growth_gom = @reaction_network begin # @species missing # @parameters missing @@ -583,7 +580,7 @@ Convert the system to a symbolic differential equation model (name it: `osys_gom """ # ╔═╡ acaf5b91-accc-4dfa-9371-6495a47c8736 -# osys_gom = missing # Uncomment and complete the instruction +# osys_gom = missing osys_gom = convert(ODESystem, growth_gom) # ╔═╡ 70df05e7-4527-4087-a762-942a33e89f74 @@ -600,7 +597,7 @@ Initialize a vector `params_gom` with the parameter values: D_gom = 0.04; # ╔═╡ bf9291d5-6fc5-4359-874e-167237a22147 -# params_gom = missing # Uncomment and complete the instruction +# params_gom = missing params_gom = [:μ => μ_gom, :D => D_gom]; # ╔═╡ 9e3158c5-d4b2-456c-b920-8431e2047770 @@ -609,7 +606,7 @@ Create the ODE problem and store it in `oprob_gom`, using the same initial condi """ # ╔═╡ 3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 -# oprob_gom = missing # Uncomment and complete the instruction +# oprob_gom = missing oprob_gom = ODEProblem(growth_gom, u0, tspan, params_gom); # ╔═╡ 30c35371-d106-4673-8569-d07b39edbdfe @@ -618,7 +615,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os """ # ╔═╡ c01fe052-3701-4e22-9a64-f681d52a445a -# osol_gom = missing # Uncomment and complete the instruction +# osol_gom = missing osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5); # ╔═╡ 4f959d21-954b-4116-a498-bafc1aada47a @@ -627,7 +624,7 @@ Plot the result: """ # ╔═╡ 5f427df2-1b7e-46c7-adb9-66b05deac6e8 -# missing # Uncomment and complete the instruction +# missing plot(osol_gom) # ╔═╡ 958854fc-8e23-4364-88c5-1ec060b86a16 @@ -636,7 +633,6 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ be21e269-16c4-42b0-8930-2d009bd91161 -# # Uncomment and complete the instruction # function growth_sim_gom(params) # missing # end @@ -653,7 +649,7 @@ Make a function based on the solution function that returns a single output. """ # ╔═╡ 69f96b57-7c8d-4e64-9e6a-9804d001b0b1 -# growth_sim_W_gom(params) = missing # Uncomment and complete the instruction +# growth_sim_W_gom(params) = missing growth_sim_W_gom(params) = growth_sim_gom(params)[:W]; # ╔═╡ 09eeab39-79a1-420a-8298-72b7002ec168 @@ -662,7 +658,7 @@ Compute the output for the given parameter values. """ # ╔═╡ 9e15601f-4ad2-41be-8168-db30548a1c3b -# W_gom = missing # Uncomment and complete the instruction +# W_gom = missing W_gom = growth_sim_W_gom([μ_gom, D_gom]); # ╔═╡ 2b731e4a-d6eb-48f2-91f0-27e916d85683 @@ -671,7 +667,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s """ # ╔═╡ 68ad7e05-7f79-4e8d-9031-fefb5f5c0897 -# sens_W_gom = missing # Uncomment and complete the instruction +# sens_W_gom = missing sens_W_gom = ForwardDiff.jacobian(growth_sim_W_gom, [μ_gom, D_gom]); # ╔═╡ 362a8e24-8fac-4978-bc0e-13596a05d39e @@ -680,11 +676,11 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. """ # ╔═╡ e9114d64-fe3e-420d-9fde-4b0fcbbc8527 -# sens_W_on_μ_gom = missing # Uncomment and complete the instruction +# sens_W_on_μ_gom = missing sens_W_on_μ_gom = sens_W_gom[:,1]; # ╔═╡ 7c9e1112-0a63-417f-a65a-9c2bdfa8bc47 -# sens_W_on_D_gom = missing # Uncomment and complete the instruction +# sens_W_on_D_gom = missing sens_W_on_D_gom = sens_W_gom[:,2]; # ╔═╡ 94e66042-b352-4bf2-a24b-b15789e10fe3 @@ -693,11 +689,11 @@ Compute the normalized sensitivities. """ # ╔═╡ 1711346a-c161-415a-96f9-1235d786a584 -# sens_W_on_μ_rel_gom = missing # Uncomment and complete the instruction +# sens_W_on_μ_rel_gom = missing sens_W_on_μ_rel_gom = sens_W_on_μ_gom .* μ_gom ./ W_gom; # ╔═╡ 7a87bfbb-a8c0-4071-8f86-66513ae40968 -# sens_W_on_D_rel_gom = missing # Uncomment and complete the instruction +# sens_W_on_D_rel_gom = missing sens_W_on_D_rel_gom = sens_W_on_D_gom .* D_gom ./ W_gom; # ╔═╡ 7593a4bf-ad0f-4fd3-aa3c-8b6a5073cb02 @@ -706,7 +702,7 @@ Plot both sensitivity functions (with appropriate title and labels): " # ╔═╡ ad70a1d5-9b68-4cd6-9259-a3c25cad706b -# missing # Uncomment and complete the instruction +# missing plot(t_vals, [sens_W_on_μ_rel_gom, sens_W_on_D_rel_gom], title="Normalized sensitivities", label=["W on μ" "W on D"], xlabel="Time (day)") # ╔═╡ 12f64333-36d0-4e2a-9061-e3dcc7a4ae96 From 715e1d0927ebc6b3dafe2ea4adf7ad7247bf7372 Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Tue, 20 Oct 2026 10:38:30 +0200 Subject: [PATCH 7/8] update student version of sens intro --- .../P7_sens_uncert/sens_intro.jl | 472 ++++++++---------- 1 file changed, 204 insertions(+), 268 deletions(-) diff --git a/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl b/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl index 28c2b68..a3850e4 100644 --- a/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl +++ b/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl @@ -1,5 +1,5 @@ ### A Pluto.jl notebook ### -# v0.19.46 +# v1.0.1 using Markdown using InteractiveUtils @@ -7,8 +7,10 @@ using InteractiveUtils # ╔═╡ e1e7bc8e-7264-4cbc-98d2-aa73679fa2df begin # add this cell if you want the notebook to use the environment from where the Pluto server is launched - using Pkg - Pkg.activate("..") + using Pkg; Pkg.activate("../../solved_notebooks") + using Catalyst + using OrdinaryDiffEq, StatsPlots + using PlutoUI; TableOfContents() end # ╔═╡ 5f4fea06-0632-11ef-102e-21f5606d2056 @@ -17,15 +19,12 @@ using Markdown # ╔═╡ 489b5399-fe4c-481a-834f-0101bbe28cea using InteractiveUtils -# ╔═╡ 22489bd4-ab64-4bf9-ad03-5372ea273935 -using Catalyst - -# ╔═╡ cdcebbb1-40e0-457f-a6ec-b769f6b1f2e9 -using OrdinaryDiffEq, StatsPlots - # ╔═╡ 65885bbe-7b73-4efd-a8e7-984a34b80548 using ForwardDiff +# ╔═╡ a493c9e6-e250-4b7b-8cf2-03e3ac40cf19 +using GlobalSensitivity + # ╔═╡ 427b509f-d08e-4d93-99ba-a79f9c244b28 md""" # Introduction to sensitivity analysis @@ -33,9 +32,12 @@ md""" # ╔═╡ 49e7085c-3691-4182-9630-68dc9371ad18 md""" -## Goal of this practicum +## Goal of this practical """ +# ╔═╡ faf316bb-0ccc-4471-a53b-2476c3646705 +md"### Sensitivity functions" + # ╔═╡ c41894dd-0f6c-483f-b19f-dbf8148f776f md""" Sensitivity functions indicate how sensitive the model output is to a change in parameter values. When a model output is very sensitive to a certain parameter, a small change in the value of this parameter will have a large influence on the value of the model output. Sensitivity functions thus provide important information about the model and are implicitely used to estimate parameters and explicitely in the context of optimal experimental design. @@ -48,19 +50,36 @@ The sensitivity function that measures how sensitive output $y_j$ is to changes $$\cfrac{\partial \hat{y}_i(\theta)}{\partial \theta_j}\tag{1}$$ """ -# ╔═╡ abd70f04-1bf8-454f-869a-0d8082d68453 +# ╔═╡ ec1ee750-1f35-4f51-a8c7-d42a18a7b233 md""" Try to understand why the above expression $(1)$ does indeed give us the information we were promised in the first paragraph. How will we be able to see from the value determined by the above expression $(1)$ that whether or not output $y_j$ is sensitive to a change in $\theta_i$? """ -# ╔═╡ e76aedec-69f2-4301-b015-6960e4503c42 +# ╔═╡ f4ee83d1-259e-4dea-a86c-8c909ae4ec61 md""" -Answer: missing +!!! example + Consider the function $f(x) = a \, e^{b \, x}$, with parameters $\vec{θ} = (a, b)$. + + The sensitivity to $a$ is + + ```math + \dfrac{\partial f}{\partial a} = e^{b \, x} \, , + ``` + and for $b$: + + ```math + \dfrac{\partial f}{\partial b} = a \, x \, e^{b \, x} \, . + ``` + + If we now take a set of values for the input and parameters, such as $(x = 2, a = 1, b = 1)$, we can calculate the sensitivities at this point: + - For $a$: $\frac{\partial f}{\partial a} = e^2$ + - For $b$: $\frac{\partial f}{\partial b} = 2 e^2$ + from which we learn the function is twice as sensitive to $b$ compared to $a$ here. """ # ╔═╡ a6a97279-0b64-45c3-8312-c22b1a8425d0 md""" -Sometimes expression $(1)$ can be evaluated analytically. Usually, however, we will have to approximate the partial derivative numerically. Expression $(1)$ can be made more specific: +For simple functions, $(1)$ can be evaluated analytically. Usually, however, we will have to approximate the partial derivative numerically. Expression $(1)$ can be made more specific: $$\cfrac{\partial \hat{y}_i(\theta)}{\partial \theta_j} \approx \cfrac{\hat{y}_i(\theta_j+\Delta\theta_j)-\hat{y}_i(\theta_j)}{\Delta\theta_j}$$ """ @@ -70,6 +89,9 @@ md""" Thus, to calculate the sensitivity function numerically, the model is evaluated for the parameter values $\theta_i$ and $\theta_i+\Delta\theta_i$ and the difference between these evaluations is taken. """ +# ╔═╡ 33a71f15-dc39-45ed-b26c-855f42c55cae +md"### Normalized sensitivity" + # ╔═╡ 2a544bef-d7fe-4300-bc59-69e6f7429304 md""" Since quantity $(1)$ is dependent on the units, a normalized variant is often used: @@ -81,9 +103,14 @@ The interpretation of $(2)$ is how much the output changes per cent if the param Using the normalized variant allows you to compare all possible sensitivity functions with each other. """ -# ╔═╡ 7107a3a9-15ef-488e-8709-52a6444d0e1c +# ╔═╡ ec590df1-91be-40a2-844e-4ec3ea5b1077 +md"### Global sensitivity analysis" + +# ╔═╡ 3e2b5d93-4385-4261-9486-d93511ee7c26 md""" -We now calculate and interpret sensitivity functions for some given models. To illustrate the concepts, we first consider three simple models modelling the growth of grass. +The sensitivity function evaluates the partial derivative of your function **given some values for your parameters** (see also the given example). Therefore, the results of your sensitivity analysis will change based on the value of θ. As these methods are dependent on the current location of θ, they are called **local sensitivity analysis** methods. + +**Global sensitivity methods** try to circumvent this problem by quantifying sensitivity using **multiple** values of the output $y$. We consider two methods in this course: the Sobol method and the Morris method. For more information regarding these methods, we refer to the theory. """ # ╔═╡ 04e95855-7c4a-4d2c-b836-c5dde291adad @@ -107,9 +134,17 @@ with output $W$ the grass yield, and $W_f$, $\mu$ and $D$ parameters. The table | Exponential | 0.02 | 10.0 | | | Gompertz | 0.09 | | 0.04 | -We will use an initial condition of $W_0 = 2.0$ for each and a simulation time of $100$ days. +For all simulations (regardless of model), we will use an initial condition of $W_0 = 2.0$, a simulation time of $100$ days and save our solution every $0.5$ days. """ +# ╔═╡ 0b3d35bb-c5b0-44b7-94b3-06fa571d339e +begin + const u0 = [:W => 2.0] + const tspan = (0.0, 100.0) + const Δt = 0.5 + const t_vals = tspan[1]:Δt:tspan[2] # time values of the simulation +end; + # ╔═╡ 8d55bc42-1f8f-4bc1-b9fa-e8152fc125ce md""" We will illustrate how to compute the local sensitivity functions for the logistic model. The same will be left as exercises below for the exponential and gompertz models. @@ -145,32 +180,21 @@ Convert the *reaction model* to check that we work with the correct differential """ # ╔═╡ ba181db8-d176-4b83-9168-d5939ffe9661 -osys_log = convert(ODESystem, growth_log) - -# ╔═╡ cd8a7ba1-194a-4b42-8c87-2c9b1fe6b475 -md""" -Setting initial conditions, timespan and parameter values: -""" - -# ╔═╡ cbb2ca49-b019-495b-9310-83fcc00cad26 -u0_log = [:W => 2.0] - -# ╔═╡ 0b3d35bb-c5b0-44b7-94b3-06fa571d339e -tspan = (0.0, 100.0) # this will be the same for the three models +osys_log = convert(ODESystem, growth_log) # ╔═╡ be565a3c-31b6-4df1-b73b-08f308a8c09b md""" For the sake of clarity, we will use the variables `μ_log` and `Wf_log` to store the parameter values. """ -# ╔═╡ f62806d1-77e1-470b-9711-33a924c788cc -μ_log = 0.07 +# ╔═╡ dca01fd2-28f9-4d36-9b8f-b209ab97d3db +μ_log = 0.07; -# ╔═╡ 546ed163-26a7-4235-982f-7568ed609488 -Wf_log = 10.0 +# ╔═╡ 292aff06-7dfa-424b-8a90-173c37b9e67b +Wf_log = 10.0; # ╔═╡ 0da53fa2-5a42-46e6-8bd8-45d6aa903d46 -params_log = [:μ => μ_log, :Wf => Wf_log] +params_log = [:μ => μ_log, :Wf => Wf_log]; # ╔═╡ b1298f40-4696-49d0-ac94-896e0cdbc996 md""" @@ -178,12 +202,12 @@ Creating and solving the ODEProblem and plotting results: """ # ╔═╡ 270647d2-1371-4272-8bc1-3a6ad77bc716 -oprob_log = ODEProblem(growth_log, u0_log, tspan, params_log) +oprob_log = ODEProblem(growth_log, u0, tspan, params_log); # Also possible: # oprob_log = ODEProblem(growth_mod_log, [], tspan, []) # ╔═╡ ac235d86-1d93-4944-aa89-1b4fd38f0e6e -osol_log = solve(oprob_log, Tsit5(), saveat=0.5) +osol_log = solve(oprob_log, Tsit5(), saveat = Δt); # ╔═╡ 3e13efa1-9bc6-456f-8e62-ecd3165e2a65 plot(osol_log) @@ -206,30 +230,18 @@ We need to write a solution function with as argument a vector of the parameters # ╔═╡ bc201b61-0f58-49d2-a20b-4f18c42fcc96 function growth_sim_log(params) μ, Wf = params - u0_log = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_log = ODEProblem(growth_log, u0_log, tspan, [:μ=>μ, :Wf=>Wf]) - osol_log = solve(oprob_log, Tsit5(), saveat=0.5) + oprob_log = ODEProblem(growth_log, u0, tspan, [:μ=>μ, :Wf=>Wf]) + osol_log = solve(oprob_log, Tsit5(), saveat = Δt) return osol_log end # ╔═╡ 11f15a0e-2958-4075-928c-5a24bcc00c69 md""" -Next we will need to make a function based on the solution function that returns a single output. +Next we will need to make a function based on the solution function that returns the output of a single variable. """ # ╔═╡ d0f3197f-3094-4f6e-b33c-e5b74e0947f6 -growth_sim_W_log(params) = growth_sim_log(params)[:W] - -# ╔═╡ ba4eaa90-b97a-4b71-b3e1-968340c6def5 -md""" -Now make a time vector that is the same as the time vector from the solution. -""" - -# ╔═╡ 71f897a6-1e91-412b-ad58-c5f1e7cd1adb -t_vals_log = 0:0.5:100.0 -# Alternative: -# t_vals_log = tspan[1]:0.5:tspan[2] +growth_sim_W_log(params) = growth_sim_log(params)[:W]; # ╔═╡ b3b97b57-a4f7-4cb6-80c8-6ad893ade75d md""" @@ -237,7 +249,7 @@ Compute the single output with the given parameter values. """ # ╔═╡ 99864b21-3db0-4d60-a84e-8e96db4de4ae -W_log = growth_sim_W_log([μ_log, Wf_log]) +W_log = growth_sim_W_log([μ_log, Wf_log]); # ╔═╡ 344eb3a7-8614-4349-a536-62d9acef6bda md""" @@ -245,7 +257,7 @@ Use the function `ForwardDiff.jacobian` to compute the sensitivities. This funct """ # ╔═╡ 670ccfa4-24fa-4167-a46a-2a48dc19538b -sens_W_log = ForwardDiff.jacobian(growth_sim_W_log, [μ_log, Wf_log]) +sens_W_log = ForwardDiff.jacobian(growth_sim_W_log, [μ_log, Wf_log]); # ╔═╡ c262d3ed-9f62-45ad-a7ee-ec78ab16f35c md""" @@ -255,10 +267,10 @@ To get the sensitivities of $W$ on $\mu$, and of $W$ on $W_f$, you need to use i """ # ╔═╡ 6600d28e-2522-4069-a5e0-643be43f6117 -sens_W_on_μ_log = sens_W_log[:,1] # sensitivity of W on μ +sens_W_on_μ_log = sens_W_log[:,1]; # sensitivity of W on μ # ╔═╡ d9e3a4ac-c138-4b93-a19c-737166a0f0ea -sens_W_on_Wf_log = sens_W_log[:,2] # sensitivity of W on Wf +sens_W_on_Wf_log = sens_W_log[:,2]; # sensitivity of W on Wf # ╔═╡ f07d8205-e172-45fa-b955-792bd95f3023 md""" @@ -266,50 +278,81 @@ We now calculate the normalized sensitivities. For that we need to multiply by t """ # ╔═╡ 1c142797-ada6-4d77-a22f-b981ffd38956 -sens_W_on_μ_rel_log = sens_W_on_μ_log .* μ_log ./ W_log +sens_W_on_μ_rel_log = sens_W_on_μ_log .* μ_log ./ W_log; # ╔═╡ cb397e0b-56f9-420f-a86e-bfab35286b44 -sens_W_on_Wf_rel_log = sens_W_on_Wf_log .* Wf_log ./ W_log +sens_W_on_Wf_rel_log = sens_W_on_Wf_log .* Wf_log ./ W_log; # ╔═╡ 9e6946a4-f207-4bfb-9e19-aa7b77c2a05b md""" -We are now ready to plot the two sensitivity functions. We provide the time vector (first argument) and a vector of the two sensitivity functions (second argument). Additionally, you can provide a title, (legend) labels and a x- and/or y-label. +We can now plot the two sensitivity functions as follows: """ -# ╔═╡ f58ce914-f366-4087-a7d1-8cfe69ac623b -plot(t_vals_log, [sens_W_on_μ_rel_log, sens_W_on_Wf_rel_log], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") +# ╔═╡ dc0507ee-0f67-4556-bdc3-1294177c86b7 +begin + plot( + t_vals, sens_W_on_μ_rel_log, title = "Normalized sensitivities", + label = "W on μ", xlabel = "Time (day)" + ) + plot!(t_vals, sens_W_on_Wf_rel_log, label="W on Wf") +end -# ╔═╡ 05972c7f-f64f-4865-b0a7-f33029d0a6fa +# ╔═╡ 04a1d1ad-53a9-4146-be09-65aa422e3730 md""" -Notice that in the `label` option there is no comma separating the labels. +Conclusions: +- From the sensitivity plot of $W$ on $\mu$ it can be seen that $W$ is most sensitive to $\mu$ in the time region $[0, 30]$. The latter corresponds to the region where the yield rate is largest (i.e., when the growth is largest). This makes sense because when looking at the differential equation, $\mu$ is approximately the growth rate for relatively small $W$ values. +- From the sensitivity plot of $W$ on $W_f$ it can be seen that $W$ is most sensitive to $W_f$ in the region where time values are large (cf. operating point). The latter corresponds to the region where the yield rate stagnates (i.e., when the yield reaches a steady value). This makes sense because when looking at the differential equation, $W_f$ is the steady state value. """ -# ╔═╡ 88aa13f6-e16c-40cd-ac10-3ede8b9fb429 -# ╠═╡ disabled = true -#=╠═╡ +# ╔═╡ 5152abfc-98fa-4e10-aa28-6455b996a256 +md"### Global sensitivity analysis (GSA)" + +# ╔═╡ 8af97e64-6a14-4c7d-b935-002af339c94a md""" -#### Alternative +To perform global sensitivity analysis, we use the `gsa` function from the `GlobalSensitivity` package. This function requires 4 inputs: +1. Similarly to local sensitivity analysis, we require a function that takes the model parameters as its input and returns a model output of interest. We can just reuse `growth_sim_log` here to get the model outputs at all timepoints. +1. A GSA method. The methods considered in this course are the Sobol method `Sobol()` and the Morris method `Morris()`. We will use Sobol here. +1. A range for every parameter from where values will be sampled. For our grass growth models, we will assume realistic values for $μ$ are $[0, 0.2]$ and for $W_f$ are $[1, 20]$. +1. The amount of samples. As the Sobol method is computationally expensive, you typically don't want to use too many for an exercise. For the Morris method on the other hand, you typically want a large amount of samples for stable results (> `1_000_000`). +""" + +# ╔═╡ b86d3fe4-c3e5-4c39-8118-7970d609c17a +p_range = [[0.0, 0.2], [1.0, 20.0]]; -For each plot command we need to provide the time vector (first argument) and the column vector with the local sensitivity (second argument). Additionally, you can provide a title, a (legend) label and a x- and/or y-label. +# ╔═╡ 66083965-8402-47b6-9ccd-9ced9c7412e7 +gsa_res = gsa(growth_sim_log, Sobol(), p_range, samples = 100) -Beware that if you want to execute multiple commands (here: `plot`) in a single cell, you need to put them in a `begin`-and-`end` block. Also if you want to visualize subsequent graphs in the same plot, the forthcoming plot command names should be followed by a `!`, as in `plot!(...)`. +# ╔═╡ 97eed321-4dfe-468b-8b4e-6748608a781d +md""" +!!! tip + Does a function return some weird type you don't know what to do with? You can find out the fields of a variable `x` using `fieldnames(typeof(x))`. """ - ╠═╡ =# -# ╔═╡ dc0507ee-0f67-4556-bdc3-1294177c86b7 -# ╠═╡ disabled = true -#=╠═╡ +# ╔═╡ 8c923fea-2e75-42c1-ae73-0bdf8d1bee91 +fieldnames(typeof(gsa_res)) + +# ╔═╡ 3016e5f5-2f13-4877-9ed8-e6a355e63d5e +md""" +The resulting variable contains the first, second and total order Sobol indices, which we can access using respectively `gsa_res.S1`, `gsa_res.S2` and `gsa_res.ST`. We will plot the total effects here: +""" + +# ╔═╡ e2ec4060-342e-41c6-8bb2-1a20192ceeff begin -plot(t_vals_log, sens_W_on_μ_rel_log, title="Normalized sensitivities", label="W on μ", xlabel="Time (day)") -plot!(t_vals_log, sens_W_on_Wf_rel_log, label="W on Wf") + plot(t_vals, gsa_res.ST[1][1, :], label = "μ", xlabel = "Time (day)", + title = "Global sensitivities", ylims = (0, 1) + ) + plot!(t_vals, gsa_res.ST[2][1, :], label = "Wf") end - ╠═╡ =# -# ╔═╡ 04a1d1ad-53a9-4146-be09-65aa422e3730 +# ╔═╡ 73990cc1-77d5-4f09-ab2e-fb67c7ff81a4 md""" -Conclusions: -- From the sensitivity plot of $W$ on $\mu$ it can be seen that $W$ is most sensitive to $\mu$ in the time region $[0, 30]$. The latter corresponds to the region where the yield rate is largest (i.e., when the growth is largest). This makes sense because when looking at the differential equation, $\mu$ is approximately the growth rate for relatively small $W$ values. -- From the sensitivity plot of $W$ on $W_f$ it can be seen that $W$ is most sensitive to $W_f$ in the region where time values are large (cf. operating point). The latter corresponds to the region where the yield rate stagnates (i.e., when the yield reaches a steady value). This makes sense because when looking at the differential equation, $W_f$ is the steady state value. +!!! question + GSA returns an appreciably positive value for the effect of growth rate μ on final grass yield, while it was approximately $0$ in LSA. Why is this effect larger when we perform a global analysis where we consider multiple possible values for $μ$ and $W_f$? +""" + +# ╔═╡ 51488d93-48ad-4a02-9b51-00ef26f64699 +md""" +**Answer**: missing """ # ╔═╡ 79f0f1dd-850f-4dfd-b895-ff41a2d8adb8 @@ -328,19 +371,10 @@ Create a *reaction network object* for the exponential growth model. Name it `gr """ # ╔═╡ 2836f231-3b90-4829-9012-3ed9a09239ff -# Uncomment and complete the instruction -# growth_exp = @reaction_network begin -# @species missing -# @parameters missing -# missing -# end growth_exp = @reaction_network begin - @species W(t)=2.0 - @parameters μ=0.02 Wf=10.0 - (μ*Wf, μ), 0 <--> W - # Alternative: - # μ*Wf, 0 --> W - # μ, W --> 0 + @species missing + @parameters missing + missing end # ╔═╡ e8db84be-31a0-415d-a194-064c4c87a293 @@ -349,22 +383,7 @@ Convert the system to a symbolic differential equation model (name it: `osys_exp """ # ╔═╡ 78d175b4-a9b2-49b5-bbc8-eb348799985b -# osys_exp = missing # Uncomment and complete the instruction -osys_exp = convert(ODESystem, growth_exp) - -# ╔═╡ 63952c54-3304-4500-9536-b375c5c8f280 -md""" -Initialize a vector `u0_exp` with the initial condition: -""" - -# ╔═╡ 7dc013f7-7b35-42e4-aa77-0dca3755f389 -# u0_exp = missing # Uncomment and complete the instruction -u0_exp = [:W => 2.0] - -# ╔═╡ 86785c55-b567-449c-ae76-bc15a16223bc -md""" -We will use the same timespan as before, so no need to redefine it. -""" +osys_exp = missing # ╔═╡ 05e96232-ce95-452b-98e4-e79817d45ae2 md""" @@ -372,10 +391,10 @@ For the sake of clarity, we will use the variables `μ_exp` and `Wf_exp` to stor """ # ╔═╡ 9a705b3f-3c45-4896-9b89-e9c37d47fd89 -μ_exp = 0.02 +μ_exp = 0.02; # ╔═╡ 55a13ebc-09e7-433c-83a8-9b6805e2fbec -Wf_exp = 10.0 +Wf_exp = 10.0; # ╔═╡ f1293675-b4d0-4827-8714-65a593ecc00e md""" @@ -383,17 +402,15 @@ Initialize a vector `params_exp` with the parameter values: """ # ╔═╡ 7005d88f-9e28-47fe-9a32-dec023903dc3 -# params_exp = missing # Uncomment and complete the instruction -params_exp = [:μ => μ_exp, :Wf => Wf_exp] +params_exp = missing # ╔═╡ c8853554-26c5-491c-a536-4c07e8c6a986 md""" -Create the ODE problem and store it in `oprob_exp`: +Create the ODE problem and store it in `oprob_exp`, using the same initial condition and timespan as before. """ # ╔═╡ 025c1154-b2ae-4e1c-af1b-277b24d648a4 -# oprob_exp = missing # Uncomment and complete the instruction -oprob_exp = ODEProblem(growth_exp, u0_exp, tspan, params_exp) +oprob_exp = missing # ╔═╡ f58825d2-55fc-43f1-b164-9555bf9f5b84 md""" @@ -401,8 +418,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os """ # ╔═╡ 8044dd5f-9996-4062-a54a-8fb0eee30b46 -# osol_exp = missing # Uncomment and complete the instruction -osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5) +osol_exp = missing # ╔═╡ a1752ad4-2c11-416c-8321-9da058f9aaea md""" @@ -410,8 +426,7 @@ Plot the result: """ # ╔═╡ 0f58fd3d-1495-46c7-9d92-5539395cf12e -# missing # Uncomment and complete the instruction -plot(osol_exp) +missing # ╔═╡ 6321774c-5bbe-42a7-bdad-0168c891b5ce md""" @@ -419,37 +434,17 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ eb1e5f97-3e09-43f0-b6ce-cc6cbc42ec2f -# Uncomment and complete the instruction -# function growth_sim_exp(params) -# missing -# ... -# end function growth_sim_exp(params) - μ, Wf = params - u0_exp = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_exp = ODEProblem(growth_exp, u0_exp, tspan, [:μ=>μ, :Wf=>Wf]) - osol_exp = solve(oprob_exp, Tsit5(), saveat=0.5) - return osol_exp + missing end # ╔═╡ b31b4d24-ae30-4cff-b6b7-348200feaa5d md""" -Make a function based on the solution function that returns a single output. +Make a function based on the solution function that returns the output of a single variable. """ # ╔═╡ f2fa716f-d7b5-4c95-a8c6-6b16fbfe1499 -# growth_sim_W_exp(params) = missing # Uncomment and complete the instruction -growth_sim_W_exp(params) = growth_sim_exp(params)[:W] - -# ╔═╡ 0b413b98-a800-44e9-a917-ba61467bd613 -md""" -Make the time vector. -""" - -# ╔═╡ 6be6fea7-561b-4eec-b249-fa4522a5b039 -# t_vals = missing # Uncomment and complete the instruction -t_vals_exp = 0:0.5:100.0 +growth_sim_W_exp(params) = missing # ╔═╡ 98070d41-e8c5-49e6-9eea-04c718cbff65 md""" @@ -457,8 +452,7 @@ Compute the output for the given parameter values. """ # ╔═╡ cef751aa-f9c8-46da-9436-9cc2e5d7515c -# W_exp = missing # Uncomment and complete the instruction -W_exp = growth_sim_W_exp([μ_exp, Wf_exp]) +W_exp = missing # ╔═╡ 9e7bb7ec-8fe3-422e-920f-13c2ef055feb md""" @@ -466,8 +460,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s """ # ╔═╡ b13b5e77-ca60-465b-9f63-be9e5da0482f -# sens_W_exp = missing # Uncomment and complete the instruction -sens_W_exp = ForwardDiff.jacobian(growth_sim_W_exp, [μ_exp, Wf_exp]) +sens_W_exp = missing # ╔═╡ 735f5191-259e-440f-a379-20a2a70ec72c md""" @@ -475,12 +468,10 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. """ # ╔═╡ 689556c9-a690-4e9e-9052-7ac66e999d4d -# sens_W_on_μ_exp = missing # Uncomment and complete the instruction -sens_W_on_μ_exp = sens_W_exp[:,1] +sens_W_on_μ_exp = missing # ╔═╡ b85f0d09-a102-435b-a898-149ba4b29caa -# sens_W_on_Wf_exp = missing # Uncomment and complete the instruction -sens_W_on_Wf_exp = sens_W_exp[:,2] +sens_W_on_Wf_exp = missing # ╔═╡ d067f814-1dd5-4e59-8cfd-3bc0b2d97612 md""" @@ -488,12 +479,10 @@ Compute the normalized sensitivities. """ # ╔═╡ e01657f3-33ac-4eb8-b39e-e2441a165a8d -# sens_W_on_μ_rel_exp = missing # Uncomment and complete the instruction -sens_W_on_μ_rel_exp = sens_W_on_μ_exp .* μ_exp ./ W_exp +sens_W_on_μ_rel_exp = missing # ╔═╡ 7231e54e-b293-48fc-b1f0-273f45f51539 -# sens_W_on_Wf_rel_exp = missing # Uncomment and complete the instruction -sens_W_on_Wf_rel_exp = sens_W_on_Wf_exp .* Wf_exp ./ W_exp +sens_W_on_Wf_rel_exp = missing # ╔═╡ 0524a976-9d4d-4ef9-b088-bdd9102415e6 md""" @@ -501,8 +490,7 @@ Plot both normalized sensitivity functions (with appropriate title and labels): """ # ╔═╡ 33031b53-6f2d-4256-a0aa-0ba97beaa7af -# missing # Uncomment and complete the instruction -plot(t_vals_exp, [sens_W_on_μ_rel_exp, sens_W_on_Wf_rel_exp], title="Normalized sensitivities", label=["W on μ" "W on Wf"], xlabel="Time (day)") +missing # ╔═╡ 22b58c97-1aa2-49c4-8a0a-488c63014a90 md" @@ -512,6 +500,17 @@ Draw your conclusions: - missing " +# ╔═╡ ee17e112-a422-4b0e-be4c-c9e6abb040cb +md""" +Perform a global sensitivity analysis on this model, now using the **Morris method** and plot the results. You can use the same parameter range as before. +""" + +# ╔═╡ a8c930a2-a119-4681-a063-7a6c99615c5b +gsa_res_exp = missing + +# ╔═╡ 1f9d7f84-fae7-4f0e-8783-7a8ce44051a3 +missing + # ╔═╡ 41abf8c4-e67e-4f66-a38a-7204a878d98d md""" ### Exercise 2 - Sensitivity analysis of the Gompertz growth model @@ -523,19 +522,10 @@ Create a *reaction network object* for the Gompertz growth model. Name it `growt """ # ╔═╡ 1ae76036-42ef-46e9-88bf-d66d4267addc -# Uncomment and complete the instruction -# growth_gom = @reaction_network begin -# @species missing -# @parameters missing -# missing -# end growth_gom = @reaction_network begin - @species W(t)=2.0 - @parameters μ=0.09 D=0.04 - μ-D*log(W), W --> 2*W - # Alternative: - # μ, W --> 2*W - # D*log(W), W --> 0 + @species missing + @parameters missing + missing end # ╔═╡ e190ceca-30b1-49e2-baf5-ceb62929f4c0 @@ -544,51 +534,29 @@ Convert the system to a symbolic differential equation model (name it: `osys_gom """ # ╔═╡ acaf5b91-accc-4dfa-9371-6495a47c8736 -# osys_gom = missing # Uncomment and complete the instruction -osys_gom = convert(ODESystem, growth_gom) - -# ╔═╡ ed1150f1-9c75-4867-b7f5-535d605810f4 -md""" -Initialize a vector `u0_gom` with the initial condition: -""" - -# ╔═╡ d35ef155-61aa-4b19-872f-d6e621f96572 -# u0_gom = missing # Uncomment and complete the instruction -u0_gom = [:W => 2.0] - -# ╔═╡ 3f546d58-15ff-442d-a4e8-7c0783ed2fe0 -md""" -We will use the same timespan as before, so no need to redefine it. -""" - -# ╔═╡ b778f786-87bc-4103-a792-2ae5cf004d50 -md""" -For the sake of clarity, we will use the variables `μ_gom` and `D_gom` to store the parameter values. -""" - -# ╔═╡ edd76a4e-e541-4a83-923f-6f593d8fadad -μ_gom = 0.09 - -# ╔═╡ 05f1ae6c-1229-4924-9c41-9554d6ca4261 -D_gom = 0.04 +osys_gom = missing # ╔═╡ 70df05e7-4527-4087-a762-942a33e89f74 md""" Initialize a vector `params_gom` with the parameter values: """ +# ╔═╡ e9408c0f-cbcb-4703-8b85-e349e37b2025 +μ_gom = missing + +# ╔═╡ 24c4e348-7322-4c70-b8da-ec70d25356a9 +D_gom = missing + # ╔═╡ bf9291d5-6fc5-4359-874e-167237a22147 -# params_gom = missing # Uncomment and complete the instruction -params_gom = [:μ => μ_gom, :D => D_gom] +params_gom = missing # ╔═╡ 9e3158c5-d4b2-456c-b920-8431e2047770 md""" -Create the ODE problem and store it in `oprob_gom`: +Create the ODE problem and store it in `oprob_gom`, using the same initial conditions and timespan as before: """ # ╔═╡ 3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 -# oprob_gom = missing # Uncomment and complete the instruction -oprob_gom = ODEProblem(growth_gom, u0_gom, tspan, params_gom) +oprob_gom = missing # ╔═╡ 30c35371-d106-4673-8569-d07b39edbdfe md""" @@ -596,8 +564,7 @@ Solve the ODE problem. Use `Tsit5()` and `saveat=0.5`. Store the solution in `os """ # ╔═╡ c01fe052-3701-4e22-9a64-f681d52a445a -# osol_gom = missing # Uncomment and complete the instruction -osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5) +osol_gom = missing # ╔═╡ 4f959d21-954b-4116-a498-bafc1aada47a md""" @@ -605,8 +572,7 @@ Plot the result: """ # ╔═╡ 5f427df2-1b7e-46c7-adb9-66b05deac6e8 -# missing # Uncomment and complete the instruction -plot(osol_gom) +missing # ╔═╡ 958854fc-8e23-4364-88c5-1ec060b86a16 md""" @@ -614,18 +580,8 @@ Write a solution function with as argument a vector of the parameters (that you """ # ╔═╡ be21e269-16c4-42b0-8930-2d009bd91161 -# Uncomment and complete the instruction -# function growth_sim_gom(params) -# missing -# ... -# end function growth_sim_gom(params) - μ, D = params - u0_gom = [:W => 2.0] - tspan = (0.0, 100.0) - oprob_gom = ODEProblem(growth_gom, u0_gom, tspan, [:μ=>μ, :D=>D]) - osol_gom = solve(oprob_gom, Tsit5(), saveat=0.5) - return osol_gom + missing end # ╔═╡ e57f5916-b753-429b-8434-8a85a67bd1fb @@ -634,17 +590,7 @@ Make a function based on the solution function that returns a single output. """ # ╔═╡ 69f96b57-7c8d-4e64-9e6a-9804d001b0b1 -# growth_sim_W_gom(params) = missing # Uncomment and complete the instruction -growth_sim_W_gom(params) = growth_sim_gom(params)[:W] - -# ╔═╡ f7f6c56b-4910-4641-a462-eacfa4b4d034 -md""" -Make the time vector. -""" - -# ╔═╡ 84c1bb38-23c3-4f70-be67-7da2300a737b -# t_vals_gom = missing # Uncomment and complete the instruction -t_vals_gom = 0:0.5:100.0 +growth_sim_W_gom(params) = missing # ╔═╡ 09eeab39-79a1-420a-8298-72b7002ec168 md""" @@ -652,8 +598,7 @@ Compute the output for the given parameter values. """ # ╔═╡ 9e15601f-4ad2-41be-8168-db30548a1c3b -# W_gom = missing # Uncomment and complete the instruction -W_gom = growth_sim_W_gom([μ_gom, D_gom]) +W_gom = missing # ╔═╡ 2b731e4a-d6eb-48f2-91f0-27e916d85683 md""" @@ -661,8 +606,7 @@ Using `ForwardDiff.jacobian` to compute the sensitivities for the single ouput(s """ # ╔═╡ 68ad7e05-7f79-4e8d-9031-fefb5f5c0897 -# sens_W_gom = missing # Uncomment and complete the instruction -sens_W_gom = ForwardDiff.jacobian(growth_sim_W_gom, [μ_gom, D_gom]) +sens_W_gom = missing # ╔═╡ 362a8e24-8fac-4978-bc0e-13596a05d39e md""" @@ -670,12 +614,10 @@ Extract the (absolute) sensitivities of the outputs on the different parameters. """ # ╔═╡ e9114d64-fe3e-420d-9fde-4b0fcbbc8527 -# sens_W_on_μ_gom = missing # Uncomment and complete the instruction -sens_W_on_μ_gom = sens_W_gom[:,1] +sens_W_on_μ_gom = missing # ╔═╡ 7c9e1112-0a63-417f-a65a-9c2bdfa8bc47 -# sens_W_on_D_gom = missing # Uncomment and complete the instruction -sens_W_on_D_gom = sens_W_gom[:,2] +sens_W_on_D_gom = missing # ╔═╡ 94e66042-b352-4bf2-a24b-b15789e10fe3 md""" @@ -683,12 +625,10 @@ Compute the normalized sensitivities. """ # ╔═╡ 1711346a-c161-415a-96f9-1235d786a584 -# sens_W_on_μ_rel_gom = missing # Uncomment and complete the instruction -sens_W_on_μ_rel_gom = sens_W_on_μ_gom .* μ_gom ./ W_gom +sens_W_on_μ_rel_gom = missing # ╔═╡ 7a87bfbb-a8c0-4071-8f86-66513ae40968 -# sens_W_on_D_rel_gom = missing # Uncomment and complete the instruction -sens_W_on_D_rel_gom = sens_W_on_D_gom .* D_gom ./ W_gom +sens_W_on_D_rel_gom = missing # ╔═╡ 7593a4bf-ad0f-4fd3-aa3c-8b6a5073cb02 md" @@ -696,8 +636,7 @@ Plot both sensitivity functions (with appropriate title and labels): " # ╔═╡ ad70a1d5-9b68-4cd6-9259-a3c25cad706b -# missing # Uncomment and complete the instruction -plot(t_vals_gom, [sens_W_on_μ_rel_gom, sens_W_on_D_rel_gom], title="Normalized sensitivities", label=["W on μ" "W on D"], xlabel="Time (day)") +missing # ╔═╡ 12f64333-36d0-4e2a-9061-e3dcc7a4ae96 md" @@ -713,30 +652,29 @@ Draw your conclusions: # ╠═e1e7bc8e-7264-4cbc-98d2-aa73679fa2df # ╟─427b509f-d08e-4d93-99ba-a79f9c244b28 # ╟─49e7085c-3691-4182-9630-68dc9371ad18 +# ╟─faf316bb-0ccc-4471-a53b-2476c3646705 # ╟─c41894dd-0f6c-483f-b19f-dbf8148f776f # ╟─7b58f2f9-cb93-4a20-bf90-62a9006b57d6 -# ╟─abd70f04-1bf8-454f-869a-0d8082d68453 -# ╟─e76aedec-69f2-4301-b015-6960e4503c42 +# ╟─ec1ee750-1f35-4f51-a8c7-d42a18a7b233 +# ╟─f4ee83d1-259e-4dea-a86c-8c909ae4ec61 # ╟─a6a97279-0b64-45c3-8312-c22b1a8425d0 # ╟─25d3ef62-382c-455c-92bb-fadfc650c5a6 +# ╟─33a71f15-dc39-45ed-b26c-855f42c55cae # ╟─2a544bef-d7fe-4300-bc59-69e6f7429304 -# ╟─7107a3a9-15ef-488e-8709-52a6444d0e1c +# ╟─ec590df1-91be-40a2-844e-4ec3ea5b1077 +# ╟─3e2b5d93-4385-4261-9486-d93511ee7c26 # ╟─04e95855-7c4a-4d2c-b836-c5dde291adad # ╟─4673bfdf-8f4a-42bd-a026-21a66d800f2b +# ╠═0b3d35bb-c5b0-44b7-94b3-06fa571d339e # ╟─8d55bc42-1f8f-4bc1-b9fa-e8152fc125ce -# ╠═22489bd4-ab64-4bf9-ad03-5372ea273935 -# ╠═cdcebbb1-40e0-457f-a6ec-b769f6b1f2e9 # ╟─676119ce-f4fe-41f9-8121-b2c21f0dd28c # ╟─5492210b-6ee6-4baf-b37e-d21358cdeb60 # ╠═18eece66-46fa-4458-aaed-f4c8fa002c20 # ╟─8bb0b24b-6cce-49cd-a625-5f375b92d9b7 # ╠═ba181db8-d176-4b83-9168-d5939ffe9661 -# ╟─cd8a7ba1-194a-4b42-8c87-2c9b1fe6b475 -# ╠═cbb2ca49-b019-495b-9310-83fcc00cad26 -# ╠═0b3d35bb-c5b0-44b7-94b3-06fa571d339e # ╟─be565a3c-31b6-4df1-b73b-08f308a8c09b -# ╠═f62806d1-77e1-470b-9711-33a924c788cc -# ╠═546ed163-26a7-4235-982f-7568ed609488 +# ╠═dca01fd2-28f9-4d36-9b8f-b209ab97d3db +# ╠═292aff06-7dfa-424b-8a90-173c37b9e67b # ╠═0da53fa2-5a42-46e6-8bd8-45d6aa903d46 # ╟─b1298f40-4696-49d0-ac94-896e0cdbc996 # ╠═270647d2-1371-4272-8bc1-3a6ad77bc716 @@ -749,8 +687,6 @@ Draw your conclusions: # ╠═bc201b61-0f58-49d2-a20b-4f18c42fcc96 # ╟─11f15a0e-2958-4075-928c-5a24bcc00c69 # ╠═d0f3197f-3094-4f6e-b33c-e5b74e0947f6 -# ╟─ba4eaa90-b97a-4b71-b3e1-968340c6def5 -# ╠═71f897a6-1e91-412b-ad58-c5f1e7cd1adb # ╟─b3b97b57-a4f7-4cb6-80c8-6ad893ade75d # ╠═99864b21-3db0-4d60-a84e-8e96db4de4ae # ╟─344eb3a7-8614-4349-a536-62d9acef6bda @@ -762,20 +698,25 @@ Draw your conclusions: # ╠═1c142797-ada6-4d77-a22f-b981ffd38956 # ╠═cb397e0b-56f9-420f-a86e-bfab35286b44 # ╟─9e6946a4-f207-4bfb-9e19-aa7b77c2a05b -# ╠═f58ce914-f366-4087-a7d1-8cfe69ac623b -# ╟─05972c7f-f64f-4865-b0a7-f33029d0a6fa -# ╠═88aa13f6-e16c-40cd-ac10-3ede8b9fb429 # ╠═dc0507ee-0f67-4556-bdc3-1294177c86b7 # ╟─04a1d1ad-53a9-4146-be09-65aa422e3730 +# ╟─5152abfc-98fa-4e10-aa28-6455b996a256 +# ╠═a493c9e6-e250-4b7b-8cf2-03e3ac40cf19 +# ╟─8af97e64-6a14-4c7d-b935-002af339c94a +# ╠═b86d3fe4-c3e5-4c39-8118-7970d609c17a +# ╠═66083965-8402-47b6-9ccd-9ced9c7412e7 +# ╟─97eed321-4dfe-468b-8b4e-6748608a781d +# ╠═8c923fea-2e75-42c1-ae73-0bdf8d1bee91 +# ╟─3016e5f5-2f13-4877-9ed8-e6a355e63d5e +# ╠═e2ec4060-342e-41c6-8bb2-1a20192ceeff +# ╟─73990cc1-77d5-4f09-ab2e-fb67c7ff81a4 +# ╠═51488d93-48ad-4a02-9b51-00ef26f64699 # ╟─79f0f1dd-850f-4dfd-b895-ff41a2d8adb8 # ╟─fa2270d0-2548-409c-a23f-4369d8bce8ec # ╟─dca50b37-b06c-4efe-881e-cf966ebc8fd7 # ╠═2836f231-3b90-4829-9012-3ed9a09239ff # ╟─e8db84be-31a0-415d-a194-064c4c87a293 # ╠═78d175b4-a9b2-49b5-bbc8-eb348799985b -# ╟─63952c54-3304-4500-9536-b375c5c8f280 -# ╠═7dc013f7-7b35-42e4-aa77-0dca3755f389 -# ╟─86785c55-b567-449c-ae76-bc15a16223bc # ╟─05e96232-ce95-452b-98e4-e79817d45ae2 # ╠═9a705b3f-3c45-4896-9b89-e9c37d47fd89 # ╠═55a13ebc-09e7-433c-83a8-9b6805e2fbec @@ -791,8 +732,6 @@ Draw your conclusions: # ╠═eb1e5f97-3e09-43f0-b6ce-cc6cbc42ec2f # ╟─b31b4d24-ae30-4cff-b6b7-348200feaa5d # ╠═f2fa716f-d7b5-4c95-a8c6-6b16fbfe1499 -# ╟─0b413b98-a800-44e9-a917-ba61467bd613 -# ╠═6be6fea7-561b-4eec-b249-fa4522a5b039 # ╟─98070d41-e8c5-49e6-9eea-04c718cbff65 # ╠═cef751aa-f9c8-46da-9436-9cc2e5d7515c # ╟─9e7bb7ec-8fe3-422e-920f-13c2ef055feb @@ -806,18 +745,17 @@ Draw your conclusions: # ╟─0524a976-9d4d-4ef9-b088-bdd9102415e6 # ╠═33031b53-6f2d-4256-a0aa-0ba97beaa7af # ╟─22b58c97-1aa2-49c4-8a0a-488c63014a90 +# ╟─ee17e112-a422-4b0e-be4c-c9e6abb040cb +# ╠═a8c930a2-a119-4681-a063-7a6c99615c5b +# ╠═1f9d7f84-fae7-4f0e-8783-7a8ce44051a3 # ╟─41abf8c4-e67e-4f66-a38a-7204a878d98d # ╟─f311e943-c96b-4af3-8757-c662fd302a88 # ╠═1ae76036-42ef-46e9-88bf-d66d4267addc # ╟─e190ceca-30b1-49e2-baf5-ceb62929f4c0 # ╠═acaf5b91-accc-4dfa-9371-6495a47c8736 -# ╟─ed1150f1-9c75-4867-b7f5-535d605810f4 -# ╠═d35ef155-61aa-4b19-872f-d6e621f96572 -# ╟─3f546d58-15ff-442d-a4e8-7c0783ed2fe0 -# ╟─b778f786-87bc-4103-a792-2ae5cf004d50 -# ╠═edd76a4e-e541-4a83-923f-6f593d8fadad -# ╠═05f1ae6c-1229-4924-9c41-9554d6ca4261 # ╟─70df05e7-4527-4087-a762-942a33e89f74 +# ╠═e9408c0f-cbcb-4703-8b85-e349e37b2025 +# ╠═24c4e348-7322-4c70-b8da-ec70d25356a9 # ╠═bf9291d5-6fc5-4359-874e-167237a22147 # ╟─9e3158c5-d4b2-456c-b920-8431e2047770 # ╠═3aab9073-cf4d-4ee5-9e8a-e3db586d9f68 @@ -829,8 +767,6 @@ Draw your conclusions: # ╠═be21e269-16c4-42b0-8930-2d009bd91161 # ╟─e57f5916-b753-429b-8434-8a85a67bd1fb # ╠═69f96b57-7c8d-4e64-9e6a-9804d001b0b1 -# ╟─f7f6c56b-4910-4641-a462-eacfa4b4d034 -# ╠═84c1bb38-23c3-4f70-be67-7da2300a737b # ╟─09eeab39-79a1-420a-8298-72b7002ec168 # ╠═9e15601f-4ad2-41be-8168-db30548a1c3b # ╟─2b731e4a-d6eb-48f2-91f0-27e916d85683 From 17f900a9ca437b41178ffdff2b98395fe31e416f Mon Sep 17 00:00:00 2001 From: Bram Spanoghe <70662421+bspanoghe@users.noreply.github.com> Date: Tue, 20 Oct 2026 10:47:55 +0200 Subject: [PATCH 8/8] correct Pkg activation in student intro --- exercises/student_notebooks/P7_sens_uncert/sens_intro.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl b/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl index a3850e4..274334a 100644 --- a/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl +++ b/exercises/student_notebooks/P7_sens_uncert/sens_intro.jl @@ -7,7 +7,7 @@ using InteractiveUtils # ╔═╡ e1e7bc8e-7264-4cbc-98d2-aa73679fa2df begin # add this cell if you want the notebook to use the environment from where the Pluto server is launched - using Pkg; Pkg.activate("../../solved_notebooks") + using Pkg; Pkg.activate("..") using Catalyst using OrdinaryDiffEq, StatsPlots using PlutoUI; TableOfContents()