Skip to content

Commit

Permalink
Merge pull request #3028 from schmitfe/add_EI_clustered_network_example
Browse files Browse the repository at this point in the history
Add EI_clustered_network by Rostami et al to PyNEST examples
  • Loading branch information
heplesser committed Apr 10, 2024
2 parents 799ab2d + 196784f commit 56c720c
Show file tree
Hide file tree
Showing 9 changed files with 1,346 additions and 18 deletions.
44 changes: 26 additions & 18 deletions doc/htmldoc/examples/index.rst
Expand Up @@ -80,31 +80,31 @@ PyNEST examples

* :doc:`../auto_examples/Potjans_2014/index`



.. grid-item-card:: EI clustered network (Rostami et al)
:img-top: ../static/img/pynest/EI_clustered_network_schematic.png

:doc:`../auto_examples/EI_clustered_network/index`



.. grid:: 1 1 2 3

.. grid-item-card:: GLIF (from Allen institute)
:img-top: ../static/img/pynest/glif_cond.png

* :doc:`../auto_examples/glif_cond_neuron`
* :doc:`../auto_examples/glif_psc_neuron`
* :doc:`../auto_examples/glif_psc_double_alpha_neuron`



.. grid:: 1 1 2 3

.. grid-item-card:: Compartmental neurons
:img-top: ../static/img/pynest/dendritic_synapse_conductances.png

* :doc:`../auto_examples/compartmental_model/receptors_and_current`
* :doc:`../auto_examples/compartmental_model/two_comps`


.. grid-item-card:: Rate neurons
:img-top: ../static/img/pynest/rate_neuron.png

* :doc:`../auto_examples/lin_rate_ipn_network`
* :doc:`../auto_examples/rate_neuron_dm`



.. grid-item-card:: GIF (from Gerstner lab)
:img-top: ../static/img/pynest/gif_pop.png
Expand All @@ -116,6 +116,13 @@ PyNEST examples

.. grid:: 1 1 2 3

.. grid-item-card:: Rate neurons
:img-top: ../static/img/pynest/rate_neuron.png

* :doc:`../auto_examples/lin_rate_ipn_network`
* :doc:`../auto_examples/rate_neuron_dm`


.. grid-item-card:: Hodgkin-Huxley
:img-top: ../static/img/pynest/hh_phase.png

Expand All @@ -127,14 +134,14 @@ PyNEST examples

* :doc:`../auto_examples/BrodyHopfield`

.. grid:: 1 1 2 3

.. grid-item-card:: Brette and Gerstner
:img-top: ../static/img/pynest/brette_gerstner2c.png

* :doc:`../auto_examples/brette_gerstner_fig_2c`
* :doc:`../auto_examples/brette_gerstner_fig_3d`

.. grid:: 1 1 2 3


.. grid-item-card:: Precise spiking
:img-top: ../static/img/pynest/precisespiking.png
Expand All @@ -146,11 +153,6 @@ PyNEST examples

* :doc:`../auto_examples/CampbellSiegert`

.. grid-item-card:: SONATA network
:img-top: ../static/img/300_pointneurons.png

* :doc:`../auto_examples/sonata_example/sonata_network`


.. grid:: 1 1 2 3

Expand Down Expand Up @@ -242,6 +244,11 @@ PyNEST examples

.. grid:: 1 1 2 3

.. grid-item-card:: SONATA network
:img-top: ../static/img/300_pointneurons.png

* :doc:`../auto_examples/sonata_example/sonata_network`

.. grid-item-card:: HPC benchmark
:img-top: ../static/img/nest_logo-faded.png

Expand Down Expand Up @@ -342,6 +349,7 @@ PyNEST examples
../auto_examples/astrocytes/astrocyte_interaction
../auto_examples/astrocytes/astrocyte_small_network
../auto_examples/astrocytes/astrocyte_brunel
../auto_examples/EI_clustered_network/index
../auto_examples/eprop_plasticity/index

.. toctree::
Expand Down
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
69 changes: 69 additions & 0 deletions pynest/examples/EI_clustered_network/README.rst
@@ -0,0 +1,69 @@
EI-clustered circuit model
==========================

This is PyNEST implementation of the EI-clustered circuit model described by Rostami et al. [1]_.

.. figure:: /static/img/pynest/EI_clustered_network_schematic.png
:alt: EI-clustered circuit model.

Schematic of the EI-clustered circuit model. The network consists of `n_clusters` with one excitatory and one inhibitory population each.

Citing this code
----------------

