Skip to content

Commit d094948

Browse files
author
Alexander Ororbia
committed
clean up and integrated hodgkin-huxley mini lesson in neurocog tutorials
1 parent 389c7aa commit d094948

File tree

8 files changed

+252
-8
lines changed

8 files changed

+252
-8
lines changed

docs/source/ngclearn.components.neurons.graded.rst

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,14 @@ ngclearn.components.neurons.graded.rateCell module
3636
:undoc-members:
3737
:show-inheritance:
3838

39+
ngclearn.components.neurons.graded.rateCellOld module
40+
-----------------------------------------------------
41+
42+
.. automodule:: ngclearn.components.neurons.graded.rateCellOld
43+
:members:
44+
:undoc-members:
45+
:show-inheritance:
46+
3947
ngclearn.components.neurons.graded.rewardErrorCell module
4048
---------------------------------------------------------
4149

docs/source/ngclearn.components.synapses.hebbian.rst

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,14 @@ ngclearn.components.synapses.hebbian.hebbianSynapse module
3636
:undoc-members:
3737
:show-inheritance:
3838

39+
ngclearn.components.synapses.hebbian.hebbianSynapseOld module
40+
-------------------------------------------------------------
41+
42+
.. automodule:: ngclearn.components.synapses.hebbian.hebbianSynapseOld
43+
:members:
44+
:undoc-members:
45+
:show-inheritance:
46+
3947
ngclearn.components.synapses.hebbian.traceSTDPSynapse module
4048
------------------------------------------------------------
4149

docs/source/ngclearn.components.synapses.modulated.rst

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,14 @@ ngclearn.components.synapses.modulated.MSTDPETSynapse module
1212
:undoc-members:
1313
:show-inheritance:
1414

15+
ngclearn.components.synapses.modulated.REINFORCESynapse module
16+
--------------------------------------------------------------
17+
18+
.. automodule:: ngclearn.components.synapses.modulated.REINFORCESynapse
19+
:members:
20+
:undoc-members:
21+
:show-inheritance:
22+
1523
Module contents
1624
---------------
1725

docs/source/ngclearn.utils.rst

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,14 @@ ngclearn.utils.io\_utils module
3232
:undoc-members:
3333
:show-inheritance:
3434

35+
ngclearn.utils.jaxProcess module
36+
--------------------------------
37+
38+
.. automodule:: ngclearn.utils.jaxProcess
39+
:members:
40+
:undoc-members:
41+
:show-inheritance:
42+
3543
ngclearn.utils.metric\_utils module
3644
-----------------------------------
3745

