From abbf6bbdf6df8d2a3cadd0b3ffe171c685d53b59 Mon Sep 17 00:00:00 2001
From: boeschf <48126478+boeschf@users.noreply.github.com>
Date: Mon, 12 Dec 2022 13:20:12 +0100
Subject: [PATCH] Sde Tutorial (#2044)
Implementation of stochastic calcium plasticity curves, originally
authored by Sebastian Schmitt @schmitts, and adapted here for Arbor
inclusion.
- python example plasticity_stdp.py
- stochastic mechanism calcium_based_synapse.mod
- tutorial
Closes #1987
---
doc/dev/index.rst | 2 +-
doc/dev/numerics.rst | 2 +
doc/tutorial/calcium_stdp.svg | 1204 +++++++++++++++++
doc/tutorial/calcium_stdp_curve.rst | 192 +++
doc/tutorial/index.rst | 8 +
mechanisms/CMakeLists.txt | 2 +-
.../stochastic/calcium_based_synapse.mod | 99 ++
python/example/calcium_stdp.py | 234 ++++
scripts/run_python_examples.sh | 1 +
9 files changed, 1742 insertions(+), 2 deletions(-)
create mode 100644 doc/tutorial/calcium_stdp.svg
create mode 100644 doc/tutorial/calcium_stdp_curve.rst
create mode 100644 mechanisms/stochastic/calcium_based_synapse.mod
create mode 100644 python/example/calcium_stdp.py
diff --git a/doc/dev/index.rst b/doc/dev/index.rst
index 17d4db34b9..25fbb8f140 100644
--- a/doc/dev/index.rst
+++ b/doc/dev/index.rst
@@ -27,4 +27,4 @@ Here we document internal components of Arbor. These pages can be useful if you'
mechanism_abi
util
version
-.. numerics
+ numerics
diff --git a/doc/dev/numerics.rst b/doc/dev/numerics.rst
index 1523aadda6..c4cd5bd3bd 100644
--- a/doc/dev/numerics.rst
+++ b/doc/dev/numerics.rst
@@ -18,3 +18,5 @@ Mechanisms
Exponential Euler `cnexp`.
Implicit Euler `sparse`.
+
+Euler-Maruyama (explicit Euler) `stochastic`
diff --git a/doc/tutorial/calcium_stdp.svg b/doc/tutorial/calcium_stdp.svg
new file mode 100644
index 0000000000..d9bea02d3e
--- /dev/null
+++ b/doc/tutorial/calcium_stdp.svg
@@ -0,0 +1,1204 @@
+
+
+
diff --git a/doc/tutorial/calcium_stdp_curve.rst b/doc/tutorial/calcium_stdp_curve.rst
new file mode 100644
index 0000000000..899925cc87
--- /dev/null
+++ b/doc/tutorial/calcium_stdp_curve.rst
@@ -0,0 +1,192 @@
+.. _tutorial_calcium_stpd_curve:
+
+Spike Timing-dependent Plasticity Curve
+=======================================
+
+This tutorial uses a single cell and reproduces `this Brian2 example
+`_. We aim
+to reproduce a spike timing-dependent plastivity curve which arises from stochastic calcium-based
+synapse dynamics described in Graupner and Brunel [1]_.
+
+The synapse is modeled as synaptic efficacy variable, :math:`\rho`, which is a function of the
+calcium concentration, :math:`c(t)`. There are two stable states at :math:`\rho=0` (DOWN) and
+:math:`\rho=1` (UP), while :math:`\rho=\rho^\ast=0.5` represents a third unstable state between the
+two stable states. The calcium concentration dynamics are represented by a simplified model which
+uses a linear sum of individual calcium transients elicited by trains of pre- and postsynaptic
+action potentials:
+
+.. math::
+
+ \begin{align*}
+ c^\prime (t) &= - \frac{1}{\tau_{Ca}}c
+ + C_{pre} \sum_i \delta(t-t_i-D)
+ + C_{post} \sum_j \delta(t-t_j), \\
+ \rho^\prime(t) &= - \frac{1}{\tau}\left [
+ \rho (1 - \rho) (\rho^\ast - \rho)
+ -\gamma_p (1-\rho) H\left(c - \theta_p \right)
+ + \gamma_d \rho H\left(c - \theta_d \right) \right ]
+ + N, \\
+ N &= \frac{\sigma}{\sqrt{\tau}} \sqrt{H\left( c - \theta_p \right)
+ + H\left( c - \theta_d \right)} W.
+ \end{align*}
+
+Here, the sums over :math:`i` and :math:`j` represent the contributions from all pre and
+postsynaptic spikes, respectively, with :math:`C_{pre}` and :math:`C_{pre}` denoting the jumps in
+concentration after a spike. The jump after the presynaptic spike is delayed by :math:`D`. The
+calcium decay time is assumed to be much faster than the synaptic time scale,
+:math:`\tau_{Ca} \ll \tau`. The subscripts :math:`p` and :math:`d` represent potentiation (increase
+in synaptic efficacy) and depression (decrease in synaptic efficacy), respectively, with
+:math:`\gamma` and :math:`\theta` being the corresponding rates and thresholds. :math:`H(x)` is the
+right-continuous heaviside step function (:math:`H(0)=1`).
+
+This mechanism is stochastic, :math:`W` represents a white noise process, and therefore our
+simulation needs to
+
+- use a stochastic synapse mechanism,
+- accumulate statistics over a large enough ensemble of initial states.
+
+
+Implementation of a Stochastic Mechanism
+----------------------------------------
+
+Implementing a stochastic mechanism which is given by a stochastic differential equation (SDE) as
+above is straightforward to implement in :ref:`Arbor's NMODL dialect `. Let's examine
+the mechanism code in the `Arbor repository
+`_.
+
+The main difference compared to a deterministic (ODE) description is the additional `WHITE_NOISE`
+block,
+
+.. code:: none
+
+ WHITE_NOISE {
+ W
+ }
+
+which declares the white noise process :math:`W`, and the specification of the `stochastic` solver
+method,
+
+.. code:: none
+
+ BREAKPOINT {
+ SOLVE state METHOD stochastic
+ }
+
+This is sufficient to inform Arbor about the stochasticity of the mechanism. For more information
+about Arbor's strategy to solve SDEs, please consult :ref:`this overview `, while
+details about the numerical solver can be found in the :ref:`developers guide `.
+
+
+The Model
+---------
+
+In this tutorial, the neuron model itself is simple with only
+passive (leaky) membrane dynamics, and it receives regular synaptic current
+input in one arbitrary chosen control volume (CV) to trigger regular spikes.
+
+First we import some required modules:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 13-18
+
+Next we set the simulation parameters in order to reproduce the plasticity curve:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 20-41
+
+The time lag resolution, together with the maximum time lag, determine the number of cases we want
+to simulate. For each such case, however, we need to run many simulations in order to get a
+statistically meaningful result. The number of simulations per case is given by the ensemble size
+and the initial conditions. In our case, we have two inital states, :math:`\rho(0)=0` and
+:math:`\rho(0)=1`, and for each initial state we want to run :math:`100` simulations. We note, that
+the stochastic synapse mechanism does not alter the state of the cell, but couples one-way only by
+reacting to spikes. Therefore, we are allowed to simply place :math:`100` synapses per initial state
+onto the cell without worrying about interference. Moreover, this has the benefit of exposing
+parallelism that Arbor can take advantage of.
+
+Thus, we create a simple cell with a midpoint at which we place our mechanisms:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 44-67
+
+Since our stochastic mechanism `calcium_based_synapse` is not within Arbor's default set of
+mechanism, we need to extend the mechanism catalogue within the cable cell properties:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 70-74
+
+Our cell and cell properties can then later be used to create a simple recipe:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 77-103
+
+Note, that the recipe takes a cell, cell properties and a list of event generators as constructor
+arguments and returns them with its corresponding methods. Furthermore, the recipe also returns a
+list of probes which contains only one item: A query for our mechanism's state variable
+:math:`\rho`. Since we placed a number of these mechanisms on our cell, we will receive a vector of
+values when probing.
+
+Next we set up the simulation logic:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 106-134
+
+The pre- and postsynaptic events are generated at explicit schedules, where the presynaptic event
+is shifted in time by :math:`D -\text{time lag}` with respect to the presynaptic event, which in
+turn is generated regularly with the frequency :math:`f`. The postsynaptic events are driven by the
+deterministic synapse with weight `1.0`, while the presynaptic events are generated at the
+stochastic calcium synapses. The postsynaptic weight can be set arbitrarily as long as it is large
+enough to trigger the spikes.
+
+Thus, we have all ingredients to create the recipe
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 136-137
+
+Now, we need to initialize the simulation, register a probe and run the simulation:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 139-154
+
+Since we are interested in the long-term average value, we only query the probe at the end of the
+simulation.
+
+After the simulation is finished, we calculate the change in synaptic strength by evaluating the
+transition probabilies from initial DOWN state to final UP state and vice versa.
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 156-174
+
+Since we need to run our simulation for each time lag case anew, we spawn a bunch of threads to
+carry out the work in parallel:
+
+.. literalinclude:: ../../python/example/calcium_stdp.py
+ :language: python
+ :lines: 177-178
+
+The collected results can then be plotted:
+
+.. figure:: calcium_stdp.svg
+ :width: 1600
+ :align: center
+
+ Comparison of this simulation with reference simulation [1]_; for a simulation duration
+ of 60 spikes at 1 Hertz, ensemble size of 2000 per initial state and time step dt=0.01 ms.
+ The shaded region indicates the 95\% confidence interval.
+
+The full code
+-------------
+You can find the full code of the example at ``python/examples/calcium_stdp.py``.
+
+References
+----------
+.. [1] Graupner and Brunel, PNAS 109 (10): 3991-3996 (2012); ``_, ``_.
diff --git a/doc/tutorial/index.rst b/doc/tutorial/index.rst
index 629e81c95e..9bbcaeac5b 100644
--- a/doc/tutorial/index.rst
+++ b/doc/tutorial/index.rst
@@ -52,6 +52,14 @@ Probes
probe_lfpykit
+Stochastic Mechanisms
+---------------------
+
+.. toctree::
+ :maxdepth: 1
+
+ calcium_stdp_curve
+
Hardware
--------
diff --git a/mechanisms/CMakeLists.txt b/mechanisms/CMakeLists.txt
index 320fd85236..de0718e52c 100644
--- a/mechanisms/CMakeLists.txt
+++ b/mechanisms/CMakeLists.txt
@@ -23,7 +23,7 @@ make_catalogue(
make_catalogue(
NAME stochastic
- MOD ou_input
+ MOD ou_input calcium_based_synapse
VERBOSE ${ARB_CAT_VERBOSE}
ADD_DEPS ON)
diff --git a/mechanisms/stochastic/calcium_based_synapse.mod b/mechanisms/stochastic/calcium_based_synapse.mod
new file mode 100644
index 0000000000..4d45a19e1c
--- /dev/null
+++ b/mechanisms/stochastic/calcium_based_synapse.mod
@@ -0,0 +1,99 @@
+: Calcium-based plasticity model
+: Based on the work of Graupner and Brunel, PNAS 109 (10): 3991-3996 (2012)
+: https://doi.org/10.1073/pnas.1109359109, https://www.pnas.org/doi/10.1073/pnas.1220044110
+:
+: Author: Sebastian Schmitt
+:
+: The synapse is modeled as synaptic efficacy variable, ρ, which is a function of the calcium
+: concentration, c(t). The synapse model features two stable states at ρ=0 (DOWN) and ρ=1 (UP),
+: while ρ=ρ_star=0.5 represents a third unstable state between the two stable states.
+: The calcium concentration dynamics are represented by a simplified model which ueses a linear sum
+: of individual calcium transients elicited by trains of pre- and postsynaptic action potentials.
+:
+: drho/dt = -(1/τ)ρ(1-ρ)(ρ_star-ρ) + (γ_p/τ)(1-ρ) H[c(t)-θ_p] - (γ_d/τ)ρ H[c(t)-θ_d] + N
+: N = (σ/√τ) √(H[c(t)-θ_p] + H[c(t)-θ_d]) W
+:
+: dc/dt = -(1/τ_Ca)c + C_pre Σ_i δ(t-t_i-D) + C_post Σ_j δ(t-t_j)
+:
+: rho synaptic efficacy variable (unit-less)
+: rho_star second root of cubic polynomial (unit-less), rho_star=0.5
+: rho_0 initial value (unit-less)
+: tau synaptic time constant (ms), order of seconds to minutes
+: gamma_p rate of synaptic increase (unit-less)
+: theta_p potentiaton threshold (concentration)
+: gamma_d rate of synaptic decrease (unit-less)
+: theta_d depression threshold (concentration)
+: sigma noise amplitude
+: W white noise
+: c calcium concentration (concentration)
+: C_pre concentration jump after pre-synaptic spike (concentration)
+: C_post concentration jump after post-synaptic spike (concentration)
+: tau_Ca Calcium decay time constant (ms), order of milliseconds
+: H right-continuous heaviside step function ( H[x]=1 for x>=0; H[x]=0 otherwise )
+: t_i presynaptic spike times
+: t_j postsynaptic spike times
+: D time delay
+
+NEURON {
+ POINT_PROCESS calcium_based_synapse
+ RANGE rho_0, tau, theta_p, gamma_p, theta_d, gamma_d, C_pre, C_post, tau_Ca, sigma
+}
+
+STATE {
+ c
+ rho
+}
+
+PARAMETER {
+ rho_star = 0.5
+ rho_0 = 1
+ tau = 150000 (ms)
+ gamma_p = 321.808
+ theta_p = 1.3
+ gamma_d = 200
+ theta_d = 1
+ sigma = 2.8248
+ C_pre = 1
+ C_post = 2
+ tau_Ca = 20 (ms)
+}
+
+ASSIGNED {
+ one_over_tau
+ one_over_tau_Ca
+ sigma_over_sqrt_tau
+}
+
+INITIAL {
+ c = 0
+ rho = rho_0
+
+ one_over_tau = 1/tau
+ one_over_tau_Ca = 1/tau_Ca
+ sigma_over_sqrt_tau = sigma/(tau^0.5)
+}
+
+BREAKPOINT {
+ SOLVE state METHOD stochastic
+}
+
+WHITE_NOISE {
+ W
+}
+
+DERIVATIVE state {
+ LOCAL hsp
+ LOCAL hsd
+ hsp = step_right(c - theta_p)
+ hsd = step_right(c - theta_d)
+ rho' = (-rho*(1-rho)*(rho_star-rho) + gamma_p*(1-rho)*hsp - gamma_d*rho*hsd)*one_over_tau + (hsp + hsd)^0.5*sigma_over_sqrt_tau*W
+ c' = -c*one_over_tau_Ca
+}
+
+NET_RECEIVE(weight) {
+ c = c + C_pre
+}
+
+POST_EVENT(time) {
+ c = c + C_post
+}
diff --git a/python/example/calcium_stdp.py b/python/example/calcium_stdp.py
new file mode 100644
index 0000000000..ba86050f8b
--- /dev/null
+++ b/python/example/calcium_stdp.py
@@ -0,0 +1,234 @@
+#!/usr/bin/env python3
+# This script is included in documentation. Adapt line numbers if touched.
+#
+# Authors: Sebastian Schmitt
+# Fabian Bösch
+#
+# Single-cell simulation: Calcium-based synapse which models synaptic efficacy as proposed by
+# Graupner and Brunel, PNAS 109 (10): 3991-3996 (2012); https://doi.org/10.1073/pnas.1109359109,
+# https://www.pnas.org/doi/10.1073/pnas.1220044110.
+# The synapse dynamics is affected by additive white noise. The results reproduce the spike
+# timing-dependent plasticity curve for the DP case described in Table S1 (supplemental material).
+
+import arbor
+import random
+import multiprocessing
+import numpy # You may have to pip install these.
+import pandas # You may have to pip install these.
+import seaborn # You may have to pip install these.
+
+# (1) Set simulation paramters
+
+# Spike response delay (ms)
+D = 13.7
+# Spike frequency in Hertz
+f = 1.0
+# Number of spike pairs
+num_spikes = 30
+# time lag resolution
+stdp_dt_step = 20.0
+# Maximum time lag
+stdp_max_dt = 100.0
+# Ensemble size per initial value
+ensemble_per_rho_0 = 100
+# Simulation time step
+dt = 0.1
+# List of initial values for 2 states
+rho_0 = [0] * ensemble_per_rho_0 + [1] * ensemble_per_rho_0
+# We need a synapse for each sample path
+num_synapses = len(rho_0)
+# Time lags between spike pairs (post-pre: < 0, pre-post: > 0)
+stdp_dt = numpy.arange(-stdp_max_dt, stdp_max_dt + stdp_dt_step, stdp_dt_step)
+
+
+# (2) Make the cell
+
+# Create a morphology with a single (cylindrical) segment of length=diameter=6 μm
+tree = arbor.segment_tree()
+tree.append(arbor.mnpos, arbor.mpoint(-3, 0, 0, 3), arbor.mpoint(3, 0, 0, 3), tag=1)
+
+# Define the soma and its midpoint
+labels = arbor.label_dict({"soma": "(tag 1)", "midpoint": "(location 0 0.5)"})
+
+# Create and set up a decor object
+decor = (
+ arbor.decor()
+ .set_property(Vm=-40)
+ .paint('"soma"', arbor.density("pas"))
+ .place('"midpoint"', arbor.synapse("expsyn"), "driving_synapse")
+ .place('"midpoint"', arbor.threshold_detector(-10), "detector")
+)
+for i in range(num_synapses):
+ mech = arbor.mechanism("calcium_based_synapse")
+ mech.set("rho_0", rho_0[i])
+ decor.place('"midpoint"', arbor.synapse(mech), f"calcium_synapse_{i}")
+
+# Create cell
+cell = arbor.cable_cell(tree, decor, labels)
+
+
+# (3) Create extended catalogue including stochastic mechanisms
+
+cable_properties = arbor.neuron_cable_properties()
+cable_properties.catalogue = arbor.default_catalogue()
+cable_properties.catalogue.extend(arbor.stochastic_catalogue(), "")
+
+
+# (4) Recipe
+
+
+class stdp_recipe(arbor.recipe):
+ def __init__(self, cell, props, gens):
+ arbor.recipe.__init__(self)
+ self.the_cell = cell
+ self.the_props = props
+ self.the_gens = gens
+
+ def num_cells(self):
+ return 1
+
+ def cell_kind(self, gid):
+ return arbor.cell_kind.cable
+
+ def cell_description(self, gid):
+ return self.the_cell
+
+ def global_properties(self, kind):
+ return self.the_props
+
+ def probes(self, gid):
+ return [arbor.cable_probe_point_state_cell("calcium_based_synapse", "rho")]
+
+ def event_generators(self, gid):
+ return self.the_gens
+
+
+# (5) run simulation for a given time lag
+
+
+def run(time_lag):
+
+ # Time between stimuli
+ T = 1000.0 / f
+
+ # Simulation duration
+ t1 = num_spikes * T
+
+ # Time difference between post and pre spike including delay
+ d = -time_lag + D
+
+ # Stimulus and sample times
+ t0_post = 0.0 if d >= 0 else -d
+ t0_pre = d if d >= 0 else 0.0
+ stimulus_times_post = numpy.arange(t0_post, t1, T)
+ stimulus_times_pre = numpy.arange(t0_pre, t1, T)
+ sched_post = arbor.explicit_schedule(stimulus_times_post)
+ sched_pre = arbor.explicit_schedule(stimulus_times_pre)
+
+ # Create strong enough driving stimulus
+ generators = [arbor.event_generator("driving_synapse", 1.0, sched_post)]
+
+ # Stimulus for calcium synapses
+ for i in range(num_synapses):
+ # Zero weight -> just modify synaptic weight via stdp
+ generators.append(arbor.event_generator(f"calcium_synapse_{i}", 0.0, sched_pre))
+
+ # Create recipe
+ recipe = stdp_recipe(cell, cable_properties, generators)
+
+ # Select one thread and no GPU
+ alloc = arbor.proc_allocation(threads=1, gpu_id=None)
+ context = arbor.context(alloc, mpi=None)
+ domains = arbor.partition_load_balance(recipe, context)
+
+ # Get random seed
+ random_seed = random.getrandbits(64)
+
+ # Create simulation
+ sim = arbor.simulation(recipe, context, domains, random_seed)
+
+ # Register prope to read out stdp curve
+ handle = sim.sample((0, 0), arbor.explicit_schedule([t1 - dt]))
+
+ # Run simulation
+ sim.run(t1, dt)
+
+ # Process sampled data
+ data, meta = sim.samples(handle)[0]
+ data_down = data[-1, 1 : ensemble_per_rho_0 + 1]
+ data_up = data[-1, ensemble_per_rho_0 + 1 :]
+ # Initial fraction of synapses in DOWN state
+ beta = 0.5
+ # Synaptic strength ratio UP to DOWN (w1/w0)
+ b = 5
+ # Transition indicator form DOWN to UP
+ P_UA = (data_down > 0.5).astype(float)
+ # Transition indicator from UP to DOWN
+ P_DA = (data_up < 0.5).astype(float)
+ # Return change in synaptic strength
+ ds_A = (
+ (1 - P_UA) * beta
+ + P_DA * (1 - beta)
+ + b * (P_UA * beta + (1 - P_DA) * (1 - beta))
+ ) / (beta + (1 - beta) * b)
+ return pandas.DataFrame({"ds": ds_A, "ms": time_lag, "type": "Arbor"})
+
+
+with multiprocessing.Pool() as p:
+ results = p.map(run, stdp_dt)
+
+ref = numpy.array(
+ [
+ [-100, 0.9793814432989691],
+ [-95, 0.981715028725338],
+ [-90, 0.9932274542583821],
+ [-85, 0.982392230227282],
+ [-80, 0.9620761851689686],
+ [-75, 0.9688482001884063],
+ [-70, 0.9512409611378684],
+ [-65, 0.940405737106768],
+ [-60, 0.9329565205853866],
+ [-55, 0.9146720800329048],
+ [-50, 0.8896156244609853],
+ [-45, 0.9024824529979171],
+ [-40, 0.8252814817763271],
+ [-35, 0.8171550637530018],
+ [-30, 0.7656877496052755],
+ [-25, 0.7176064429672677],
+ [-20, 0.7582385330838939],
+ [-15, 0.7981934216985763],
+ [-10, 0.8835208109434913],
+ [-5, 0.9390513341028807],
+ [0, 0.9927519271849183],
+ [5, 1.2354639175257733],
+ [10, 1.2255075694250952],
+ [15, 1.1760718597832],
+ [20, 1.1862298823123565],
+ [25, 1.1510154042112806],
+ [30, 1.125958948639361],
+ [35, 1.1205413366238108],
+ [40, 1.0812636495110723],
+ [45, 1.0717828284838595],
+ [50, 1.0379227533866708],
+ [55, 1.0392771563905585],
+ [60, 1.023024320343908],
+ [65, 1.046049171409996],
+ [70, 1.040631559394446],
+ [75, 1.0257331263516831],
+ [80, 1.0013538722817072],
+ [85, 1.0121890963128077],
+ [90, 1.0013538722817072],
+ [95, 1.0094802903050326],
+ [100, 0.9918730512544945],
+ ]
+)
+df_ref = pandas.DataFrame({"ds": ref[:, 1], "ms": ref[:, 0], "type": "Reference"})
+
+seaborn.set_theme()
+df = pandas.concat(results)
+df = pandas.concat([df, df_ref])
+plt = seaborn.relplot(kind="line", data=df, x="ms", y="ds", hue="type")
+plt.set_xlabels("lag time difference (ms)")
+plt.set_ylabels("change in synaptic strenght (after/before)")
+plt._legend.set_title("")
+plt.savefig("calcium_stdp.svg")
diff --git a/scripts/run_python_examples.sh b/scripts/run_python_examples.sh
index 3bdc5a7d4c..e561ca58b1 100755
--- a/scripts/run_python_examples.sh
+++ b/scripts/run_python_examples.sh
@@ -32,3 +32,4 @@ $PREFIX python python/example/network_two_cells_gap_junctions.py
$PREFIX python python/example/diffusion.py
$PREFIX python python/example/plasticity.py
$PREFIX python python/example/v-clamp.py
+$PREFIX python python/example/calcium_stdp.py