FitzHugh-Nagumo ODE
The FitzHugh-Nagumo neuron, the two-variable reduction of Hodgkin and Huxley's excitable membrane.
- Name:
fitzhugh_nagumo_ode - Field:
biology/neuroscience - Formalism: ordinary differential equations
- References: FitzHugh 1961; Nagumo, Arimoto, and Yoshizawa 1962; Trzaska 2023; Datseris and Parlitz 2022
Factory
@register_model(
field='biology/neuroscience',
tags=('neuron', 'excitability', 'relaxation oscillations'),
references=('FitzHugh 1961', 'Nagumo, Arimoto, and Yoshizawa 1962', 'Trzaska 2023', 'Datseris and Parlitz 2022'),
defaults={'until': 100.0, 'dt': 0.2},
)
def fitzhugh_nagumo_ode(
v: float = -1.2,
w: float = -0.6,
*,
a: float = 0.7,
b: float = 0.8,
epsilon: float = 0.08,
I: float = 0.5,
) -> System:
"""Return the FitzHugh-Nagumo neuron, the two-variable reduction of Hodgkin and Huxley's excitable membrane.
Without current the neuron rests on a stable fixed point and answers a large enough kick with one spike. A current
`I` between about 0.33 and 1.42 destabilizes the rest and the neuron fires periodically.
Args:
v:
The initial potential.
w:
The initial recovery.
a:
The offset of the recovery.
b:
The damping of the recovery.
epsilon:
The time scale of the recovery.
I:
The applied current.
Returns:
system:
The model, with the potential `neuron.v` and the recovery `neuron.w`.
Examples:
>>> from skcomplex.models.biology.neuroscience import fitzhugh_nagumo_ode
>>> fitzhugh_nagumo_ode(I=0.0).interactions[0].I
0.0
"""
kinetics = FitzHughNagumoInteraction('kinetics', target='neuron', a=a, b=b, epsilon=epsilon, I=I)
return System(elements=[FitzHughNagumoNeuronElement('neuron', v=v, w=w)], interactions=[kinetics])
Components
class FitzHughNagumoNeuronElement(Element):
"""FitzHugh-Nagumo neuron, with a fast potential `v` and a slow recovery `w`.
The neuron has no laws of its own. FitzHugh-Nagumo kinetics change it.
Args:
label:
The label of the neuron.
visual:
How the element is drawn, or `None` for `default_visual_`.
v:
The initial potential.
w:
The initial recovery.
Examples:
>>> from skcomplex.models.biology.neuroscience import FitzHughNagumoNeuronElement
>>> FitzHughNagumoNeuronElement('neuron').w
-0.6
"""
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
v: Annotated[float, StateSpec()] = -1.2,
w: Annotated[float, StateSpec()] = -0.6,
) -> None:
super().__init__(label, visual=visual)
self.v = v
self.w = w
class FitzHughNagumoInteraction(Interaction):
"""FitzHugh-Nagumo kinetics of a fast variable `v` and a slow recovery `w` of a target.
The fast variable follows `v' = v - v³/3 - w + I` and the recovery `w' = epsilon (v + a - b w)`. On a field the
kinetics act on every cell, and a diffusion of the fast variable makes the medium excitable.
Args:
label:
The label of the interaction.
target:
The path of the element or field whose states follow the kinetics.
fast:
The name of the fast state, such as a membrane potential.
slow:
The name of the slow recovery state.
a:
The offset of the recovery.
b:
The damping of the recovery.
epsilon:
The time scale of the recovery.
I:
The applied current.
Raises:
ValidationError:
A state name is not a name, or the fast and the slow state are the same.
Examples:
>>> from typing import Annotated
>>> from skcomplex.core.elements import Element
>>> from skcomplex.core.system import System
>>> from skcomplex.engine.runtime import Simulation
>>> from skcomplex.models.biology.common.excitable import FitzHughNagumoInteraction
>>> from skcomplex.specs import StateSpec
>>> class CellElement(Element):
... def __init__(
... self, label, visual=None, *, v: Annotated[float, StateSpec()] = -1.2,
... w: Annotated[float, StateSpec()] = -0.6,
... ):
... super().__init__(label, visual=visual)
... self.v = v
... self.w = w
>>> kinetics = FitzHughNagumoInteraction('kinetics', target='cell', I=0.5)
>>> system = System(elements=[CellElement('cell')], interactions=[kinetics])
>>> trajectory = Simulation(system, dt=0.1).check().prepare().run(until=100.0).trajectory_
>>> bool(trajectory.values('cell.v').max() > 1.5)
True
"""
def __init__(
self,
label: str,
*,
target: Annotated[str, RoleSpec()],
fast: str = 'v',
slow: str = 'w',
a: Annotated[float, ParameterSpec()] = 0.7,
b: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.8,
epsilon: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.08,
I: Annotated[float, ParameterSpec()] = 0.0,
) -> None:
super().__init__(label)
self.target = target
self.fast = fast
self.slow = slow
self.a = a
self.b = b
self.epsilon = epsilon
self.I = I
def _check(self) -> None:
"""Raise an error for state names that are not names or that coincide.
Raises:
ValidationError:
A state name is not a name, or the fast and the slow state are the same.
"""
super()._check()
names = (self.fast, self.slow)
if not all(isinstance(name, str) and name.isidentifier() for name in names) or self.fast == self.slow:
error_msg = f'{self.label!r} reads the states {names!r}. Give two different names of states.'
raise ValidationError(error_msg)
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the cubic fast kinetics and the slow linear recovery.
Args:
s:
The namespace of the symbols the laws can use, with the namespace of the target.
Returns:
laws:
The rate laws of the fast and the slow state.
"""
v, w = getattr(s.target, self.fast), getattr(s.target, self.slow)
return [RateLaw(v, v - v**3 / 3 - w + s.I), RateLaw(w, s.epsilon * (v + s.a - s.b * w))]
Equations
neuron.v: diff(v_neuron(t), t) == -1/3*v_neuron(t)^3 + I_kinetics + v_neuron(t) - w_neuron(t)
neuron.w: diff(w_neuron(t), t) == -(b_kinetics*w_neuron(t) - a_kinetics - v_neuron(t))*epsilon_kinetics