If you use this code, we ask you to cite the paper by Rostami et al. [1]_ and the NEST release on Zenodo.

File structure
--------------

* :doc:`run_simulation.py <run_simulation>`: an example script to try out the EI-clustered circuit model
* :doc:`network.py <network>`: the main ``Network`` class with functions to build and simulate the network
* :doc:`helper.py <helper>`: helper functions for calculation of synaptic weights and currents and plot function for raster plots
* :doc:`network_params.py <network_params>`: network and neuron parameters
* :doc:`stimulus_params.py <stimulus_params>`: parameters for optional external stimulation
* :doc:`sim_params.py <sim_params>`: simulation parameters

Running the simulation
----------------------

.. code-block:: bash
python run_simulation.py
A raster plot of the network activity is saved as ``clustered_ei_raster.png``.

The code can be parallelized by using multiple threads during the NEST simulation.
This can be done by setting the parameter ``n_vp`` in the ``run_simulation.py`` script.

Contributions to this PyNEST model implementation
-------------------------------------------------

2023: initial version of code and documentation by Felix J. Schmitt, Vahid Rostami and Martin Nawrot.

Acknowledgments
---------------

Funding for the study by Rostami et al. [1]_: This work was supported by the German Research Foundation (DFG),
in parts through the Collaborative Research Center ’Motor Control in Health and Disease’
(DFG-SFB 1451, Project-ID 431549029) and under the Institutional Strategy of the University of Cologne within the
German Excellence Initiative (DFG-ZUK 81/1) and in parts through the DFG graduate school
’Neural Circuit Analysis’ (DFG-RTG 1960, ID 365082554) and through the European Union’s Horizon 2020 Framework
Programme for Research and Innovation under grant agreement number 945539 (Human Brain Project SGA3).
The figure is created with BioRender.com.

Other implementations of the EI-clustered model
-----------------------------------------------

A `GeNN version <https://github.com/nawrotlab/SNN_GeNN_Nest>`__ by Felix J. Schmitt, Vahid Rostami and Martin Nawrot [2]_.

References
----------

.. [1] Rostami V, Rost T, Riehle A, van Albada SJ and Nawrot MP. 2020.
Excitatory and inhibitory motor cortical clusters account for balance, variability, and task performance.
bioRxiv 2020.02.27.968339. DOI: `10.1101/2020.02.27.968339 <https://doi.org/10.1101/2020.02.27.968339>`__.
.. [2] Schmitt FJ, Rostami V and Nawrot MP. 2023.
Efficient parameter calibration and real-time simulation of large-scale spiking neural networks with GeNN
and NEST. Front. Neuroinform. 17:941696. DOI: `10.3389/fninf.2023.941696 <https://doi.org/10.3389/fninf.2023.941696>`__.
190 changes: 190 additions & 0 deletions pynest/examples/EI_clustered_network/helper.py
@@ -0,0 +1,190 @@
# -*- coding: utf-8 -*-
#
# helper.py
#
# This file is part of NEST.
#
# Copyright (C) 2004 The NEST Initiative
#
# NEST is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 2 of the License, or
# (at your option) any later version.
#
# NEST is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with NEST. If not, see <http://www.gnu.org/licenses/>.

"""PyNEST EI-clustered network: Helper Functions
------------------------------------------------
Helper functions to calculate synaptic weights to construct
random balanced networks and plot raster plot with color groups.
"""

import matplotlib.pyplot as plt
import numpy as np


def postsynaptic_current_to_potential(tau_m, tau_syn, c_m=1.0, e_l=0.0):
"""Maximum post-synaptic potential amplitude
for exponential synapses and a synaptic efficacy J of 1 pA.
Parameters
----------
tau_m: float
Membrane time constant [ms]
tau_syn: float
Synapse time constant [ms]
c_m: float (optional)
Membrane capacity [pF] (default: 1.0)
e_l: float (optional)
Resting potential [mV] (default: 0.0)
Returns
-------
float
maximum psp amplitude (for a 1 pA current spike) [mV]
"""
# time of maximum deflection of the psp [ms]
tmax = np.log(tau_syn / tau_m) / (1 / tau_m - 1 / tau_syn)
# we assume here the current spike is 1 pA, otherwise [mV/pA]
pre = tau_m * tau_syn / c_m / (tau_syn - tau_m)
return (e_l - pre) * np.exp(-tmax / tau_m) + pre * np.exp(-tmax / tau_syn)


