example_6_distributed_synaptic_input.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example: active cell model adapted from Mainen and Sejnowski (1996)\n",
    "This is an example scripts using LFPy with an active cell model adapted from\n",
    "Mainen and Sejnowski, Nature 1996, for the original files, see\n",
    "http://senselab.med.yale.edu/modeldb/ShowModel.asp?model=2488\n",
    "\n",
    "This scripts is set up to use the model, where the active conductances and\n",
    "spine corrections are set in the function `active_declarations()`, and uses\n",
    "the mechanisms from the .mod-files provided here.\n",
    "For this example to work, run `nrnivmodl` in this folder to compile these\n",
    "mechanisms (e.g . `/$PATHTONEURON/nrn/x86_64/bin/nrnivmodl`).\n",
    "\n",
    "The active parameters are similar to the file `active_declarations_example3.py`\n",
    "\n",
    "Here, excitatory and inhibitory neurons are distributed on different parts of\n",
    "the morphology, with stochastic spike times produced by the\n",
    "LFPy.inputgenerators.stationary_gamma() function.\n",
    "\n",
    "\n",
    "Copyright (C) 2017 Computational Neuroscience Group, NMBU.\n",
    "\n",
    "This program is free software: you can redistribute it and/or modify\n",
    "it under the terms of the GNU General Public License as published by\n",
    "the Free Software Foundation, either version 3 of the License, or\n",
    "(at your option) any later version.\n",
    "\n",
    "This program is distributed in the hope that it will be useful,\n",
    "but WITHOUT ANY WARRANTY; without even the implied warranty of\n",
    "MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the\n",
    "GNU General Public License for more details."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# importing some modules, including LFPy\n",
    "import LFPy\n",
    "import os\n",
    "from os.path import join\n",
    "import sys\n",
    "from urllib.request import urlopen\n",
    "import zipfile\n",
    "import ssl\n",
    "from warnings import warn\n",
    "import numpy as np\n",
    "import scipy.stats\n",
    "import neuron\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams.update({'font.size' : 12,\n",
    "                     'figure.facecolor' : '1',\n",
    "                     'figure.subplot.wspace' : 0.5,\n",
    "                     'figure.subplot.hspace' : 0.5})\n",
    "\n",
    "#seed for random generation\n",
    "np.random.seed(1234)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Fetch Mainen&Sejnowski 1996 model files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if not os.path.isfile(join('cells', 'cells', 'j4a.hoc')):\n",
    "    # get the model files:\n",
    "    u = urlopen('http://senselab.med.yale.edu/ModelDB/eavBinDown.asp?o=2488&a=23&mime=application/zip',\n",
    "                context=ssl._create_unverified_context())\n",
    "    localFile = open('patdemo.zip', 'wb')\n",
    "    localFile.write(u.read())\n",
    "    localFile.close()\n",
    "    # unzip:\n",
    "    myzip = zipfile.ZipFile('patdemo.zip', 'r')\n",
    "    myzip.extractall('.')\n",
    "    myzip.close()\n",
    "\n",
    "# compile mod files every time, because of incompatibility with Hay2011 files:\n",
    "if \"win32\" in sys.platform:\n",
    "    pth = \"cells\"\n",
    "    warn(\"no autompile of NMODL (.mod) files on Windows. \" \n",
    "         + \"Run mknrndll from NEURON bash in the folder cells and rerun example script\")\n",
    "    if not pth in neuron.nrn_dll_loaded:\n",
    "        neuron.h.nrn_load_dll(pth+\"/nrnmech.dll\")\n",
    "    neuron.nrn_dll_loaded.append(pth)\n",
    "else:\n",
    "    os.system('''\n",
    "              cd cells\n",
    "              nrnivmodl\n",
    "              ''')\n",
    "    neuron.load_mechanisms('cells')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A couple of function declarations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def active_declarations(cell):\n",
    "    '''set active conductances and correct for spines,\n",
    "    see file active_declarations_example3.hoc'''\n",
    "    spine_dens = 1\n",
    "    spine_area = 0.83 # // um^2  -- K Harris\n",
    "        \n",
    "    cm_myelin = 0.04\n",
    "    g_pas_node = 0.02\n",
    "    \n",
    "    celsius   = 37.\n",
    "    \n",
    "    Ek = -85.\n",
    "    Ena = 60.\n",
    "    \n",
    "    gna_dend = 20.\n",
    "    gna_node = 30000.\n",
    "    gna_soma = gna_dend * 10\n",
    "    \n",
    "    gkv_axon = 2000.\n",
    "    gkv_soma = 200.\n",
    "    \n",
    "    gca = .3\n",
    "    gkm = .1\n",
    "    gkca = 3\n",
    "    \n",
    "    gca_soma = gca\n",
    "    gkm_soma = gkm\n",
    "    gkca_soma = gkca\n",
    "    \n",
    "    dendritic = neuron.h.SectionList()\n",
    "    for sec in neuron.h.allsec():\n",
    "        if sec.name()[:4] == 'soma':\n",
    "            dendritic.append(sec)\n",
    "        if sec.name()[:4] == 'dend':\n",
    "            dendritic.append(sec)\n",
    "        if sec.name()[:4] == 'apic':\n",
    "            dendritic.append(sec)\n",
    "    \n",
    "    def add_spines(section):\n",
    "        '''add spines as a modification of the compartment area'''\n",
    "        is_spiny = 1\n",
    "        if section == \"dend\":\n",
    "            print(\"adding spines\")\n",
    "            for sec in neuron.h.dendlist:\n",
    "                a = 0\n",
    "                for seg in sec:\n",
    "                    a += neuron.h.area(seg.x)\n",
    "                \n",
    "                F = (sec.L*spine_area*spine_dens + a)/a\n",
    "                sec.L = sec.L * F**(2./3.)\n",
    "                for seg in sec:\n",
    "                    seg.diam = seg.diam * F**(1./3.)\n",
    "                \n",
    "        neuron.h.define_shape()\n",
    "\n",
    "    # Insert active channels\n",
    "    def set_active(cell):\n",
    "        '''set the channel densities'''\n",
    "        # exceptions along the axon\n",
    "        for sec in neuron.h.myelin:\n",
    "            sec.cm = cm_myelin\n",
    "        for sec in neuron.h.node: \n",
    "            sec.g_pas = g_pas_node\n",
    "            \n",
    "        # na+ channels\n",
    "        for sec in neuron.h.allsec():\n",
    "            sec.insert('na')\n",
    "        for sec in dendritic:\n",
    "            sec.gbar_na = gna_dend            \n",
    "        for sec in neuron.h.myelin:\n",
    "            sec.gbar_na = gna_dend\n",
    "        for sec in neuron.h.hill:\n",
    "            sec.gbar_na = gna_node\n",
    "        for sec in neuron.h.iseg:\n",
    "            sec.gbar_na = gna_node\n",
    "        for sec in neuron.h.node:\n",
    "            sec.gbar_na = gna_node\n",
    "                    \n",
    "        # kv delayed rectifier channels\n",
    "        neuron.h.iseg.insert('kv')\n",
    "        neuron.h.iseg.gbar_kv = gkv_axon\n",
    "        neuron.h.hill.insert('kv')\n",
    "        neuron.h.hill.gbar_kv = gkv_axon\n",
    "        for sec in neuron.h.soma:\n",
    "            sec.insert('kv')\n",
    "            sec.gbar_kv = gkv_soma\n",
    "        \n",
    "        # dendritic channels\n",
    "        for sec in dendritic:\n",
    "            sec.insert('km')\n",
    "            sec.gbar_km  = gkm\n",
    "            sec.insert('kca')\n",
    "            sec.gbar_kca = gkca\n",
    "            sec.insert('ca')\n",
    "            sec.gbar_ca = gca\n",
    "            sec.insert('cad')\n",
    "        \n",
    "        # somatic channels\n",
    "        for sec in neuron.h.soma:\n",
    "            sec.gbar_na = gna_soma\n",
    "            sec.gbar_km = gkm_soma\n",
    "            sec.gbar_kca = gkca_soma\n",
    "            sec.gbar_ca = gca_soma\n",
    "        \n",
    "        for sec in neuron.h.allsec():\n",
    "            if neuron.h.ismembrane('k_ion', sec=sec):\n",
    "                sec.ek = Ek\n",
    "        \n",
    "        for sec in neuron.h.allsec():\n",
    "            if neuron.h.ismembrane('na_ion', sec=sec):\n",
    "                sec.ena = Ena\n",
    "                neuron.h.vshift_na = -5\n",
    "                \n",
    "        for sec in neuron.h.allsec():\n",
    "            if neuron.h.ismembrane('ca_ion', sec=sec):\n",
    "                sec.eca = 140\n",
    "                neuron.h.ion_style(\"ca_ion\", 0, 1, 0, 0, 0)\n",
    "                neuron.h.vshift_ca = 0\n",
    "        \n",
    "        #set the temperature for the neuron dynamics\n",
    "        neuron.h.celsius = celsius\n",
    "        \n",
    "        print(\"active ion-channels inserted.\")\n",
    "        \n",
    "        \n",
    "    #// Insert spines\n",
    "    add_spines('dend')\n",
    "    \n",
    "    #// Insert active channels\n",
    "    set_active(cell)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def insert_synapses(synparams, section, n, spTimesFun, args):\n",
    "    '''find n compartments to insert synapses onto'''\n",
    "    idx = cell.get_rand_idx_area_norm(section=section, nidx=n)\n",
    "\n",
    "    #Insert synapses in an iterative fashion\n",
    "    for i in idx:\n",
    "        synparams.update({'idx' : int(i)})\n",
    "\n",
    "        # Some input spike train using the function call\n",
    "        [spiketimes] = spTimesFun(**args)\n",
    "        \n",
    "        # Create synapse(s) and setting times using the Synapse class in LFPy\n",
    "        s = LFPy.Synapse(cell, **synparams)\n",
    "        s.set_spike_times(spiketimes)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parameters etc.:\n",
    "Define parameters, using dictionaries. It is possible to set a few more parameters for each class or functions, but\n",
    "we chose to show only the most important ones here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define cell parameters used as input to cell-class\n",
    "cellParameters = {\n",
    "    'morphology' : join('morphologies', 'L5_Mainen96_wAxon_LFPy.hoc'),\n",
    "    'cm' : 1.0,                 # membrane capacitance\n",
    "    'Ra' : 150,                 # axial resistance\n",
    "    'v_init' : -65,             # initial crossmembrane potential\n",
    "    'passive' : True,           # switch on passive mechs\n",
    "    'passive_parameters' : {'g_pas' : 1./30000, 'e_pas' : -65}, # passive params\n",
    "    'nsegs_method' : 'lambda_f',# method for setting number of segments,\n",
    "    'lambda_f' : 100,           # segments are isopotential at this frequency\n",
    "    'dt' : 2**-4,               # dt of LFP and NEURON simulation.\n",
    "    'tstart' : -100,          #start time, recorders start at t=0\n",
    "    'tstop' : 200,           #stop time of simulation\n",
    "    'custom_fun'  : [active_declarations], # will execute this function\n",
    "    'custom_fun_args' : [{}],\n",
    "}\n",
    "\n",
    "# Synaptic parameters taken from Hendrickson et al 2011\n",
    "# Excitatory synapse parameters:\n",
    "synapseParameters_AMPA = {\n",
    "    'e' : 0,                    #reversal potential\n",
    "    'syntype' : 'Exp2Syn',      #conductance based exponential synapse\n",
    "    'tau1' : 1.,                #Time constant, rise\n",
    "    'tau2' : 3.,                #Time constant, decay\n",
    "    'weight' : 0.005,           #Synaptic weight\n",
    "    'record_current' : True,    #record synaptic currents\n",
    "}\n",
    "# Excitatory synapse parameters\n",
    "synapseParameters_NMDA = {         \n",
    "    'e' : 0,\n",
    "    'syntype' : 'Exp2Syn',\n",
    "    'tau1' : 10.,\n",
    "    'tau2' : 30.,\n",
    "    'weight' : 0.005,\n",
    "    'record_current' : True,\n",
    "}\n",
    "# Inhibitory synapse parameters\n",
    "synapseParameters_GABA_A = {         \n",
    "    'e' : -80,\n",
    "    'syntype' : 'Exp2Syn',\n",
    "    'tau1' : 1.,\n",
    "    'tau2' : 12.,\n",
    "    'weight' : 0.005,\n",
    "    'record_current' : True\n",
    "}\n",
    "# where to insert, how many, and which input statistics\n",
    "insert_synapses_AMPA_args = {\n",
    "    'section' : 'apic',\n",
    "    'n' : 100,\n",
    "    'spTimesFun' : LFPy.inputgenerators.get_activation_times_from_distribution,\n",
    "    'args' : dict(n=1, tstart=0, tstop=cellParameters['tstop'],\n",
    "                  distribution=scipy.stats.gamma,\n",
    "                  rvs_args=dict(a=0.5, loc=0., scale=40)\n",
    "                  )\n",
    "}\n",
    "insert_synapses_NMDA_args = {\n",
    "    'section' : ['dend', 'apic'],\n",
    "    'n' : 15,\n",
    "    'spTimesFun' : LFPy.inputgenerators.get_activation_times_from_distribution,\n",
    "    'args' : dict(n=1, tstart=0, tstop=cellParameters['tstop'],\n",
    "                  distribution=scipy.stats.gamma,\n",
    "                  rvs_args=dict(a=2, loc=0, scale=50)\n",
    "                  )\n",
    "}\n",
    "insert_synapses_GABA_A_args = {\n",
    "    'section' : 'dend',\n",
    "    'n' : 100,\n",
    "    'spTimesFun' : LFPy.inputgenerators.get_activation_times_from_distribution,\n",
    "    'args' : dict(n=1, tstart=0, tstop=cellParameters['tstop'],\n",
    "                  distribution=scipy.stats.gamma,\n",
    "                  rvs_args=dict(a=0.5, loc=0., scale=40)\n",
    "                  )\n",
    "}\n",
    "\n",
    "# Define electrode geometry corresponding to a laminar electrode, where contact\n",
    "# points have a radius r, surface normal vectors N, and LFP calculated as the\n",
    "# average LFP in n random points on each contact:\n",
    "N = np.empty((16, 3))\n",
    "for i in range(N.shape[0]): N[i,] = [1, 0, 0] #normal unit vec. to contacts\n",
    "# put parameters in dictionary\n",
    "electrodeParameters = {\n",
    "    'sigma' : 0.3,              # Extracellular potential\n",
    "    'x' : np.zeros(16) + 25,      # x,y,z-coordinates of electrode contacts\n",
    "    'y' : np.zeros(16),\n",
    "    'z' : np.linspace(-500, 1000, 16),\n",
    "    'n' : 20,\n",
    "    'r' : 10,\n",
    "    'N' : N,\n",
    "}\n",
    "\n",
    "# Parameters for the cell.simulate() call, recording membrane- and syn.-currents\n",
    "simulationParameters = {\n",
    "    'rec_imem' : True,  # Record Membrane currents during simulation\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Main simulation procedure:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Initialize cell instance, using the LFPy.Cell class\n",
    "cell = LFPy.Cell(**cellParameters)\n",
    "\n",
    "# Insert synapses using the function defined earlier\n",
    "insert_synapses(synapseParameters_AMPA, **insert_synapses_AMPA_args)\n",
    "insert_synapses(synapseParameters_NMDA, **insert_synapses_NMDA_args)\n",
    "insert_synapses(synapseParameters_GABA_A, **insert_synapses_GABA_A_args)\n",
    "\n",
    "# perform NEURON simulation, results saved as attributes in the cell instance\n",
    "cell.simulate(**simulationParameters)\n",
    "\n",
    "# Initialize electrode geometry, then calculate the LFP, using the\n",
    "# LFPy.RecExtElectrode class. Note that now cell is given as input to electrode\n",
    "# and created after the NEURON simulations are finished\n",
    "electrode = LFPy.RecExtElectrode(cell, **electrodeParameters)\n",
    "electrode.data = electrode.get_transformation_matrix() @ cell.imem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## plotting some variables and geometry"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from example_suppl import plot_ex3\n",
    "fig = plot_ex3(cell, electrode)\n",
    "# Optional save\n",
    "# fig.savefig('LFPy-example-7.pdf', dpi=300)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}