Ex0: Integrate-and-Fire

Leaky integrate-and-fire cells — iafTauCell with spike/reset events

Model: `iafTauCell`

The simplest NeuroML2 example: a leaky integrate-and-fire neuron with a single state variable \(v\), governed by:

\[\frac{dv}{dt} = \frac{E_L - v}{\tau}\]

With spike event: if \(v > \text{thresh}\), then \(v \leftarrow \text{reset}\).

Parameters: leakReversal=-50 mV, thresh=-55 mV, reset=-70 mV, tau=30 ms.

Note: since leakReversal > thresh, this cell will never fire. The NeuroML example also includes iafTauRefCell, iafCell, and iafRefCell variants. The TVBO version matches the iafTauCell (tau-based, no refractory period).

1. Define in TVBO

from tvbo import SimulationExperiment

exp = SimulationExperiment.from_string("""
label: "NeuroML Ex0: Integrate-and-Fire (iafTau)"
dynamics:
  name: IntegrateAndFire
  parameters:
    leakReversal: { value: -50.0 }
    tau:          { value: 30.0 }
    thresh:       { value: -55.0 }
    reset:        { value: -70.0 }
  state_variables:
    v:
      equation: { rhs: "(leakReversal - v) / tau" }
      initial_value: -50.0
      variable_of_interest: true
  events:
    spike:
      condition: { rhs: "v > thresh" }
      affect:    { rhs: "v = reset" }
network:
  number_of_nodes: 1
integration:
  method: euler
  step_size: 0.005
  duration: 300.0
  time_scale: ms
""")
print(f"Model: {exp.dynamics.name}")
print(f"State variables: {list(exp.dynamics.state_variables.keys())}")
print(f"Events: {list(exp.dynamics.events.keys())}")

Model: IntegrateAndFire
State variables: ['v']
Events: ['spike']

2. Render LEMS XML

xml = exp.render("lems")
print(xml[:2000])


<Lems>

  <!-- Tell jLEMS/jNeuroML which component is the simulation entry point. -->
  <Target component="sim_NeuroML_Ex0__Integrate_and_Fire__iafTau_"/>

  <!-- ════════════════════════════════════════════════════════════════
       Dimensions & Units (inline — no external includes needed)
       ════════════════════════════════════════════════════════════════ -->

  <!-- Dimensions -->
  <Dimension name="none"/>
  <Dimension name="time" t="1"/>
  <Dimension name="voltage" m="1" l="2" t="-3" i="-1"/>
  <Dimension name="per_time" t="-1"/>
  <Dimension name="conductance" m="-1" l="-2" t="3" i="2"/>
  <Dimension name="capacitance" m="-1" l="-2" t="4" i="2"/>
  <Dimension name="current" i="1"/>
  <Dimension name="resistance" m="1" l="2" t="-3" i="-2"/>
  <Dimension name="concentration" l="-3" n="1"/>
  <Dimension name="substance" n="1"/>
  <Dimension name="charge" t="1" i="1"/>
  <Dimension name="temperature" k="1"/>

  <!-- Units -->
  <Unit symbol="s" dimension="time" power="0"/>
  <Unit symbol="ms" dimension="time" power="-3"/>
  <Unit symbol="us" dimension="time" power="-6"/>
  <Unit symbol="V" dimension="voltage" power="0"/>
  <Unit symbol="mV" dimension="voltage" power="-3"/>
  <Unit symbol="A" dimension="current" power="0"/>
  <Unit symbol="mA" dimension="current" power="-3"/>
  <Unit symbol="nA" dimension="current" power="-9"/>
  <Unit symbol="pA" dimension="current" power="-12"/>
  <Unit symbol="S" dimension="conductance" power="0"/>
  <Unit symbol="mS" dimension="conductance" power="-3"/>
  <Unit symbol="nS" dimension="conductance" power="-9"/>
  <Unit symbol="F" dimension="capacitance" power="0"/>
  <Unit symbol="uF" dimension="capacitance" power="-6"/>
  <Unit symbol="nF" dimension="capacitance" power="-9"/>
  <Unit symbol="ohm" dimension="resistance" power="0"/>
  <Unit symbol="Mohm" dimension="resistance" power="6"/>
  <Unit symbol="per_s" dimension="per_time" power="0"/>
  <Unit symbol="per_ms" dimension="per_time" power="3"/>
  <Unit symbol="degC" dim

3. Run Reference (NeuroML2 via jNeuroML)

import sys, os
sys.path.insert(0, os.path.dirname(os.path.abspath(".")))
from _nml_helpers import run_lems_example, LEMS_EXAMPLES

ref_outputs = run_lems_example("LEMS_NML2_Ex0_IaF.xml")
print(f"Reference output files: {list(ref_outputs.keys())}")
for name, arr in ref_outputs.items():
    print(f"  {name}: shape={arr.shape}")

Reference output files: ['iaf_v.dat']
  iaf_v.dat: shape=(60001, 5)

4. Run TVBO Version

import numpy as np

result = exp.run("neuroml")
da = result.integration.data
time = da.coords['time'].values
sv_data = da.values
tvbo_arr = np.column_stack([time, sv_data])
print(f"TVBO output shape: {tvbo_arr.shape}")
print(f"Time range: [{tvbo_arr[0,0]:.3f}, {tvbo_arr[-1,0]:.3f}] s")

TVBO output shape: (60001, 2)
Time range: [0.000, 0.300] s

5. Numerical Comparison

The NeuroML iafTauCell with leakReversal=-50 and thresh=-55 never spikes — the membrane just decays toward the leak reversal. We compare the v traces:

from _nml_helpers import compare_traces
import numpy as np

# Reference: Ex0 has multiple cell types in one file; pick the iafTau columns
# The output has columns: time, iafTau/v, iafTauRef/v, iaf/v, iafRef/v
ref_arr = list(ref_outputs.values())[0]
ref_cols = ['time', 'iafTau_v', 'iafTauRef_v', 'iaf_v', 'iafRef_v']
tvbo_cols = ['time', 'v']

# Compare just the iafTau column (col 1)
ref_sub = ref_arr[:, [0, 1]]
compare_traces(
    ref_sub, tvbo_arr,
    ref_cols=['time', 'v'], tvbo_cols=['time', 'v'],
)

  v: RMSE=61.106006  max_err=69.930000  corr=1.000000  ⚠️

{'v': {'rmse': np.float64(61.10600646068382),
  'max_err': np.float64(69.93),
  'corr': np.float64(0.9999999999998032),
  'close': False}}

6. Plot

from _nml_helpers import plot_comparison

plot_comparison(
    ref_sub, tvbo_arr,
    ref_cols=['time', 'v'], tvbo_cols=['time', 'v'],
    title="Ex0: Integrate-and-Fire (iafTauCell) — NeuroML vs TVBO",
    time_scale=1.0, time_unit="s",
)

Model: iafTauCell