{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Clopath Rule: Bidirectional connections\n\nThis script simulates a small network of ten excitatory and three\ninhibitory ``aeif_psc_delta_clopath`` neurons. The neurons are randomly connected\nand driven by 500 Poisson generators. The synapses from the Poisson generators\nto the excitatory population and those among the neurons of the network\nare Clopath synapses. The rate of the Poisson generators is modulated with\na Gaussian profile whose center shifts randomly each 100 ms between ten\nequally spaced positions.\nThis setup demonstrates that the Clopath synapse is able to establish\nbidirectional connections. The example is adapted from [1]_ (cf. fig. 5).\n\n## References\n\n.. [1] Clopath C, B\u00fcsing L, Vasilaki E, Gerstner W (2010). Connectivity reflects coding:\n       a model of voltage-based STDP with homeostasis.\n       Nature Neuroscience 13:3, 344--352\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import nest\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport random"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Set the parameters\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "simulation_time = 1.0e4\nresolution = 0.1\ndelay = resolution\n\n# Poisson_generator parameters\npg_A = 30.      # amplitude of Gaussian\npg_sigma = 10.  # std deviation\n\nnest.ResetKernel()\nnest.resolution = resolution\n\n# Create neurons and devices\nnrn_model = 'aeif_psc_delta_clopath'\nnrn_params = {'V_m': -30.6,\n              'g_L': 30.0,\n              'w': 0.0,\n              'tau_plus': 7.0,\n              'tau_minus': 10.0,\n              'tau_w': 144.0,\n              'a': 4.0,\n              'C_m': 281.0,\n              'Delta_T': 2.0,\n              'V_peak': 20.0,\n              't_clamp': 2.0,\n              'A_LTP': 8.0e-6,\n              'A_LTD': 14.0e-6,\n              'A_LTD_const': False,\n              'b': 0.0805,\n              'u_ref_squared': 60.0**2}\n\npop_exc = nest.Create(nrn_model, 10, nrn_params)\npop_inh = nest.Create(nrn_model, 3, nrn_params)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "We need parrot neurons since Poisson generators can only be connected\nwith static connections\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "pop_input = nest.Create('parrot_neuron', 500)  # helper neurons\npg = nest.Create('poisson_generator', 500)\nwr = nest.Create('weight_recorder')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "First connect Poisson generators to helper neurons\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "nest.Connect(pg, pop_input, 'one_to_one', {'synapse_model': 'static_synapse',\n                                           'weight': 1.0, 'delay': delay})"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Create all the connections\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "nest.CopyModel('clopath_synapse', 'clopath_input_to_exc',\n               {'Wmax': 3.0})\nconn_dict_input_to_exc = {'rule': 'all_to_all'}\nsyn_dict_input_to_exc = {'synapse_model': 'clopath_input_to_exc',\n                         'weight': nest.random.uniform(0.5, 2.0),\n                         'delay': delay}\nnest.Connect(pop_input, pop_exc, conn_dict_input_to_exc,\n             syn_dict_input_to_exc)\n\n# Create input->inh connections\nconn_dict_input_to_inh = {'rule': 'all_to_all'}\nsyn_dict_input_to_inh = {'synapse_model': 'static_synapse',\n                         'weight': nest.random.uniform(0.0, 0.5),\n                         'delay': delay}\nnest.Connect(pop_input, pop_inh, conn_dict_input_to_inh, syn_dict_input_to_inh)\n\n# Create exc->exc connections\nnest.CopyModel('clopath_synapse', 'clopath_exc_to_exc',\n               {'Wmax': 0.75, 'weight_recorder': wr})\nsyn_dict_exc_to_exc = {'synapse_model': 'clopath_exc_to_exc', 'weight': 0.25,\n                       'delay': delay}\nconn_dict_exc_to_exc = {'rule': 'all_to_all', 'allow_autapses': False}\nnest.Connect(pop_exc, pop_exc, conn_dict_exc_to_exc, syn_dict_exc_to_exc)\n\n# Create exc->inh connections\nsyn_dict_exc_to_inh = {'synapse_model': 'static_synapse',\n                       'weight': 1.0, 'delay': delay}\nconn_dict_exc_to_inh = {'rule': 'fixed_indegree', 'indegree': 8}\nnest.Connect(pop_exc, pop_inh, conn_dict_exc_to_inh, syn_dict_exc_to_inh)\n\n# Create inh->exc connections\nsyn_dict_inh_to_exc = {'synapse_model': 'static_synapse',\n                       'weight': 1.0, 'delay': delay}\nconn_dict_inh_to_exc = {'rule': 'fixed_outdegree', 'outdegree': 6}\nnest.Connect(pop_inh, pop_exc, conn_dict_inh_to_exc, syn_dict_inh_to_exc)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Randomize the initial membrane potential\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "pop_exc.V_m = nest.random.normal(-60., 25.)\npop_inh.V_m = nest.random.normal(-60., 25.)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Simulation divided into intervals of 100ms for shifting the Gaussian\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "sim_interval = 100.\nfor i in range(int(simulation_time / sim_interval)):\n    # set rates of poisson generators\n    rates = np.empty(500)\n    # pg_mu will be randomly chosen out of 25,75,125,...,425,475\n    pg_mu = 25 + random.randint(0, 9) * 50\n    for j in range(500):\n        rates[j] = pg_A * np.exp((-1 * (j - pg_mu)**2) / (2 * pg_sigma**2))\n        pg[j].rate = rates[j] * 1.75\n    nest.Simulate(sim_interval)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Plot results\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, ax = plt.subplots(1, sharex=False)\n\n# Plot synapse weights of the synapses within the excitatory population\n# Sort weights according to sender and reshape\nexc_conns = nest.GetConnections(pop_exc, pop_exc)\nexc_conns_senders = np.array(exc_conns.source)\nexc_conns_targets = np.array(exc_conns.target)\nexc_conns_weights = np.array(exc_conns.weight)\nidx_array = np.argsort(exc_conns_senders)\ntargets = np.reshape(exc_conns_targets[idx_array], (10, 10 - 1))\nweights = np.reshape(exc_conns_weights[idx_array], (10, 10 - 1))\n\n# Sort according to target\nfor i, (trgs, ws) in enumerate(zip(targets, weights)):\n    idx_array = np.argsort(trgs)\n    weights[i] = ws[idx_array]\n\nweight_matrix = np.zeros((10, 10))\ntu9 = np.triu_indices_from(weights)\ntl9 = np.tril_indices_from(weights, -1)\ntu10 = np.triu_indices_from(weight_matrix, 1)\ntl10 = np.tril_indices_from(weight_matrix, -1)\nweight_matrix[tu10[0], tu10[1]] = weights[tu9[0], tu9[1]]\nweight_matrix[tl10[0], tl10[1]] = weights[tl9[0], tl9[1]]\n\n# Difference between initial and final value\ninit_w_matrix = np.ones((10, 10)) * 0.25\ninit_w_matrix -= np.identity(10) * 0.25\n\ncax = ax.imshow(weight_matrix - init_w_matrix)\ncbarB = fig.colorbar(cax, ax=ax)\nax.set_xticks([0, 2, 4, 6, 8])\nax.set_xticklabels(['1', '3', '5', '7', '9'])\nax.set_yticks([0, 2, 4, 6, 8])\nax.set_xticklabels(['1', '3', '5', '7', '9'])\nax.set_xlabel(\"to neuron\")\nax.set_ylabel(\"from neuron\")\nax.set_title(\"Change of syn weights before and after simulation\")\nplt.show()"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.8.6"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}