1# -*- coding: utf-8 -*-
2#
3# intrinsic_currents_spiking.py
4#
5# This file is part of NEST.
6#
7# Copyright (C) 2004 The NEST Initiative
8#
9# NEST is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 2 of the License, or
12# (at your option) any later version.
13#
14# NEST is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with NEST.  If not, see <http://www.gnu.org/licenses/>.
21#
22
23"""
24Intrinsic currents spiking
25--------------------------
26
27This example illustrates a neuron receiving spiking input through
28several different receptors (AMPA, NMDA, GABA_A, GABA_B), provoking
29spike output. The model, ``ht_neuron``, also has intrinsic currents
30(``I_NaP``, ``I_KNa``, ``I_T``, and ``I_h``). It is a slightly simplified implementation of
31neuron model proposed in [1]_.
32
33The neuron is bombarded with spike trains from four Poisson generators,
34which are connected to the AMPA, NMDA, GABA_A, and GABA_B receptors,
35respectively.
36
37References
38~~~~~~~~~~
39
40.. [1] Hill and Tononi (2005) Modeling sleep and wakefulness in the
41       thalamocortical system. J Neurophysiol 93:1671
42       http://dx.doi.org/10.1152/jn.00915.2004.
43
44See Also
45~~~~~~~~
46
47:doc:`intrinsic_currents_subthreshold`
48
49"""
50
51###############################################################################
52# We imported all necessary modules for simulation, analysis and plotting.
53
54import nest
55import matplotlib.pyplot as plt
56
57###############################################################################
58# Additionally, we set the verbosity using ``set_verbosity`` to suppress info
59# messages. We also reset the kernel to be sure to start with a clean NEST.
60
61nest.set_verbosity("M_WARNING")
62nest.ResetKernel()
63
64###############################################################################
65# We define the simulation parameters:
66#
67# - The rate of the input spike trains
68# - The weights of the different receptors (names must match receptor types)
69# - The time to simulate
70#
71# Note that all parameter values should be doubles, since NEST expects doubles.
72
73rate_in = 100.
74w_recep = {'AMPA': 30., 'NMDA': 30., 'GABA_A': 5., 'GABA_B': 10.}
75t_sim = 250.
76
77num_recep = len(w_recep)
78
79###############################################################################
80# We create
81#
82# - one neuron instance
83# - one Poisson generator instance for each synapse type
84# - one multimeter to record from the neuron:
85
86#   - membrane potential
87#   - threshold potential
88#   - synaptic conductances
89#   - intrinsic currents
90#
91# See :doc:`intrinsic_currents_subthreshold` for more details on ``multimeter``
92# configuration.
93
94nrn = nest.Create('ht_neuron')
95p_gens = nest.Create('poisson_generator', 4, params={'rate': rate_in})
96mm = nest.Create('multimeter',
97                 params={'interval': 0.1,
98                         'record_from': ['V_m', 'theta',
99                                         'g_AMPA', 'g_NMDA',
100                                         'g_GABA_A', 'g_GABA_B',
101                                         'I_NaP', 'I_KNa', 'I_T', 'I_h']})
102
103###############################################################################
104# We now connect each Poisson generator with the neuron through a different
105# receptor type.
106#
107# First, we need to obtain the numerical codes for the receptor types from
108# the model. The ``receptor_types`` entry of the default dictionary for the
109# ``ht_neuron`` model is a dictionary mapping receptor names to codes.
110#
111# In the loop, we use Python's tuple unpacking mechanism to unpack
112# dictionary entries from our `w_recep` dictionary.
113#
114# Note that we need to pack the `pg` variable into a list before
115# passing it to ``Connect``, because iterating over the `p_gens` list
116# makes `pg` a "naked" node ID.
117
118receptors = nest.GetDefaults('ht_neuron')['receptor_types']
119for index, (rec_name, rec_wgt) in enumerate(w_recep.items()):
120    nest.Connect(p_gens[index], nrn, syn_spec={'receptor_type': receptors[rec_name], 'weight': rec_wgt})
121
122###############################################################################
123# We then connect the ``multimeter``. Note that the multimeter is connected to
124# the neuron, not the other way around.
125
126nest.Connect(mm, nrn)
127
128###############################################################################
129# We are now ready to simulate.
130
131nest.Simulate(t_sim)
132
133###############################################################################
134# We now fetch the data recorded by the multimeter. The data are returned as
135# a dictionary with entry ``times`` containing timestamps for all
136# recorded data, plus one entry per recorded quantity.
137# All data is contained in the ``events`` entry of the status dictionary
138# returned by the multimeter.
139
140data = mm.events
141t = data['times']
142
143###############################################################################
144# The following function turns a name such as ``I_NaP`` into proper TeX code
145# :math:`I_{\mathrm{NaP}}` for a pretty label.
146
147
148def texify_name(name):
149    return r'${}_{{\mathrm{{{}}}}}$'.format(*name.split('_'))
150
151###############################################################################
152# The next step is to plot the results. We create a new figure, and add one
153# subplot each for membrane and threshold potential, synaptic conductances,
154# and intrinsic currents.
155
156
157fig = plt.figure()
158
159Vax = fig.add_subplot(311)
160Vax.plot(t, data['V_m'], 'b', lw=2, label=r'$V_m$')
161Vax.plot(t, data['theta'], 'g', lw=2, label=r'$\Theta$')
162Vax.set_ylabel('Potential [mV]')
163
164try:
165    Vax.legend(fontsize='small')
166except TypeError:
167    Vax.legend()  # work-around for older Matplotlib versions
168Vax.set_title('ht_neuron driven by Poisson processes')
169
170Gax = fig.add_subplot(312)
171for gname in ('g_AMPA', 'g_NMDA', 'g_GABA_A', 'g_GABA_B'):
172    Gax.plot(t, data[gname], lw=2, label=texify_name(gname))
173
174try:
175    Gax.legend(fontsize='small')
176except TypeError:
177    Gax.legend()  # work-around for older Matplotlib versions
178Gax.set_ylabel('Conductance [nS]')
179
180Iax = fig.add_subplot(313)
181for iname, color in (('I_h', 'maroon'), ('I_T', 'orange'),
182                     ('I_NaP', 'crimson'), ('I_KNa', 'aqua')):
183    Iax.plot(t, data[iname], color=color, lw=2, label=texify_name(iname))
184
185try:
186    Iax.legend(fontsize='small')
187except TypeError:
188    Iax.legend()  # work-around for older Matplotlib versions
189Iax.set_ylabel('Current [pA]')
190Iax.set_xlabel('Time [ms]')
191