1# -*- coding: utf-8 -*-
2#
3# clopath_synapse_small_network.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"""
23Clopath Rule: Bidirectional connections
24---------------------------------------
25
26This script simulates a small network of ten excitatory and three
27inhibitory ``aeif_psc_delta_clopath`` neurons. The neurons are randomly connected
28and driven by 500 Poisson generators. The synapses from the Poisson generators
29to the excitatory population and those among the neurons of the network
30are Clopath synapses. The rate of the Poisson generators is modulated with
31a Gaussian profile whose center shifts randomly each 100 ms between ten
32equally spaced positions.
33This setup demonstrates that the Clopath synapse is able to establish
34bidirectional connections. The example is adapted from [1]_ (cf. fig. 5).
35
36References
37~~~~~~~~~~
38
39.. [1] Clopath C, Büsing L, Vasilaki E, Gerstner W (2010). Connectivity reflects coding:
40       a model of voltage-based STDP with homeostasis.
41       Nature Neuroscience 13:3, 344--352
42"""
43
44import nest
45import numpy as np
46import matplotlib.pyplot as plt
47import random
48
49##############################################################################
50# Set the parameters
51
52simulation_time = 1.0e4
53resolution = 0.1
54delay = resolution
55
56# Poisson_generator parameters
57pg_A = 30.      # amplitude of Gaussian
58pg_sigma = 10.  # std deviation
59
60nest.ResetKernel()
61nest.resolution = resolution
62
63# Create neurons and devices
64nrn_model = 'aeif_psc_delta_clopath'
65nrn_params = {'V_m': -30.6,
66              'g_L': 30.0,
67              'w': 0.0,
68              'tau_plus': 7.0,
69              'tau_minus': 10.0,
70              'tau_w': 144.0,
71              'a': 4.0,
72              'C_m': 281.0,
73              'Delta_T': 2.0,
74              'V_peak': 20.0,
75              't_clamp': 2.0,
76              'A_LTP': 8.0e-6,
77              'A_LTD': 14.0e-6,
78              'A_LTD_const': False,
79              'b': 0.0805,
80              'u_ref_squared': 60.0**2}
81
82pop_exc = nest.Create(nrn_model, 10, nrn_params)
83pop_inh = nest.Create(nrn_model, 3, nrn_params)
84
85##############################################################################
86# We need parrot neurons since Poisson generators can only be connected
87# with static connections
88
89pop_input = nest.Create('parrot_neuron', 500)  # helper neurons
90pg = nest.Create('poisson_generator', 500)
91wr = nest.Create('weight_recorder')
92
93##############################################################################
94# First connect Poisson generators to helper neurons
95nest.Connect(pg, pop_input, 'one_to_one', {'synapse_model': 'static_synapse',
96                                           'weight': 1.0, 'delay': delay})
97
98##############################################################################
99# Create all the connections
100
101nest.CopyModel('clopath_synapse', 'clopath_input_to_exc',
102               {'Wmax': 3.0})
103conn_dict_input_to_exc = {'rule': 'all_to_all'}
104syn_dict_input_to_exc = {'synapse_model': 'clopath_input_to_exc',
105                         'weight': nest.random.uniform(0.5, 2.0),
106                         'delay': delay}
107nest.Connect(pop_input, pop_exc, conn_dict_input_to_exc,
108             syn_dict_input_to_exc)
109
110# Create input->inh connections
111conn_dict_input_to_inh = {'rule': 'all_to_all'}
112syn_dict_input_to_inh = {'synapse_model': 'static_synapse',
113                         'weight': nest.random.uniform(0.0, 0.5),
114                         'delay': delay}
115nest.Connect(pop_input, pop_inh, conn_dict_input_to_inh, syn_dict_input_to_inh)
116
117# Create exc->exc connections
118nest.CopyModel('clopath_synapse', 'clopath_exc_to_exc',
119               {'Wmax': 0.75, 'weight_recorder': wr})
120syn_dict_exc_to_exc = {'synapse_model': 'clopath_exc_to_exc', 'weight': 0.25,
121                       'delay': delay}
122conn_dict_exc_to_exc = {'rule': 'all_to_all', 'allow_autapses': False}
123nest.Connect(pop_exc, pop_exc, conn_dict_exc_to_exc, syn_dict_exc_to_exc)
124
125# Create exc->inh connections
126syn_dict_exc_to_inh = {'synapse_model': 'static_synapse',
127                       'weight': 1.0, 'delay': delay}
128conn_dict_exc_to_inh = {'rule': 'fixed_indegree', 'indegree': 8}
129nest.Connect(pop_exc, pop_inh, conn_dict_exc_to_inh, syn_dict_exc_to_inh)
130
131# Create inh->exc connections
132syn_dict_inh_to_exc = {'synapse_model': 'static_synapse',
133                       'weight': 1.0, 'delay': delay}
134conn_dict_inh_to_exc = {'rule': 'fixed_outdegree', 'outdegree': 6}
135nest.Connect(pop_inh, pop_exc, conn_dict_inh_to_exc, syn_dict_inh_to_exc)
136
137##############################################################################
138# Randomize the initial membrane potential
139
140pop_exc.V_m = nest.random.normal(-60., 25.)
141pop_inh.V_m = nest.random.normal(-60., 25.)
142
143##############################################################################
144# Simulation divided into intervals of 100ms for shifting the Gaussian
145
146sim_interval = 100.
147for i in range(int(simulation_time / sim_interval)):
148    # set rates of poisson generators
149    rates = np.empty(500)
150    # pg_mu will be randomly chosen out of 25,75,125,...,425,475
151    pg_mu = 25 + random.randint(0, 9) * 50
152    for j in range(500):
153        rates[j] = pg_A * np.exp((-1 * (j - pg_mu)**2) / (2 * pg_sigma**2))
154        pg[j].rate = rates[j] * 1.75
155    nest.Simulate(sim_interval)
156
157##############################################################################
158# Plot results
159
160fig, ax = plt.subplots(1, sharex=False)
161
162# Plot synapse weights of the synapses within the excitatory population
163# Sort weights according to sender and reshape
164exc_conns = nest.GetConnections(pop_exc, pop_exc)
165exc_conns_senders = np.array(exc_conns.source)
166exc_conns_targets = np.array(exc_conns.target)
167exc_conns_weights = np.array(exc_conns.weight)
168idx_array = np.argsort(exc_conns_senders)
169targets = np.reshape(exc_conns_targets[idx_array], (10, 10 - 1))
170weights = np.reshape(exc_conns_weights[idx_array], (10, 10 - 1))
171
172# Sort according to target
173for i, (trgs, ws) in enumerate(zip(targets, weights)):
174    idx_array = np.argsort(trgs)
175    weights[i] = ws[idx_array]
176
177weight_matrix = np.zeros((10, 10))
178tu9 = np.triu_indices_from(weights)
179tl9 = np.tril_indices_from(weights, -1)
180tu10 = np.triu_indices_from(weight_matrix, 1)
181tl10 = np.tril_indices_from(weight_matrix, -1)
182weight_matrix[tu10[0], tu10[1]] = weights[tu9[0], tu9[1]]
183weight_matrix[tl10[0], tl10[1]] = weights[tl9[0], tl9[1]]
184
185# Difference between initial and final value
186init_w_matrix = np.ones((10, 10)) * 0.25
187init_w_matrix -= np.identity(10) * 0.25
188
189cax = ax.imshow(weight_matrix - init_w_matrix)
190cbarB = fig.colorbar(cax, ax=ax)
191ax.set_xticks([0, 2, 4, 6, 8])
192ax.set_xticklabels(['1', '3', '5', '7', '9'])
193ax.set_yticks([0, 2, 4, 6, 8])
194ax.set_xticklabels(['1', '3', '5', '7', '9'])
195ax.set_xlabel("to neuron")
196ax.set_ylabel("from neuron")
197ax.set_title("Change of syn weights before and after simulation")
198plt.show()
199