def calculate_RBN_weights(params):
"""Calculate synaptic weights for a random balanced network.
The synaptic weights are calculated according to the method
described in Eqs. 7-10 [1]_.
References
----------
.. [1] Rostami V, Rost T, Riehle A, van Albada SJ and Nawrot MP. 2020.
Excitatory and inhibitory motor cortical clusters account for balance,
variability, and task performance. bioRxiv 2020.02.27.968339.
DOI: `10.1101/2020.02.27.968339
<https://doi.org/10.1101/2020.02.27.968339>`__.
Parameters
----------
params: dict
Dictionary of network parameters
Returns
-------
ndarray
synaptic weights 2x2 matrix [[EE, EI], [IE, II]]
"""

N_E = params.get("N_E") # excitatory units
N_I = params.get("N_I") # inhibitory units
N = N_E + N_I # total units

E_L = params.get("E_L")
V_th_E = params.get("V_th_E") # threshold voltage
V_th_I = params.get("V_th_I")

tau_E = params.get("tau_E")
tau_I = params.get("tau_I")

tau_syn_ex = params.get("tau_syn_ex")
tau_syn_in = params.get("tau_syn_in")

gei = params.get("gei")
gii = params.get("gii")
gie = params.get("gie")

amp_EE = postsynaptic_current_to_potential(tau_E, tau_syn_ex)
amp_EI = postsynaptic_current_to_potential(tau_E, tau_syn_in)
amp_IE = postsynaptic_current_to_potential(tau_I, tau_syn_ex)
amp_II = postsynaptic_current_to_potential(tau_I, tau_syn_in)

baseline_conn_prob = params.get("baseline_conn_prob") # connection probs

js = np.zeros((2, 2))
K_EE = N_E * baseline_conn_prob[0, 0]
js[0, 0] = (V_th_E - E_L) * (K_EE**-0.5) * N**0.5 / amp_EE
js[0, 1] = -gei * js[0, 0] * baseline_conn_prob[0, 0] * N_E * amp_EE / (baseline_conn_prob[0, 1] * N_I * amp_EI)
K_IE = N_E * baseline_conn_prob[1, 0]
js[1, 0] = gie * (V_th_I - E_L) * (K_IE**-0.5) * N**0.5 / amp_IE
js[1, 1] = -gii * js[1, 0] * baseline_conn_prob[1, 0] * N_E * amp_IE / (baseline_conn_prob[1, 1] * N_I * amp_II)
return js


def rheobase_current(tau_m, e_l, v_th, c_m):
"""Rheobase current for membrane time constant and resting potential.
Parameters
----------
tau_m: float
Membrane time constant [ms]
e_l: float
Resting potential [mV]
v_th: float
Threshold potential [mV]
c_m: float
Membrane capacity [pF]
Returns
-------
float
rheobase current [pA]
"""
return (v_th - e_l) * c_m / tau_m


def raster_plot(spiketimes, tlim=None, colorgroups=None, ax=None, markersize=0.5):
"""Raster plot of spiketimes.
Plots raster plot of spiketimes withing given time limits and
colors neurons according to colorgroups.
Parameters
----------
spiketimes: ndarray
2D array [2xN_Spikes]
of spiketimes with spiketimes in row 0 and neuron IDs in row 1.
tlim: list of floats (optional)
Time limits of plot: [tmin, tmax],
if None: [min(spiketimes), max(spiketimes)]
colorgroups: list of tuples (optional)
Each element is a tuple in the format
(color, start_neuron_ID, stop_neuron_ID]) for coloring neurons in
given range, if None is given all neurons are black.
ax: axis object (optional)
If None a new figure is created,
else the plot is added to the given axis.
markersize: float (optional)
Size of markers. (default: 0.5)
Returns
-------
axis object
Axis object with raster plot.
"""
if ax is None:
fig, ax = plt.subplots()
if tlim is None:
tlim = [min(spiketimes[0]), max(spiketimes[0])]
if colorgroups is None:
colorgroups = [("k", 0, max(spiketimes[1]))]
for color, start, stop in colorgroups:
ax.plot(
spiketimes[0][np.logical_and(spiketimes[1] >= start, spiketimes[1] < stop)],
spiketimes[1][np.logical_and(spiketimes[1] >= start, spiketimes[1] < stop)],
color=color,
marker=".",
linestyle="None",
markersize=markersize,
)
ax.set_xlim(tlim)
ax.set_ylim([0, max(spiketimes[1])])
ax.set_xlabel("Time [ms]")
ax.set_ylabel("Neuron ID")
return ax

0 comments on commit 56c720c

Please sign in to comment.