Update spiking_neural_systems.md

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"

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)
-    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 = 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.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 `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 `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 = [10.0, -75.0, 5.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 = [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 = 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)
@@ -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