Lines changed: 213 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,213 @@
1+
# Lecture 2E: The Hodgkin-Huxley Cell
2+
3+
In this tutorial, we will study/setup one of the most important biophysical
4+
neuronal models in computational neuroscience -- the Hodgkin-Huxley (H-H) spiking
5+
cell model.
6+
7+
## Using and Probing the H-H Cell
8+
9+
Go ahead and make a new folder for this study and create a Python script,
10+
i.e., `run_hhcell.py`, to write your code for this part of the tutorial.
11+
12+
Now let's set up the controller for this lesson's simulation and construct a
13+
single component system made up of an H-H cell.
14+
15+
16+
### Instantiating the H-H Neuronal Cell
17+
18+
Constructing/setting up a dynamical system made up of only a single
19+
H-H cell amounts to the following:
20+
21+
```python
22+
from jax import numpy as jnp, random, jit
23+
import numpy as np
24+
25+
from ngclearn.utils.model_utils import scanner
26+
from ngcsimlib.context import Context
27+
from ngcsimlib.compilers.process import Process
28+
## import model-specific mechanisms
29+
from ngclearn.components.neurons.spiking.hodgkinHuxleyCell import HodgkinHuxleyCell
30+
31+
## create seeding keys (JAX-style)
32+
dkey = random.PRNGKey(1234)
33+
dkey, *subkeys = random.split(dkey, 6)
34+
35+
T = 20000 ## number of simulation steps to run
36+
dt = 0.01 # ms ## compute integration time constant
37+
38+
## H-H cell hyperparameters
39+
v_Na=115. ## sodium reversal potential
40+
v_K=-35. ## potassium reversal potential
41+
v_L=10.6 ## leak reversal potential
42+
g_Na=100. ## sodium conductance (per unit area)
43+
g_K=5. ## potassium conductance (per unit area)
44+
g_L=0.3 ## leak conductance (per unit area)
45+
v_thr=4. ## membrane potential threshold for binary pulse production
46+
47+
## create simple system with only one AdEx
48+
with Context("Model") as model:
49+
cell = HodgkinHuxleyCell(
50+
name="z0", n_units=1, tau_v=1., resist_m=1., v_Na=v_Na, v_K=v_K, v_L=v_L,
51+
g_Na=g_Na, g_K=g_K, g_L=g_L, thr=v_thr, integration_type="euler", key=subkeys[0]
52+
)
53+
54+
## create and compile core simulation commands
55+
advance_process = (Process()
56+
>> cell.advance_state)
57+
model.wrap_and_add_command(jit(advance_process.pure), name="advance")
58+
59+
reset_process = (Process()
60+
>> cell.reset)
61+
model.wrap_and_add_command(jit(reset_process.pure), name="reset")
62+
63+
## set up non-compiled utility commands
64+
@Context.dynamicCommand
65+
def clamp(x):
66+
cell.j.set(x)
67+
```
68+
69+
Notably, the H-H model is a four-dimensional differential equation system, invented in 1952
70+
by Alan Hodgkin and Andrew Huxley to describe the ionic mechanisms that underwrite the
71+
initiation and propagation of action potentials within squid giant axons (the axons responsible for
72+
controlling the squid's water jet propulsion system in squid) <b>[1]</b>.
73+
Notably, the H-H cell is one of the more complex cells supported by ngc-learn since it
74+
models membrane potential `v` ($\mathbf{v}_t$) as well as three gates or channels. The three channels are
75+
essentially probability values:
76+
`n` ($\mathbf{n}_t$) for the probability of potassium channel subunit activation,
77+
`m` ($\mathbf{m}_t$) for the probability of sodium channel subunit activation, and
78+
`h` ($\mathbf{h}_t$) for the probability of sodium channel subunit inactivation.
79+
80+
neurons and muscle cells. It is a continuous-time dynamical system.
81+
82+
Formally, the core dynamics of the H-H cell can be written out as follows:
83+
84+
$$
85+
\tau_v \frac{\partial \mathbf{v}_t}{\partial t} &= \mathbf{j}_t - g_Na * \mathbf{m}^3_t * \mathbf{h}_t * (\mathbf{v}_t - v_Na) - g_K * \mathbf{n}^4_t * (\mathbf{v}_t - v_K) - g_L * (\mathbf{v}_t - v_L) \\
86+
\frac{\partial \mathbf{n}_t}{\partial t} &= alpha_n(\mathbf{v}_t) * (1 - \mathbf{n}_t) - beta_n(\mathbf{v}_t) * \mathbf{n}_t \\
87+
\frac{\partial \mathbf{m}_t}{\partial t} &= alpha_m(\mathbf{v}_t) * (1 - \mathbf{m}_t) - beta_m(\mathbf{v}_t) * \mathbf{m}_t \\
88+
\frac{\partial \mathbf{h}_t}{\partial t} &= alpha_h(\mathbf{v}_t) * (1 - \mathbf{h}_t) - beta_h(\mathbf{v}_t) * \mathbf{h}_t
89+
$$
90+
91+
where we observe that the above four-dimensional set of dynamics is composed of nonlinear ODEs. Notice that, in each gate or channel probability ODE, there are two generator functions (each of which is a function of the membrane potential $\mathbf{v}_t$) that produces the necessary dynamic coefficients at time $t$; $\alpha_x(\mathbf{v}_t)$ and $\beta_x(\mathbf{v}_t)$ produce different biopphysical weighting values depending on which channel $x = \{n, m, h\}$ they are related to.
92+
Note that, in ngc-learn's implementation of the H-H cell model, most of the core coefficients have been generally set according to Hodgkin and Huxley's 1952 work but can be configured by the experimenter to obtain different kinds of behavior/dynamics.
93+
94+
### Simulating the H-H Neuronal Cell
95+
96+
To see how the H-H cell works, we next write some code for visualizing how
97+
the node's membrane potential and core related gates/channels evolve with time
98+
(over a period of about `200` milliseconds). We will inject a square input pulse current
99+
into our H-H cell (specifically into its `j` compartment) and observe how the cell behaves in response.
100+
Specifically, we simulate the injection of this kind of current via the code below:
101+
102+
```python
103+
## create input stimulation feed
104+
stim = np.zeros((1, T))
105+
stim[0, 7000:13000] = 50 ## inject square input pulse
106+
x_seq = jnp.array(stim)
107+
108+
## now simulate the model
109+
time_span = []
110+
outs = []
111+
v = []
112+
n = []
113+
m = []
114+
h = []
115+
model.reset()
116+
for ts in range(x_seq.shape[1]):
117+
x_t = jnp.array([[x_seq[0, ts]]]) ## get data at time t
118+
model.clamp(x_t)
119+
model.run(t=ts * dt, dt=dt)
120+
outs.append(a.s.value)
121+
n.append(cell.n.value[0, 0])
122+
m.append(cell.m.value[0, 0])
123+
h.append(cell.h.value[0, 0])
124+
v.append(cell.v.value[0, 0])
125+
print(f"\r {ts} v = {cell.v.value}", end="")
126+
time_span.append(ts*dt)
127+
outs = jnp.concatenate(outs, axis=1)
128+
v = jnp.array(v)
129+
time_span = jnp.array(time_span)
130+
outs = jnp.squeeze(outs)
131+
```
132+
133+
and we can plot the dynamics of the neuron's voltage `v` and its three gate/channel
134+
variables, `h`, `m`, and `n`, with the following:
135+
136+
```python
137+
import matplotlib.pyplot as plt
138+
139+
plt.figure(figsize=(10, 8))
140+
## create subplot 1 -- show membrane voltage potential
141+
ax1 = plt.subplot(311)
142+
ax1.plot(time_span, v - 70., color='b')
143+
ax1.set_ylabel("Potential (mV)")
144+
ax1.set_title("Hodgkin-Huxley Spike Model")
145+
146+
## create subplot 2 -- show electrical input current
147+
ax2 = plt.subplot(312)
148+
ax2.plot(time_span, jnp.squeeze(x_seq), color='r')
149+
ax2.set_ylabel("Stimulation (µA/cm^2)")
150+
151+
## create subplot 3 -- show activation of gate/channel variables
152+
ax3 = plt.subplot(313, sharex=ax1)
153+
ax3.plot(time_span, h, label='h')
154+
ax3.plot(time_span, m, label='m')
155+
ax3.plot(time_span, n, label='n')
156+
ax3.set_ylabel("Activation (frac)")
157+
ax3.legend()
158+
159+
plt.tight_layout()
160+
plt.savefig("{0}".format("hh_plot.jpg"))
161+
plt.close()
162+
```
163+
164+
You should get a compound plot that depict the evolution of the H-H cell's voltage
165+
and channel/gate variables, i.e., saved as `hh_plot.jpg` locally to
166+
disk, like the one below:
167+
168+
```{eval-rst}
169+
.. table::
170+
:align: center
171+
172+
+--------------------------------------------------------+
173+
| .. image:: ../../images/tutorials/neurocog/hh_plot.jpg |
174+
| :scale: 75% |
175+
| :align: center |
176+
+--------------------------------------------------------+
177+
```
178+
179+
A useful note is that the H-H cell above used Euler integration to step through its
180+
dynamics (this is the default/base routine for all cell components in ngc-learn).
181+
However, one could configure the cell to use the midpoint method for integration
182+
by setting its argument `integration_type = rk2` or the Runge-Kutta fourth-order
183+
routine via `integration_type=rk4` for cases where, at the cost of increased
184+
compute time, more accurate dynamics are possible.
185+
186+
## Optional: Setting Up The Components with a JSON Configuration
187+
188+
While you are not required to create a JSON configuration file for ngc-learn,
189+
to get rid of the warning that ngc-learn will throw at the start of your
190+
program's execution (indicating that you do not have a configuration set up yet),
191+
all you need to do is create a sub-directory for your JSON configuration
192+
inside of your project code's directory, i.e., `json_files/modules.json`.
193+
Inside the JSON file, you would write the following:
194+
195+
```json
196+
[
197+
{"absolute_path": "ngclearn.components",
198+
"attributes": [
199+
{"name": "HodgkinHuxleyCell"}]
200+
},
201+
{"absolute_path": "ngcsimlib.operations",
202+
"attributes": [
203+
{"name": "overwrite"}]
204+
}
205+
]
206+
```
207+
208+
## References
209+
210+
<b>[1]</b> Hodgkin, Alan L., and Andrew F. Huxley. "A quantitative description
211+
of membrane current and its application to conduction and excitation in nerve."
212+
The Journal of physiology 117.4 (1952): 500.
213+

docs/tutorials/neurocog/index.rst

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,7 @@ work towards more advanced concepts.
4747
fitzhugh_nagumo_cell
4848
izhikevich_cell
4949
adex_cell
50+
hodgkin_huxley_cell
5051

5152
.. toctree::
5253
:maxdepth: 1

ngclearn/components/neurons/spiking/hodgkinHuxleyCell.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -166,7 +166,7 @@ def advance_state(
166166
_, _m = step_rk2(0., m, dx_dt, dt, (alpha_m_of_v, beta_m_of_v))
167167
_, _h = step_rk2(0., h, dx_dt, dt, (alpha_h_of_v, beta_h_of_v))
168168
elif intgFlag == 4: ## Runge-Kutta 4th order
169-
_, _v = step_rk4(0., v, dv_dt, dt, (_j, m + 0., n + 0., h + 0., tau_v, g_Na, g_K, g_L, v_Na, v_K, v_L))
169+
_, _v = step_rk4(0., v, dv_dt, dt, (_j, m + 0., n + 0., h + 0., tau_v, g_Na, g_K, g_L, v_Na, v_K, v_L))
170170
## next, integrate different channels
171171
_, _n = step_rk4(0., n, dx_dt, dt, (alpha_n_of_v, beta_n_of_v))
172172
_, _m = step_rk4(0., m, dx_dt, dt, (alpha_m_of_v, beta_m_of_v))

ngclearn/components/synapses/modulated/MSTDPETSynapse.py

Lines changed: 5 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -94,13 +94,11 @@ def __init__(
9494
def evolve(
9595
dt, w_bound, w_eps, preTrace_target, mu, Aplus, Aminus, tau_elg, elg_decay, preSpike, postSpike, preTrace,
9696
postTrace, weights, dWeights, eta, modulator, eligibility
97-
):
98-
'''
99-
dW_dt = TraceSTDPSynapse._compute_update( ## use Hebbian/STDP rule to obtain a non-modulated update
100-
dt, w_bound, preTrace_target, mu, Aplus, Aminus, preSpike, postSpike, preTrace, postTrace, weights
101-
)
102-
dWeights = dW_dt ## can think of this as eligibility at time t
103-
'''
97+
):
98+
# dW_dt = TraceSTDPSynapse._compute_update( ## use Hebbian/STDP rule to obtain a non-modulated update
99+
# dt, w_bound, preTrace_target, mu, Aplus, Aminus, preSpike, postSpike, preTrace, postTrace, weights
100+
# )
101+
# dWeights = dW_dt ## can think of this as eligibility at time t
104102

105103
if tau_elg > 0.: ## perform dynamics of M-STDP-ET
106104
eligibility = eligibility * jnp.exp(-dt / tau_elg) * elg_decay + dWeights/tau_elg

0 commit comments

Comments
 (0)