From 1d5511be5d8b874b6a2b4ab6219a42087c5f1511 Mon Sep 17 00:00:00 2001 From: jbrea Date: Wed, 6 Mar 2024 10:12:25 +0100 Subject: [PATCH 1/3] Update spiking_neural_systems.md It looks odd, when the membrane potential in an integrate-and-fire neuron overshoots the threshold. Btw: I find it convenient to use `ComponentArrays` for parameters, such that I can write e.g. `integrator.p.Vth` instead of `integrator.p[4]`, but I didn't apply this change here, as it would introduce a new dependency. --- docs/src/examples/spiking_neural_systems.md | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/docs/src/examples/spiking_neural_systems.md b/docs/src/examples/spiking_neural_systems.md index 50d8ea418..ab80f312b 100644 --- a/docs/src/examples/spiking_neural_systems.md +++ b/docs/src/examples/spiking_neural_systems.md @@ -56,30 +56,30 @@ They can make discontinuous changes to the model when certain conditions are met ```@example spikingneural function thr(u, t, integrator) - integrator.u > integrator.p[4] + u - integrator.p[4] end function reset!(integrator) integrator.u = integrator.p[2] end -threshold = DiscreteCallback(thr, reset!) -current_step = PresetTimeCallback([2, 15], integrator -> integrator.p[5] += 210.0) +threshold = ContinuousCallback(thr, reset!, nothing) +current_step = PresetTimeCallback([2, 15], integrator -> integrator.p[5] += 150.0) cb = CallbackSet(current_step, threshold) ``` -Our condition is `thr(u,t,integrator)` and the condition kicks in when `integrator.u > integrator.p[4]` where `p[4]` is our threshold parameter `Vth`. Our effect of the condition is `reset!(integrator)`. It sets `u` back to the equilibrium potential `p[2]`. We then wrap both the condition and the effect into a `DiscreteCallback` called threshold. There is one more callback called `PresetTimeCallback` that is particularly useful. This one increases the input `p[5]` at `t=2` and `t=15` by `210.0`. Both callbacks are then combined into a `CallbackSet`. We are almost done to simulate our system, we just need to put numbers on our initial voltage and parameters. +Our condition is `thr(u,t,integrator)` and the condition kicks in when `integrator.u > integrator.p[4]` where `p[4]` is our threshold parameter `Vth`. Our effect of the condition is `reset!(integrator)`. It sets `u` back to the equilibrium potential `p[2]`. We then wrap both the condition and the effect into a `DiscreteCallback` called threshold. There is one more callback called `PresetTimeCallback` that is particularly useful. This one increases the input `p[5]` at `t=2` and `t=15` by `150.0`. Both callbacks are then combined into a `CallbackSet`. We are almost done to simulate our system, we just need to put numbers on our initial voltage and parameters. ```@example spikingneural u0 = -75 tspan = (0.0, 40.0) # p = (gL, EL, C, Vth, I) -p = [10.0, -75.0, 5.0, -55.0, 0] +p = [5.0, -75.0, 50.0, -55.0, 0] prob = ODEProblem(lif, u0, tspan, p, callback = cb) ``` -Our initial voltage is `u0 = - 75`, which will be the same as our equilibrium potential, so we start at a stable point. Then we define the timespan we want to simulate. The timescale of the LIF as it is defined conforms roughly to milliseconds. Then we define our parameters as `p = [10.0, -75.0, 5.0, -55.0, 0]`. Remember that `gL, EL, C, Vth, I = p`. Finally, we wrap everything into a call to `ODEProblem`. Can't forget the `CallbackSet`. With that, our model is defined. Now we just need to solve it with a quick call to `solve`. +Our initial voltage is `u0 = - 75`, which will be the same as our equilibrium potential, so we start at a stable point. Then we define the timespan we want to simulate. The timescale of the LIF as it is defined conforms roughly to milliseconds. Then we define our parameters as `p = [5.0, -75.0, 50.0, -55.0, 0]`. Remember that `gL, EL, C, Vth, I = p`. Finally, we wrap everything into a call to `ODEProblem`. Can't forget the `CallbackSet`. With that, our model is defined. Now we just need to solve it with a quick call to `solve`. ```@example spikingneural sol = solve(prob) @@ -91,7 +91,7 @@ First of all the `solve` output tells us if solving the system generally worked. plot(sol) ``` -We see that the model is resting at `-75` while there is no input. At `t=2` the input increases by `210` and the model starts to spike. Spiking does not start immediately because the input first has to charge the membrane capacitance. Notice how once spiking starts, it rapidly becomes extremely regular. Increasing the input again at `t=15` increases firing as we would expect, but it is still extremely regular. This is one of the features of the LIF. The firing frequency is regular for constant input and a linear function of the input strength. There are ways to make LIF models less regular. For example, we could use certain noise types at the input. We could also simulate a large number of LIF models and connect them synaptically. Instead of going into those topics, we will move on to the Izhikevich model, which is known for its ability to generate a large variety of spiking dynamics during constant inputs. +We see that the model is resting at `-75` while there is no input. At `t=2` the input increases by `150` and the model starts to spike. Spiking does not start immediately because the input first has to charge the membrane capacitance. Notice how once spiking starts, it rapidly becomes extremely regular. Increasing the input again at `t=15` increases firing as we would expect, but it is still extremely regular. This is one of the features of the LIF. The firing frequency is regular for constant input and a linear function of the input strength. There are ways to make LIF models less regular. For example, we could use certain noise types at the input. We could also simulate a large number of LIF models and connect them synaptically. Instead of going into those topics, we will move on to the Izhikevich model, which is known for its ability to generate a large variety of spiking dynamics during constant inputs. ## The Izhikevich Model From 855db8cad2ab55f0f27dde049797c484ed77049f Mon Sep 17 00:00:00 2001 From: Johanni Brea Date: Wed, 6 Mar 2024 10:51:36 +0100 Subject: [PATCH 2/3] use ComponentArrays --- docs/src/examples/spiking_neural_systems.md | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/docs/src/examples/spiking_neural_systems.md b/docs/src/examples/spiking_neural_systems.md index ab80f312b..57a8e8c4f 100644 --- a/docs/src/examples/spiking_neural_systems.md +++ b/docs/src/examples/spiking_neural_systems.md @@ -20,6 +20,7 @@ The LIF model has five parameters, `gL, EL, C, Vth, I` and we define it in the ` ```@example spikingneural using DifferentialEquations +using ComponentArrays using Plots gr() @@ -56,30 +57,29 @@ They can make discontinuous changes to the model when certain conditions are met ```@example spikingneural function thr(u, t, integrator) - u - integrator.p[4] + u - integrator.p.Vth end function reset!(integrator) - integrator.u = integrator.p[2] + integrator.u = integrator.p.EL end threshold = ContinuousCallback(thr, reset!, nothing) -current_step = PresetTimeCallback([2, 15], integrator -> integrator.p[5] += 150.0) +current_step = PresetTimeCallback([2, 15], integrator -> integrator.p.I += 150.0) cb = CallbackSet(current_step, threshold) ``` -Our condition is `thr(u,t,integrator)` and the condition kicks in when `integrator.u > integrator.p[4]` where `p[4]` is our threshold parameter `Vth`. Our effect of the condition is `reset!(integrator)`. It sets `u` back to the equilibrium potential `p[2]`. We then wrap both the condition and the effect into a `DiscreteCallback` called threshold. There is one more callback called `PresetTimeCallback` that is particularly useful. This one increases the input `p[5]` at `t=2` and `t=15` by `150.0`. Both callbacks are then combined into a `CallbackSet`. We are almost done to simulate our system, we just need to put numbers on our initial voltage and parameters. +Our condition is `thr(u,t,integrator)` and the condition kicks in when `u > integrator.p.Vth`. Our effect of the condition is `reset!(integrator)`. It sets `u` back to the equilibrium potential `p.EL`. We then wrap both the condition and the effect into a `ContinuousCallback` called threshold. There is one more callback called `PresetTimeCallback` that is particularly useful. This one increases the input `p.I` at `t=2` and `t=15` by `150.0`. Both callbacks are then combined into a `CallbackSet`. We are almost done to simulate our system, we just need to put numbers on our initial voltage and parameters. ```@example spikingneural u0 = -75 tspan = (0.0, 40.0) -# p = (gL, EL, C, Vth, I) -p = [5.0, -75.0, 50.0, -55.0, 0] +p = ComponentArray(gL = 5.0, EL = -75.0, C = 50.0, Vth = -55.0, I = 0) prob = ODEProblem(lif, u0, tspan, p, callback = cb) ``` -Our initial voltage is `u0 = - 75`, which will be the same as our equilibrium potential, so we start at a stable point. Then we define the timespan we want to simulate. The timescale of the LIF as it is defined conforms roughly to milliseconds. Then we define our parameters as `p = [5.0, -75.0, 50.0, -55.0, 0]`. Remember that `gL, EL, C, Vth, I = p`. Finally, we wrap everything into a call to `ODEProblem`. Can't forget the `CallbackSet`. With that, our model is defined. Now we just need to solve it with a quick call to `solve`. +Our initial voltage is `u0 = - 75`, which will be the same as our equilibrium potential, so we start at a stable point. Then we define the timespan we want to simulate. The timescale of the LIF as it is defined conforms roughly to milliseconds. Then we define our parameters as `p = ComponentArray(gL = 5.0, EL = -75.0, C = 50.0, Vth = -55.0, I = 0)`. Finally, we wrap everything into a call to `ODEProblem`. Can't forget the `CallbackSet`. With that, our model is defined. Now we just need to solve it with a quick call to `solve`. ```@example spikingneural sol = solve(prob) From 09faf0e9972499881b079dbad76de9b72c2b2259 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 8 Mar 2024 06:57:45 -0500 Subject: [PATCH 3/3] Update Project.toml --- docs/Project.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/docs/Project.toml b/docs/Project.toml index 484fa0be8..a59e8191d 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -5,6 +5,7 @@ BVProblemLibrary = "ded0fc24-dfea-4565-b1d9-79c027d14d84" BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" BoundaryValueDiffEq = "764a87c0-6b3e-53db-9096-fe964310641d" CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b" +ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66" DAEProblemLibrary = "dfb8ca35-80a1-48ba-a605-84916a45b4f8" DASKR = "165a45c3-f624-5814-8e85-3bdf39a7becd" DDEProblemLibrary = "f42792ee-6ffc-4e2a-ae83-8ee2f22de800"