{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# IF curve example\n\nThis example illustrates how to measure the I-F curve of a neuron.\nThe program creates a small group of neurons and injects a noisy current\n$I(t) = I_mean + I_std*W(t)$\nwhere $W(t)$ is a white noise process.\nThe programm systematically drives the current through a series of  values in\nthe two-dimensional `(I_mean, I_std)` space and measures the firing rate of\nthe neurons.\n\nIn this example, we measure the I-F curve of the adaptive exponential\nintegrate and fire neuron (``aeif_cond_exp``), but any other neuron model that\naccepts current inputs is possible. The model and its parameters are\nsupplied when the IF_curve object is created.\n\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import numpy\nimport nest\nimport shelve"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Here we define which model and the neuron parameters to use for measuring\nthe transfer function.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "model = 'aeif_cond_exp'\nparams = {'a': 4.0,\n          'b': 80.8,\n          'V_th': -50.4,\n          'Delta_T': 2.0,\n          'I_e': 0.0,\n          'C_m': 281.0,\n          'g_L': 30.0,\n          'V_reset': -70.6,\n          'tau_w': 144.0,\n          't_ref': 5.0,\n          'V_peak': -40.0,\n          'E_L': -70.6,\n          'E_ex': 0.,\n          'E_in': -70.}\n\n\nclass IF_curve():\n\n    t_inter_trial = 200.  # Interval between two successive measurement trials\n    t_sim = 1000.         # Duration of a measurement trial\n    n_neurons = 100       # Number of neurons\n    n_threads = 4         # Nubmer of threads to run the simulation\n\n    def __init__(self, model, params=False):\n        self.model = model\n        self.params = params\n        self.build()\n        self.connect()\n\n    def build(self):\n        #######################################################################\n        #  We reset NEST to delete information from previous simulations\n        # and adjust the number of threads.\n\n        nest.ResetKernel()\n        nest.SetKernelStatus({'local_num_threads': self.n_threads})\n\n        #######################################################################\n        # We set the default parameters of the neuron model to those\n        # defined above and create neurons and devices.\n\n        if self.params:\n            nest.SetDefaults(self.model, self.params)\n        self.neuron = nest.Create(self.model, self.n_neurons)\n        self.noise = nest.Create('noise_generator')\n        self.spike_recorder = nest.Create('spike_recorder')\n\n    def connect(self):\n        #######################################################################\n        # We connect the noisy current to the neurons and the neurons to\n        # the spike recorders.\n\n        nest.Connect(self.noise, self.neuron, 'all_to_all')\n        nest.Connect(self.neuron, self.spike_recorder, 'all_to_all')\n\n    def output_rate(self, mean, std):\n        self.build()\n        self.connect()\n\n        #######################################################################\n        # We adjust the parameters of the noise according to the current\n        # values.\n\n        self.noise.set(mean=mean, std=std, start=0.0, stop=1000., origin=0.)\n\n        # We simulate the network and calculate the rate.\n\n        nest.Simulate(self.t_sim)\n        rate = self.spike_recorder.n_events * 1000. / (1. * self.n_neurons * self.t_sim)\n        return rate\n\n    def compute_transfer(self, i_mean=(400.0, 900.0, 50.0),\n                         i_std=(0.0, 600.0, 50.0)):\n        #######################################################################\n        # We loop through all possible combinations of `(I_mean, I_sigma)`\n        # and measure the output rate of the neuron.\n\n        self.i_range = numpy.arange(*i_mean)\n        self.std_range = numpy.arange(*i_std)\n        self.rate = numpy.zeros((self.i_range.size, self.std_range.size))\n        nest.set_verbosity('M_WARNING')\n        for n, i in enumerate(self.i_range):\n            print('I  =  {0}'.format(i))\n            for m, std in enumerate(self.std_range):\n                self.rate[n, m] = self.output_rate(i, std)\n\n\ntransfer = IF_curve(model, params)\ntransfer.compute_transfer()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "After the simulation is finished, we store the data into a file for\nlater analysis.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "with shelve.open(model + '_transfer.dat') as dat:\n    dat['I_mean'] = transfer.i_range\n    dat['I_std'] = transfer.std_range\n    dat['rate'] = transfer.rate"
      ]
    }
  ],
  "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.6.12"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}