diff --git a/models/CMakeLists.txt b/models/CMakeLists.txt index b1c4af83a8..9468de4c94 100644 --- a/models/CMakeLists.txt +++ b/models/CMakeLists.txt @@ -88,6 +88,7 @@ set( models_sources tsodyks_connection.h tsodyks_connection_hom.h tsodyks_connection_hom.cpp volume_transmitter.h volume_transmitter.cpp + vogels_sprekeler_connection.h spike_dilutor.h spike_dilutor.cpp ) diff --git a/models/modelsmodule.cpp b/models/modelsmodule.cpp index 24a3864901..bd40742e13 100644 --- a/models/modelsmodule.cpp +++ b/models/modelsmodule.cpp @@ -115,6 +115,7 @@ #include "tsodyks2_connection.h" #include "tsodyks_connection.h" #include "tsodyks_connection_hom.h" +#include "vogels_sprekeler_connection.h" // Includes from nestkernel: #include "common_synapse_properties.h" @@ -571,6 +572,21 @@ ModelsModule::init( SLIInterpreter* ) .model_manager .register_connection_model< STDPDopaConnection< TargetIdentifierIndex > >( "stdp_dopamine_synapse_hpc" ); + + /* BeginDocumentation + Name: vogels_sprekeler_synapse_hpc - Variant of vogels_sprekeler_synapse + with low memory + consumption. + SeeAlso: synapsedict, vogels_sprekeler_synapse + */ + kernel() + .model_manager + .register_connection_model< VogelsSprekelerConnection< TargetIdentifierPtrRport > >( + "vogels_sprekeler_synapse" ); + kernel() + .model_manager + .register_connection_model< VogelsSprekelerConnection< TargetIdentifierIndex > >( + "vogels_sprekeler_synapse_hpc" ); } } // namespace nest diff --git a/models/vogels_sprekeler_connection.h b/models/vogels_sprekeler_connection.h new file mode 100644 index 0000000000..f0790bb6cf --- /dev/null +++ b/models/vogels_sprekeler_connection.h @@ -0,0 +1,313 @@ +/* + * vogels_sprekeler_connection.h + * + * This file is part of NEST. + * + * Copyright (C) 2004 The NEST Initiative + * + * NEST is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 2 of the License, or + * (at your option) any later version. + * + * NEST is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with NEST. If not, see . + * + */ + +#ifndef VOGELS_SPREKELER_CONNECTION_H +#define VOGELS_SPREKELER_CONNECTION_H + +/* BeginDocumentation + Name: vogels_sprekeler_synapse - Synapse type for symmetric spike-timing + dependent + plasticity with constant depression. + + Description: + vogels_sprekeler_synapse is a connector to create synapses with symmetric + spike time dependent plasticity and constant depression (as defined in [1]). + The learning rule is symmetric, i.e., the synapse is strengthened + irrespective of the order of the pre and post-synaptic spikes. Each + pre-synaptic spike also causes a constant depression of the synaptic weight + which differentiates this rule from other classical stdp rules. + + Parameters: + tau double - Time constant of STDP window, potentiation in ms + Wmax double - Maximum allowed weight + eta double - learning rate + alpha double - constant depression (= 2 * tau * target firing rate in + [1]) + + Transmits: SpikeEvent + + References: + [1] Vogels et al. (2011) Inhibitory Plasticity Balances Excitation and + Inhibition in Sensory Pathways and Memory Networks. + Science Vol. 334, Issue 6062, pp. 1569-1573 + DOI: 10.1126/science.1211095 + + FirstVersion: January 2016 + Author: Ankur Sinha + SeeAlso: synapsedict +*/ + +// C-header for math.h since copysign() is in C99 but not C++98 +#include +#include "connection.h" + +namespace nest +{ + +// connections are templates of target identifier type (used for pointer / +// target index addressing) +// derived from generic connection template +template < typename targetidentifierT > +class VogelsSprekelerConnection : public Connection< targetidentifierT > +{ + +public: + typedef CommonSynapseProperties CommonPropertiesType; + typedef Connection< targetidentifierT > ConnectionBase; + + /** + * Default Constructor. + * Sets default values for all parameters. Needed by GenericConnectorModel. + */ + VogelsSprekelerConnection(); + + + /** + * Copy constructor. + * Needs to be defined properly in order for GenericConnector to work. + */ + VogelsSprekelerConnection( const VogelsSprekelerConnection& ); + + // Explicitly declare all methods inherited from the dependent base + // ConnectionBase. + // This avoids explicit name prefixes in all places these functions are used. + // Since ConnectionBase depends on the template parameter, they are not + // automatically + // found in the base class. + using ConnectionBase::get_delay_steps; + using ConnectionBase::get_delay; + using ConnectionBase::get_rport; + using ConnectionBase::get_target; + + /** + * Get all properties of this connection and put them into a dictionary. + */ + void get_status( DictionaryDatum& d ) const; + + /** + * Set properties of this connection from the values given in dictionary. + */ + void set_status( const DictionaryDatum& d, ConnectorModel& cm ); + + /** + * Send an event to the receiver of this connection. + * \param e The event to send + * \param t_lastspike Point in time of last spike sent. + * \param cp common properties of all synapses (empty). + */ + void send( Event& e, + thread t, + double_t t_lastspike, + const CommonSynapseProperties& cp ); + + + class ConnTestDummyNode : public ConnTestDummyNodeBase + { + public: + // Ensure proper overriding of overloaded virtual functions. + // Return values from functions are ignored. + using ConnTestDummyNodeBase::handles_test_event; + port + handles_test_event( SpikeEvent&, rport ) + { + return invalid_port_; + } + }; + + void + check_connection( Node& s, + Node& t, + rport receptor_type, + double_t t_lastspike, + const CommonPropertiesType& ) + { + ConnTestDummyNode dummy_target; + + ConnectionBase::check_connection_( dummy_target, s, t, receptor_type ); + + t.register_stdp_connection( t_lastspike - get_delay() ); + } + + void + set_weight( double_t w ) + { + weight_ = w; + } + +private: + double_t + facilitate_( double_t w, double_t kplus ) + { + double_t new_w = std::abs( w ) + ( eta_ * kplus ); + return copysign( new_w < std::abs( Wmax_ ) ? new_w : Wmax_, Wmax_ ); + } + + double_t + depress_( double_t w ) + { + double_t new_w = std::abs( w ) - ( alpha_ * eta_ ); + return copysign( new_w > 0.0 ? new_w : 0.0, Wmax_ ); + } + + // data members of each connection + double_t weight_; + double_t tau_; + double_t alpha_; + double_t eta_; + double_t Wmax_; + double_t Kplus_; +}; + + +/** + * Send an event to the receiver of this connection. + * \param e The event to send + * \param t The thread on which this connection is stored. + * \param t_lastspike Time point of last spike emitted + * \param cp Common properties object, containing the stdp parameters. + */ +template < typename targetidentifierT > +inline void +VogelsSprekelerConnection< targetidentifierT >::send( Event& e, + thread t, + double_t t_lastspike, + const CommonSynapseProperties& ) +{ + // synapse STDP depressing/facilitation dynamics + double_t t_spike = e.get_stamp().get_ms(); + // t_lastspike_ = 0 initially + + // use accessor functions (inherited from Connection< >) to obtain delay and + // target + Node* target = get_target( t ); + double_t dendritic_delay = get_delay(); + + // get spike history in relevant range (t1, t2] from post-synaptic neuron + std::deque< histentry >::iterator start; + std::deque< histentry >::iterator finish; + target->get_history( + t_lastspike - dendritic_delay, t_spike - dendritic_delay, &start, &finish ); + + // presynaptic neuron j, post synaptic neuron i + // Facilitation for each post synaptic spike + // Wij = Wij + eta*xj + double_t minus_dt; + while ( start != finish ) + { + minus_dt = t_lastspike - ( start->t_ + dendritic_delay ); + ++start; + if ( minus_dt == 0 ) + continue; + weight_ = facilitate_( weight_, Kplus_ * std::exp( minus_dt / tau_ ) ); + } + + // For pre-synaptic spikes + // Wij = Wij + eta(xi - alpha) + // Facilitation and constant depression + // Getting kvalue at required time already for deferred processing, so no + // need to transform it to the current time, and so, no exponential required + weight_ = + facilitate_( weight_, target->get_K_value( t_spike - dendritic_delay ) ); + weight_ = depress_( weight_ ); + + e.set_receiver( *target ); + e.set_weight( weight_ ); + // use accessor functions (inherited from Connection< >) to obtain delay in + // steps and rport + e.set_delay( get_delay_steps() ); + e.set_rport( get_rport() ); + e(); + + // exponential part for the decay, addition of one for each spike + Kplus_ = Kplus_ * std::exp( ( t_lastspike - t_spike ) / tau_ ) + 1.0; +} + + +template < typename targetidentifierT > +VogelsSprekelerConnection< targetidentifierT >::VogelsSprekelerConnection() + : ConnectionBase() + , weight_( 0.5 ) + , tau_( 20.0 ) + , alpha_( 0.12 ) + , eta_( 0.001 ) + , Wmax_( 1.0 ) + , Kplus_( 0.0 ) +{ +} + +template < typename targetidentifierT > +VogelsSprekelerConnection< targetidentifierT >::VogelsSprekelerConnection( + const VogelsSprekelerConnection< targetidentifierT >& rhs ) + : ConnectionBase( rhs ) + , weight_( rhs.weight_ ) + , tau_( rhs.tau_ ) + , alpha_( rhs.alpha_ ) + , eta_( rhs.eta_ ) + , Wmax_( rhs.Wmax_ ) + , Kplus_( rhs.Kplus_ ) +{ +} + +template < typename targetidentifierT > +void +VogelsSprekelerConnection< targetidentifierT >::get_status( + DictionaryDatum& d ) const +{ + ConnectionBase::get_status( d ); + def< double_t >( d, names::weight, weight_ ); + def< double_t >( d, "tau", tau_ ); + def< double_t >( d, "alpha", alpha_ ); + def< double_t >( d, "eta", eta_ ); + def< double_t >( d, "Wmax", Wmax_ ); + def< double_t >( d, "Kplus", Kplus_ ); + def< long_t >( d, names::size_of, sizeof( *this ) ); +} + +template < typename targetidentifierT > +void +VogelsSprekelerConnection< targetidentifierT >::set_status( + const DictionaryDatum& d, + ConnectorModel& cm ) +{ + ConnectionBase::set_status( d, cm ); + updateValue< double_t >( d, names::weight, weight_ ); + updateValue< double_t >( d, "tau", tau_ ); + updateValue< double_t >( d, "alpha", alpha_ ); + updateValue< double_t >( d, "eta", eta_ ); + updateValue< double_t >( d, "Wmax", Wmax_ ); + updateValue< double_t >( d, "Kplus", Kplus_ ); + + // check if weight_ and Wmax_ has the same sign + if ( not( ( ( weight_ >= 0 ) - ( weight_ < 0 ) ) + == ( ( Wmax_ >= 0 ) - ( Wmax_ < 0 ) ) ) ) + { + throw BadProperty( "Weight and Wmax must have same sign." ); + } + + if ( not( Kplus_ >= 0 ) ) + { + throw BadProperty( "State Kplus must be positive." ); + } +} +} // of namespace nest + +#endif // of #ifndef VOGELS_SPREKELER_CONNECTION_H diff --git a/pynest/nest/tests/test_vogels_sprekeler_synapse.py b/pynest/nest/tests/test_vogels_sprekeler_synapse.py new file mode 100644 index 0000000000..9bc67c9646 --- /dev/null +++ b/pynest/nest/tests/test_vogels_sprekeler_synapse.py @@ -0,0 +1,236 @@ +# -*- coding: utf-8 -*- +# +# test_vogels_sprekeler_synapse.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +# This script tests the vogels_sprekeler_synapse in NEST. + +import nest +import unittest +from math import exp + + +@nest.check_stack +class VogelsSprekelerConnectionTestCase(unittest.TestCase): + + """Check vogels_sprekeler_synapse model properties.""" + + def setUp(self): + """Set up the test.""" + nest.set_verbosity('M_WARNING') + nest.ResetKernel() + + # settings + self.dendritic_delay = 1.0 + self.decay_duration = 5.0 + self.synapse_model = "vogels_sprekeler_synapse" + self.syn_spec = { + "model": self.synapse_model, + "delay": self.dendritic_delay, + "weight": 5.0, + "eta": 0.001, + "alpha": 0.1, + "tau": 20., + "Kplus": 0.0, + "Wmax": 15., + } + + # setup basic circuit + self.pre_neuron = nest.Create("parrot_neuron") + self.post_neuron = nest.Create("parrot_neuron") + nest.Connect(self.pre_neuron, self.post_neuron, + syn_spec=self.syn_spec) + + def generateSpikes(self, neuron, times): + """Trigger spike to given neuron at specified times.""" + delay = self.dendritic_delay + gen = nest.Create("spike_generator", 1, + {"spike_times": [t-delay for t in times]}) + nest.Connect(gen, neuron, syn_spec={"delay": delay}) + + def status(self, which): + """Get synapse parameter status.""" + stats = nest.GetConnections(self.pre_neuron, + synapse_model=self.synapse_model) + return nest.GetStatus(stats, [which])[0][0] + + def decay(self, time, Kvalue): + """Decay variables.""" + Kvalue *= exp(- time / self.syn_spec["tau"]) + return Kvalue + + def facilitate(self, w, Kplus): + """Facilitate weight.""" + new_w = w + (self.syn_spec['eta'] * Kplus) + return new_w if (new_w / self.status('Wmax') < 1.0) else \ + self.syn_spec['Wmax'] + + def depress(self, w): + """Depress weight.""" + new_w = w - (self.syn_spec['alpha'] * self.syn_spec['eta']) + return new_w if (new_w / self.status('Wmax') > 0.0) else 0 + + def assertAlmostEqualDetailed(self, expected, given, message): + """Improve assetAlmostEqual with detailed message.""" + messageWithValues = ("%s (expected: `%s` was: `%s`" % (message, + str(expected), + str(given))) + self.assertAlmostEqual(given, expected, msg=messageWithValues) + + def test_badPropertiesSetupsThrowExceptions(self): + """Check that exceptions are thrown when setting bad parameters.""" + def setupProperty(property): + bad_syn_spec = self.syn_spec.copy() + bad_syn_spec.update(property) + nest.Connect(self.pre_neuron, self.post_neuron, + syn_spec=bad_syn_spec) + + def badPropertyWith(content, parameters): + self.assertRaisesRegexp(nest.NESTError, "BadProperty(.+)" + + content, setupProperty, parameters) + + badPropertyWith("Kplus", {"Kplus": -1.0}) + + def test_varsZeroAtStart(self): + """Check that pre and post-synaptic variables are zero at start.""" + self.assertAlmostEqualDetailed(0.0, self.status("Kplus"), + "Kplus should be zero") + + def test_preVarsIncreaseWithPreSpike(self): + """ + Check that pre-synaptic variable, Kplus increase after each + pre-synaptic spike. + """ + self.generateSpikes(self.pre_neuron, [2.0]) + + Kplus = self.status("Kplus") + + nest.Simulate(20.0) + self.assertAlmostEqualDetailed(Kplus + 1.0, self.status("Kplus"), + "Kplus should have increased by 1") + + def test_preVarsDecayAfterPreSpike(self): + """ + Check that pre-synaptic variables Kplus decay after each + pre-synaptic spike. + """ + self.generateSpikes(self.pre_neuron, [2.0]) + self.generateSpikes(self.pre_neuron, + [2.0 + self.decay_duration]) # trigger computation + + Kplus = self.decay(self.decay_duration, 1.0) + Kplus += 1.0 + + nest.Simulate(20.0) + self.assertAlmostEqualDetailed(Kplus, self.status("Kplus"), + "Kplus should have decay") + + def test_preVarsDecayAfterPostSpike(self): + """ + Check that pre-synaptic variables Kplus decay after each post-synaptic + spike. + """ + self.generateSpikes(self.pre_neuron, [2.0]) + self.generateSpikes(self.post_neuron, [3.0, 4.0]) + self.generateSpikes(self.pre_neuron, + [2.0 + self.decay_duration]) # trigger computation + + Kplus = self.decay(self.decay_duration, 1.0) + Kplus += 1.0 + + nest.Simulate(20.0) + self.assertAlmostEqualDetailed(Kplus, self.status("Kplus"), + "Kplus should have decay") + + def test_weightChangeWhenPrePostSpikes(self): + """Check that weight changes whenever a pre-post spike pair happen.""" + self.generateSpikes(self.pre_neuron, [2.0]) + self.generateSpikes(self.post_neuron, [4.0]) + self.generateSpikes(self.pre_neuron, [6.0]) # trigger computation + + print("") + weight = self.status("weight") + Kplus = self.status("Kplus") + Kminus = 0. + + Kplus = self.decay(2.0, Kplus) + # first pre-synaptic spike + weight = self.facilitate(weight, Kminus) + weight = self.depress(weight) + Kplus += 1.0 + + # Resultant postspike at 3.0ms (because we're using parrot neurons, the + # prespike causes a postspike too + Kminus += 1.0 + + # Planned postspike at 4.0ms + Kminus = self.decay(1.0, Kminus) + Kminus += 1.0 + + # next pre-synaptic spike + # first postspike in history + dt = 2.0 + Kplus_temp1 = self.decay(2.0, Kplus) + weight = self.facilitate(weight, Kplus_temp1) + # second postspike in history + dt = 3.0 + Kplus_temp2 = self.decay(3.0, Kplus) + weight = self.facilitate(weight, Kplus_temp2) + + Kminus = self.decay(1.0, Kminus) + weight = self.facilitate(weight, Kminus) + weight = self.depress(weight) + + Kplus = self.decay(4.0, Kplus) + + # The simulation using Nest + nest.Simulate(20.0) + self.assertAlmostEqualDetailed(weight, self.status("weight"), + "weight should have increased") + + def test_maxWeightStaturatesWeight(self): + """Check that setting maximum weight property keep weight limited.""" + limited_weight = self.status("weight") + 1e-10 + conn = nest.GetConnections(target=self.post_neuron, + source=self.pre_neuron) + # disable depression to make it get to max weight + # increase eta to cause enough facilitation + nest.SetStatus(conn, "Wmax", limited_weight) + nest.SetStatus(conn, "eta", 5.) + nest.SetStatus(conn, "alpha", 0.) + + self.generateSpikes(self.pre_neuron, [2.0]) + self.generateSpikes(self.post_neuron, [3.0]) + self.generateSpikes(self.pre_neuron, [4.0]) # trigger computation + + self.assertAlmostEqualDetailed(limited_weight, + self.status("weight"), + "weight should have been limited") + + +def suite(): + return unittest.makeSuite(VogelsSprekelerConnectionTestCase, "test") + + +def run(): + runner = unittest.TextTestRunner(verbosity=2) + runner.run(suite()) + +if __name__ == "__main__": + run()